The 3D version of this simulator runs a kinetic Monte Carlo process: it picks random transitions for each individual surface site and lets the average emerge from noise. This 2D counterpart solves the exact same birth–death reaction network with a completely different, independently computed method — it never simulates a single site. Instead it integrates the master equation for pₙ(t), the population fraction of sites that currently hold exactly n adsorbed layers, using 4th-order Runge–Kutta:
dp₀/dt = −Cx·p₀ + 1·p₁
dp₁/dt = Cx·p₀ − (x+1)·p₁ + 1·p₂
dpₙ/dt = x·pₙ₋₁ − (x+1)·pₙ + 1·pₙ₊₁ (2 ≤ n < Nmax)
dp_Nmax/dt = x·p_Nmax₋₁ − 1·p_Nmax
x = P/P₀ is the relative pressure and C the BET affinity constant — the same two control parameters and the same four transition rates (Cx, 1, x, 1) as the 3D engine's per-site trials. Instead of counting how many sites hold n layers out of a finite grid, this model tracks the exact probability pₙ(t) directly as a continuous quantity — a deterministic relaxation with no sampling noise at all. The bars in the main panel are pₙ(t) itself; the thin flowing arrows between them show the real net probability current Jₙ = (rate n→n+1)·pₙ − (rate n+1→n)·pₙ₊₁ that the master equation computes at every step. At long times this exactly reproduces the closed-form BET isotherm:
V/Vₘ = Σ n·pₙ → C·x / [ (1−x)·(1−x+C·x) ]
- P/P₀ — relative pressure; higher values pack on more layers and slow the relaxation toward steady state.
- C — surface affinity; low C (~1) gives gradual uptake, high C (~100+) a sharp monolayer step.
- Layer capacity Nmax — the truncation of the infinite chain of layer-states, mirroring the 3D engine's own hard MAX_LAYERS=9 cap. Push P/P₀ high with a low Nmax and the solved n̄ will visibly saturate below the diverging theory value — the same real truncation effect the 3D grid shows once a site's stack hits its ceiling, not a bug in either model.
- Integration speed — how many RK4 steps run per animation frame, i.e. how fast simulated time advances toward the steady state.