Gillespie direct method with horizon and budget statuses

Inputs. A finite state set $E\subseteq {\mathbb{Z}}^d$, a finite-channel reaction propensity system from Finite-channel reaction propensity systems and their generator, an initial state $x_0\in E$, a horizon $T>0$, and an optional event budget $B\in {\mathbb{N}}_0\cup \{\mathrm{\infty }\}$, defaulting to $\mathrm{\infty }$, where ${\mathbb{N}}_0:=\{0,1,2,\dots \}$. The random source is a sequence $(E_n,U_n{)}_{n\ge 1}$ with $E_n\sim \mathrm{Exp}(1)$ and $U_n$ uniform on $[0,1]$, all coordinates mutually independent. The countable product construction for these input laws uses Countable Choice and Dependent Choice by Assuming countable and dependent choice, countable products of arbitrary probability spaces; the exponential law uses Countable Choice by Exponential and Weibull survival specifications define probability laws. For the uniform law, Assuming countable choice, $\mathcal{L}({\mathbb{R}}^n)$ is a sigma-algebra containing every elementary set and ${\lambda }_n$ is a complete measure extending elementary volume makes the set function in Lebesgue measurable sets, the family $\mathcal{L}({\mathbb{R}}^n)$, and the restricted set function ${\lambda }_n$ a measure, and A box in ${\mathbb{R}}^n$ with parameters $a_i\le b_i$ is Lebesgue measurable of measure ${\prod }_{i<n}(b_i-a_i)$, whichever of its faces are included gives its restriction to $[0,1]$ mass one and interval probabilities equal to interval lengths.

Procedure. Set $t=0$, $x=x_0$, $n=0$, and record the initial state. While $t<T$:

1. Compute $\lambda (x)={\sum }_{j=1}^R{a}_j(x)$. If $\lambda (x)=0$, keep the state equal to $x$ through $T$ and return status absorbed.
2. If $n=B$, return the path prefix through $(t,x)$ with status budget-exhausted; this status does not assert that the horizon was reached.
3. Set $\tau =E_{n+1}\mathrm{/}\lambda (x)$. If $t+\tau \ge T$, keep the state equal to $x$ through $T$ and return status horizon.
4. Set $J$ to the least index with $U_{n+1}<\lambda (x{)}^{-1}{\sum }_{j=1}^J{a}_j(x)$. If the null event $U_{n+1}=1$ occurs, use the largest index with $a_j(x)>0$. This rule never selects a zero-rate channel.
5. Set $t\leftarrow t+\tau$, $x\leftarrow x+{\nu }_J$, append the marked event $(t,J,x)$, and set $n\leftarrow n+1$.

If the loop reaches $t=T$, return the path through $T$ with status horizon. Each event uses a linear scan of the $R$ propensities and a linear cumulative channel search; no update-tree acceleration is part of this procedure. If $B<\mathrm{\infty }$ and the budget expires before $T$, the returned prefix is incomplete.

The procedure assumes an ideal independent random source with the displayed laws. It makes no claim that a finite-state pseudorandom generator supplies independent exact exponential or uniform variables.