The Model

The model

The balances follow Ruthven (1984), and the equation numbers cited on this page (Ruthven Eq. 7.11 and so on) refer to that book; the dynamic-column-breakthrough framing follows Wilkins, Rajendran & Farooq (2021). The equations below are what the solver integrates, in dimensionless form, with

$$Z = z/L, \quad \tau = t\,v_0/L, \quad P \equiv P/P_0, \quad T \equiv T/T_0, \quad v \equiv v/v_0, \quad x_i = q_i/q_{s,0}$$

Governing equations

Component 1, gas phase

$$\frac{\partial y_1}{\partial\tau}=\frac{T}{Pe_a P}\frac{\partial}{\partial Z}\left(\frac{P}{T}\frac{\partial y_1}{\partial Z}\right)-\frac{T}{P}\frac{\partial}{\partial Z}\left(\frac{y_1 P v}{T}\right)-\psi\frac{T}{P}\frac{\partial x_1}{\partial\tau}-\frac{y_1}{P}\frac{\partial P}{\partial\tau}+\frac{y_1}{T}\frac{\partial T}{\partial\tau}$$

Overall mass (total pressure)

$$\mathcal{R}_P=-T\frac{\partial}{\partial Z}\left(\frac{Pv}{T}\right)-\psi T\left(\frac{\partial x_1}{\partial\tau}+\frac{\partial x_2}{\partial\tau}\right)$$

Column energy

$$\mathcal{R}_T=\omega_1\frac{\partial^2 T}{\partial Z^2}-\omega_2\frac{\partial (Pv)}{\partial Z}-\omega_3 T\left(\frac{\partial x_1}{\partial\tau}+\frac{\partial x_2}{\partial\tau}\right)+\sigma_1\frac{\partial x_1}{\partial\tau}+\sigma_2\frac{\partial x_2}{\partial\tau}-\omega_4\left(T-T_w\right)$$

Pressure–temperature decoupling. The two balances above share the same time derivatives, so they are solved simultaneously at each step:

$$\frac{\partial P}{\partial\tau}=\frac{\mathcal{R}_P+(P/T)\mathcal{R}_T}{1+P\omega_2/T},\qquad \frac{\partial T}{\partial\tau}=\frac{\mathcal{R}_T-\omega_2\mathcal{R}_P}{1+P\omega_2/T}$$

Wall energy

$$\frac{\partial T_w}{\partial\tau}=\pi_1\frac{\partial^2 T_w}{\partial Z^2}+\pi_2\left(T-T_w\right)-\pi_3\left(T_w-T_{amb}\right)$$

Mass transfer (linear driving force)

$$\frac{\partial x_i}{\partial\tau}=K_i\frac{L}{v_0}\left(x_i^{\ast}-x_i\right)$$

Momentum (Ergun)

$$A\,v+B\,v\left|v\right|=-\frac{\partial P}{\partial Z}$$

Dimensionless groups

Every group is formed from the run's own inputs and reported per run in the Computed params tab. $\Lambda$ is the shared volumetric heat capacity of the column energy balance:

$$\Lambda=\rho_s C_{p,s}+\rho_s C_{p,ads}\,q_{s,0}\left(x_1+x_2\right)$$

Group	Definition	Meaning
$\psi$	$\dfrac{R T_0\rho_s q_{s,0}}{P_0}\dfrac{1-\varepsilon}{\varepsilon}$	Adsorbed-phase coupling: converts solid uptake into the equivalent gas-phase mole change
$Pe_a$	$\dfrac{v_0 L}{D_L}$, $D_L=\gamma_1 D_m+\tfrac{1}{2}d_p v_0$	Axial mass Péclet number; strength of the dispersion term, evaluated per cell. $\gamma_1$ depends on the selected correlation — see below
$Pe_h$	$\dfrac{\varepsilon v_0 L\,\rho_g C_{p,gas}}{k_z}$, $\rho_g=\dfrac{P}{RT}$	Axial thermal Péclet number, evaluated at the local pressure and temperature; used in the inlet half-cell temperature boundary condition
$\omega_1$	$\dfrac{k_z}{(1-\varepsilon)v_0 L\,\Lambda}$	Axial conduction through the bed
$\omega_2$	$\dfrac{C_{p,gas}P_0}{R T_0}\dfrac{\varepsilon}{1-\varepsilon}\dfrac{1}{\Lambda}$	Compression and convective enthalpy carried by the flowing gas
$\omega_3$	$\dfrac{C_{p,ads}q_{s,0}\rho_s}{\Lambda}$	Sensible heat of the adsorbed phase as loading changes
$\omega_4$	$\dfrac{2h_{in}L}{r_{in}v_0(1-\varepsilon)\Lambda}$	Bed-to-wall heat transfer
$\sigma_i$	$-\dfrac{\Delta H_i\,q_{s,0}\rho_s}{T_0\,\Lambda}$	Heat of adsorption of component $i$; negative $\Delta H_i$ gives an exotherm on adsorption
$\pi_1$	$\dfrac{k_w}{\rho_w C_{p,w}v_0 L}$	Axial conduction along the wall
$\pi_2$	$\dfrac{2r_{in}h_{in}L}{(r_{out}^2-r_{in}^2)\rho_w C_{p,w}v_0}$	Bed-to-wall transfer as seen by the wall
$\pi_3$	$\dfrac{2r_{out}h_{out}L}{(r_{out}^2-r_{in}^2)\rho_w C_{p,w}v_0}$	Wall-to-ambient loss; with $\pi_2$ this sets the thermal time constant

Gas heat capacities use Shomate polynomials per species. All heat capacities are on a molar basis ($C_{p,gas}$, $C_{p,ads}$ in J/mol·K). A mass-basis value entered in these fields gives a wrong thermal response.

Axial dispersion

Every available correlation has the same form, $D_L=\gamma_1 D_m+\tfrac{1}{2}d_p v_0$, and differs only in $\gamma_1$. The correlations are listed under Thermal & Transport, each evaluated at the current operating point, and selected by clicking one.

Correlation	$\gamma_1$	When it applies
Basic correlation, Ruthven Eq. 7.6 (default)	$0.70$	$D_L=0.7D_m+0.5\,v_0d_p$, the book's basic correlation; a lower bound for a porous, strongly adsorbing system
Wakao, Ruthven Eq. 7.11	$20/\varepsilon$	Porous particle carrying a strongly favourable isotherm. Ruthven is explicit that dispersion here is far larger than the nonporous correlations give
Edwards–Richardson, Ruthven Eq. 7.9	resolves $1/Pe'$ directly	Intermediate $Re_p$, where Ruthven Eq. 7.6 under-predicts dispersion
Wicke, Ruthven Eq. 7.8	$0.45+0.55\varepsilon$	Voidage-corrected form of Ruthven Eq. 7.6
Custom $D_L$	—	Your own value in m²/s, e.g. fitted to a tracer pulse, held constant along the bed; prefilled with the selected correlation's value

Ruthven's basic correlation is the default. Wakao (Ruthven Eq. 7.11) is derived for a porous particle with a strongly favourable isotherm in laminar flow, the regime the bundled presets are in; at typical conditions it gives $D_L$ about 2.6 times Ruthven Eq. 7.6, so the front is correspondingly broader. Select it deliberately, on a case you are prepared to re-fit.

Cases saved before 2026.09v3 carry no dispersion selection and load with the default, Ruthven Eq. 7.6 — the only correlation available then — so an old case reproduces its original numbers. The Logs tab records it.

One caveat the book attaches to the eddy term rather than to $\gamma_1$: by Ruthven Eq. 7.10, below $R_p<0.15$ cm the limiting Péclet number is particle-size dependent rather than 2.0, so $\gamma_2=0.5$ is itself outside its stated range for millimetre pellets.

Momentum balance

Ergun by default — viscous (Blake–Kozeny) plus inertial (Burke–Plummer). The inertial term can be switched off to give pure Darcy, $v = -A^{-1}\partial P/\partial Z$.

Ergun is the default. Darcy is appropriate where $Re_p$ (Computed params) is small; when you use it, compare against Ergun before trusting the result. Where the inertial term is negligible the two agree closely; a significant difference shows that the inertial term matters for your case.

Numerics

• Discretisation: finite volume, $N$ cells, WENO reconstruction at cell faces, upwinded by the sign of the local face velocity. Where $v_{j+1/2}\ge0$ the face value is $W(u_{j-1},u_j,u_{j+1})$ as usual; where the flow reverses it is the mirrored $W(u_{j+2},u_{j+1},u_j)$. Some runs briefly reverse the flow near the start, and direction-aware upwinding keeps such runs fast and stable. Where the flow never reverses, the result is identical to forward-only WENO.
• System: $6N$ coupled ODEs.
• Integrator: SUNDIALS CVODE (implicit BDF) through scikit-sundae, with a sparse SuperLU_MT linear solver and an explicitly supplied Jacobian sparsity pattern.
• Limits enforced by the backend: $N_{grid} \in [5, 100]$, $t_{run} \in [1, 10^6]$ s, two concurrent solves, 900 s wall-clock cap per solve.

The finite-volume scheme is conservative: what leaves one cell enters the next, so moles are conserved at any cell count, and the mole-conservation check measures the time integration rather than the discretisation. The default is $N=30$. It is not grid-converged for the front shape: on the default column $t_{5\%}$ comes out 2–3 % earlier than at $N=150$ (numerical dispersion), while $t_{50\%}$ moves by about 1 %. Use a finer grid when fitting to measurements.

On a coarse grid the outlet flow also shows a small periodic ripple before breakthrough, with a period equal to the time the front takes to cross one cell: each cell takes up most of its adsorbate, and releases its heat, as the front passes through it. On the default case the ripple is about 1.3 % of the flow at $N=30$, 0.7 % at $N=60$ and 0.3 % at $N=100$; solver tolerance has no effect on it. It does not affect breakthrough times or the mass balance.

How results are checked

The engine includes a built-in self-test of 139 checks covering isotherm limiting behaviour, IAST consistency, heat-of-adsorption signs, dimensionless group formation, and mass and energy closure; all pass for every published build.

In addition, every run reports two checks. The loading the bed accumulated is set against the mixture equilibrium loading the isotherm gives at feed conditions, and, separately, the integrated inlet-minus-outlet flow is set against the change in bed inventory, per component. A discretisation or integration error would have to shift both by the same amount to go unnoticed. These are consistency checks on the numerics; they do not establish that the model is appropriate for your system.

The default Zeolite 13X case uses the isotherm parameters measured by Wilkins & Rajendran (2019). To compare a run against a measurement, load the data through Experimental validation; the app reports RMSE, MAE and R².