# Finite-ring TASEP: probabilities, trajectories and a Bethe mode This standalone calculation treats continuous-time right-moving exclusion on four periodic sites with two particles and homogeneous rate r>0 per eligible jump. It tests the probability law, one reconstructed Bethe mode and sampling conventions. It does not establish general Bethe completeness, large-system scaling or any open-boundary result. ## Reproduce The saved run used Python 3.9.6 and NumPy 2.0.2. The pinned NumPy release supports Python 3.9–3.12. Download `experiment.py`, `inputs.json`, `results.json`, `requirements.txt` and this file into one directory. From there: ```sh python3 -m venv .venv .venv/bin/python -m pip install -r requirements.txt .venv/bin/python experiment.py --check ``` On Windows use `.venv\Scripts\python.exe`. With NumPy installed in the repository: ```sh python3 public/computations/tasep/experiment.py --check ``` Expected output: ```text TASEP: generator, current, spectrum, Bethe mode, uniformization, trajectories, negative controls and saved-result checks passed. ``` `--check` recomputes the experiment without changing files. It requires the complete input object to match exactly, compares integer trajectory counts exactly, and compares numerical results with absolute tolerance 2e-12 and relative tolerance 1e-9. Environment strings are informational. No flag prints fresh JSON; `--write` intentionally replaces `results.json` after all checks pass. Paths are relative to the script. Routine site builds do not run it. Explore changed inputs in a separate copy. The fixed PCG64 seed and pinned NumPy support reproduction; the statistical tests below provide a separate uncertainty interpretation. ## Generator and exact checks The basis is the lexicographic occupied-site list `12,13,14,23,24,34`. Probability columns evolve by `p'=Qp`. Each off-diagonal entry has orientation `Q(destination,source)`; diagonal entries are minus total exit rates. The code constructs integer `Q/r` by bit-mask moves, including the seam 4→1. It compares against an independently specified matrix from these transitions, each rate r: ```text 12 -> 13 13 -> 14 and 23 14 -> 24 23 -> 24 24 -> 12 and 34 34 -> 13. ``` Exact integer checks establish nonnegative off-diagonal rates, zero column sums and zero row sums. The latter proves the uniform stationary distribution in this sector. Exact Fraction arithmetic computes every characteristic coefficient by the Faddeev-LeVerrier trace recurrence, independently of numerical eigenvalue finding: ```text det(q I-Q/r) = q(q+1)^2(q+3)(q^2+3q+4) descending coefficients = [1,8,26,44,37,12,0]. ``` NumPy diagonalization is only a cross-check against the analytic spectrum `{0,-1,-1,-3,(-3+i sqrt(7))/2,(-3-i sqrt(7))/2}` in units r. Its largest saved eigenvalue discrepancy is 4.00e-15 (rounded upward). Counting allowed jumps in the six equally probable configurations gives current r/3 across each directed bond, total jump rate 4r/3 and rate 2r/3 per particle. These are exact rational results. They agree with the finite-ring expression `j/r=M(N-M)/(N(N-1))`, not the uncorrelated-density approximation rho(1-rho). ## Uniformization and its error The deterministic initial distribution is concentrated at occupied sites 1,2. Saved rate r=1 and dimensionless observation times rt=0.25,1,2 are explicit in `inputs.json`. Changing r while holding rt fixed rescales physical time. The maximal exit rate is nu=2r. Set `P=I+Q/nu`; P is column-stochastic. Uniformization evaluates ```text p(t) = exp(-nu*t) sum_(k>=0) (nu*t)^k/k! P^k p(0). ``` No eigenvector decomposition enters this calculation. With mu=nu*t, the Poisson weights satisfy `w_(k+1)=w_k*mu/(k+1)`. After retaining degree K, the omitted tail is at most `w_(K+1)/(1-mu/(K+2))` when K+2>mu. This follows by bounding the remaining weight ratios by a geometric series. The code stops below 1e-13 and does not renormalize the partial sum. In exact arithmetic its L1 error is exactly the missing Poisson mass, because every `P^k p(0)` is a normalized nonnegative vector. Binary64 roundoff is separate; deterministic checks allow 2e-12 in addition to the analytic tail. | rt | Last included degree | Analytic tail bound | Analytic mode L1 error | |---:|---:|---:|---:| | 0.25 | 12 | 1.2331e-14 | 1.2435e-14 | | 1 | 19 | 6.4470e-14 | 6.4060e-14 | | 2 | 26 | 3.5352e-14 | 3.5528e-14 | Positive entries here are rounded upward. A floating-point mass deficit can exceed the analytic tail slightly through rounding; this is not a failure of the truncation argument. The implementation is intended for this short finite-time experiment, not very large Poisson means where exp(-mu) underflows. ## Independent paths and statistical checks The trajectory routine never reads Q or calls the bit-mask move generator. It stores particle positions, finds vacant successors, draws an exponential waiting time of rate r times the number of available moves, then chooses an available move uniformly. It records the configuration at prescribed times between events. There are 50,000 independent paths, all initially 12, generated with `Generator(PCG64(20261003))`. At rt=1 the saved counts are `(19365,11816,6986,6884,3878,1071)`. Their frequency vector is `(0.38730,0.23632,0.13972,0.13768,0.07756,0.02142)`. The deterministic vector, rounded, is `(0.388883,0.240126,0.136688,0.136688,0.076612,0.021003)`. Each fixed-state count at one time is binomial under the independent-path model. `binomial_standard_errors_from_deterministic_p` reports `sqrt(p(1-p)/n)` using the deterministic probability. At rt=0.25,1,2 the largest absolute errors are 0.00140748,0.00380552,0.00492868; the largest errors divided by this standard error are 1.29179,1.99209,2.63576. No convergence rate is inferred from these three times or from one seed. `wilson_95_pointwise_intervals` are approximate 95% score intervals. Starting with `n(p_hat-p)^2 <= z^2 p(1-p)`, solve the quadratic for p. With `z=1.959963984540054`, the center and radius are ```text center = (p_hat+z^2/(2n))/(1+z^2/n) radius = z*sqrt(p_hat*(1-p_hat)/n+z^2/(4n^2))/(1+z^2/n). ``` This is Wilson's score inversion [Wilson 1927 below]. They are pointwise, approximate intervals, not a simultaneous guarantee. At rt=1 the probabilities for states 13 and 14 narrowly miss their intervals in this recorded run. These misses are retained and are not themselves the acceptance criterion. For the acceptance criterion, Hoeffding's Theorem 1, equation (2.3), applied separately to Bernoulli indicators and complements gives `Pr(abs(p_hat-p)>=delta) <= 2 exp(-2*n*delta^2)` [Hoeffding 1963 below]. A union bound over m=18 state/time comparisons chooses `delta=sqrt(log(2*m/alpha)/(2*n))=0.013190538084710716`, with alpha=1e-6. This bounds the probability of at least one larger sampling discrepancy under the independent-trajectory model. Neither state counts within one multinomial sample nor times on the same trajectory must be independent for this union bound. The code allows the deterministic tail and rounding in addition. This is deliberately more conservative than the descriptive pointwise intervals. ## Bethe mode and a separate probability initial condition For x