---
title: McNabb–Foster Equations in Hydrogen Diffusion
url: https://www.emergentmind.com/topics/mcnabb-foster-equations
type: topic
---

# McNabb–Foster Equations in Hydrogen Diffusion

The McNabb–Foster equations are a generalized set of equations for hydrogen diffusion in the presence of saturable traps. In contemporary arXiv treatments of hydrogen transport in steels, they provide the governing framework for coupled lattice diffusion, trap occupancy evolution, permeation, room-temperature desorption, and thermal desorption spectroscopy (TDS). The framework is explicitly kinetic rather than purely local-equilibrium, and recent studies emphasize both its forward-modeling power and a central inverse-problem limitation: electrochemical permeation (EP) and TDS data can often be reproduced by multiple combinations of trap density and trapping parameters, so individual parameters are not necessarily uniquely identifiable [2507.15711], [2504.17532].

## 1. Governing equations and conservation structure

A recent finite-element treatment describes the McNabb–Foster framework as the “generalised set of equations for hydrogen diffusion in the presence of saturable traps” [2507.15711]. In that formulation, the lattice hydrogen concentration \(C_L\) evolves according to  
\[
\frac{\partial C_L}{\partial t} = D_L \nabla^2 C_L - \sum_{i=1}^{n} N_i \frac{\partial \theta_i}{\partial t},
\]
while each trap family \(i\) satisfies  
\[
\frac{\partial \theta_i}{\partial t} = K_i C_L (1-\theta_i) - \lambda_i \theta_i.
\]

Here \(\theta_i=C_i/N_i\) is the fractional occupancy of trap family \(i\), \(N_i\) is the trap-site density, \(K_i\) is the trapping rate constant, and \(\lambda_i\) is the detrapping rate constant. The total hydrogen concentration may be written as  
\[
C_{\text{tot}} = C_L + \sum_{i=1}^{n} N_i\theta_i,
\]
so the lattice equation is equivalently a conservation statement,
\[
\frac{\partial}{\partial t}\left(C_L+\sum_{i=1}^{n}N_i\theta_i\right)=D_L\nabla^2 C_L.
\]

In fitting examples for ferritic steels containing \((\mathrm{Ti},\mathrm{Cr})\mathrm{C}\) particles, the general system was reduced to a single reversible saturable trap,
\[
\frac{\partial C_L}{\partial t} = D_L \nabla^2 C_L - N \frac{\partial \theta}{\partial t},
\qquad
\frac{\partial \theta}{\partial t} = K C_L (1-\theta)-\lambda\theta.
\]

This structure distinguishes the McNabb–Foster description from a purely Fickian effective-diffusivity treatment. It retains explicit exchange between mobile lattice hydrogen and trapped hydrogen, and therefore represents transient storage and release rather than absorbing all trap effects into a single apparent diffusivity [2507.15711].

## 2. Kinetic parameters, occupancies, and related formulations

The same recent treatment gives Arrhenius forms for the trapping and detrapping rate constants,
\[
K_i = \nu \, N_L \exp\!\left(-\frac{E_{K,i}}{RT}\right),
\qquad
\lambda_i = \nu \exp\!\left(-\frac{E_{\lambda,i}}{RT}\right),
\]
with trap binding energy
\[
\Delta E_{T,i} = E_{K,i} - E_{\lambda,i}.
\]
In later notation used in the fitting discussion, this becomes
\[
\Delta E_T = E_K - E_a.
\]

The lattice diffusivity was fixed as
\[
D_L = 7.32\times 10^{-8}\exp\!\left(-\frac{E_L}{RT}\right),
\qquad
E_L = 0.059\ \text{eV},
\]
with
\[
\nu = 10^{13}\ \text{s}^{-1},
\qquad
N_L = 5.2\times 10^{29}\ \text{sites·m}^{-3}.
\]
These constants were used in the finite-element model for ferritic steels [2507.15711].

A related arXiv formulation presents the same physical content in terms of occupancies of lattice and trap sites rather than directly in terms of \(\theta_i\) and \(K_i,\lambda_i\) [2504.17532]. There,
\[
C_L = y_L\,N_L,
\qquad
C_{Tj} = y_{Tj}\,N_{Tj},
\]
and
\[
C_{\text{tot}} = y_L\,N_L + \sum_j y_{Tj}\,N_{Tj}.
\]
The transport equation is written as
\[
\frac{\partial C_L}{\partial t} = D_H\,\frac{\partial^2 C_L}{\partial z^2} - \sum_j \frac{\partial C_{Tj}^{\text{in}}}{\partial t} + \sum_j \frac{\partial C_{Tj}^{\text{out}}}{\partial t}.
\]

In that specialization, adsorption and desorption are occupancy-dependent and thermally activated. The trapping term carries the barrier \(Q_D\), whereas the detrapping term carries the barrier \(Q_D+E_{Tj}\). The same paper states that diffusion and adsorption have the same activation barrier \(Q_D\), while desorption must overcome both diffusion and trap binding. This asymmetry is central to its representation of slow release from deep traps during TDS [2504.17532].

The relation to Oriani’s local-equilibrium approximation is explicit. One paper states that Oriani’s simplification retains only \(\Delta E_T\), so the separate barriers \(E_K\) and \(E_a\) cannot be extracted, and therefore deliberately uses the full McNabb–Foster kinetic model to “capture as much detail as possible regarding the trap energetics” [2507.15711]. Another presents the Oriani occupancy relation
\[
y_T = \frac{y_L \exp\!\left(\frac{E_T}{RT}\right)} {1+y_L \exp\!\left(\frac{E_T}{RT}\right)}
\]
but keeps finite trapping and detrapping kinetics in its numerical model [2504.17532]. This suggests that the McNabb–Foster framework is treated as the more general kinetic description, with Oriani recovered only in appropriate fast-equilibrium limits.

## 3. Boundary-value problems for permeation, desorption, and TDS

The McNabb–Foster equations are used as a unified framework for EP and TDS. In the finite-element implementation applied to ferritic steels, the EP problem used the initial condition
\[
C_L(t=0)=0,
\]
a fixed charging-side lattice concentration at the cathodic surface,
\[
C_{L,x=0} = \frac{i_{\max} X}{E_c D_L},
\]
and
\[
C_{L,x=L}=0
\]
at the anodic side. The observable current density at the exit surface was computed from the flux \(F\) as
\[
i_{x=L}=E_c F.
\]
The same study notes that a simple Fickian baseline,
\[
I = I_{\max}\left\{ 1 + 2 \sum_{n=1}^{\infty} (-1)^n \exp\!\left(-\frac{D_{\text{eff}}n^2\pi^2 t}{X^2}\right)\right\},
\qquad
D_{\text{eff}}=\frac{X^2}{6t_{\text{lag}}},
\]
is inadequate for the measured curves, which deviate significantly from the analytical Fickian form [2507.15711].

For TDS, the same paper imposed zero lattice concentration on both surfaces,
\[
C_L=0 \quad \text{at both surfaces, all times},
\]
took the initial total hydrogen concentration from the integrated TDS spectrum, and assumed that all initial hydrogen is trapped because any lattice hydrogen would have escaped before heating started. The measured temperature history was imposed experimentally, and the total desorption signal was taken as the sum of the fluxes from both surfaces [2507.15711].

A related sensitivity-analysis study uses the same unified sequence of charging/permeation, free desorption, and TDS. During permeation it imposes
\[
C_L(0,t)=C_0,
\qquad
C_L(L,t)=0,
\]
and during free desorption and thermal desorption both surfaces are treated as zero-concentration sinks. The surface flux is computed as
\[
J_{\text{out}}(t) = m_H\,D_H\, \frac{C_{L,i-1}(t)-C_{L,i}(t)}{\Delta z}.
\]
That paper emphasizes that a TDS peak is not a pure trap-emptying event but arises from coupled detrapping, lattice diffusion, specimen thickness, and heating rate [2504.17532].

## 4. Numerical realization and parameter-estimation workflows

One arXiv implementation solves the McNabb–Foster equations in **FESTIM**, using a 1D through-thickness membrane geometry discretized with **100 unidimensional linear elements of equal length**. Typical time steps were **100 s** for EP and **5 s** for TDS. In all fitting examples, **a single type of trap was assumed**. Least-squares fitting was performed using **lmfit**, based on the **Levenberg–Marquardt algorithm**, and the optimized parameters were usually \(N\), \(E_K\), and \(E_a\) [2507.15711].

That same study also documents what was not explicitly specified: no objective-function equation was printed, no weighting strategy was described, no parameter bounds were reported, no initial guesses were reported, and no convergence criteria were reported. Instead, the practical strategy for identifiability assessment was to fix one parameter—especially \(N\)—over broad ranges and re-optimize the others, then examine whether good fits persisted [2507.15711].

A second implementation uses an explicit finite-difference scheme rather than finite elements. The specimen thickness is divided into \(n+1\) uniformly spaced nodes,
\[
i=0,\dots,n,
\qquad
\Delta z = L/n,
\]
and the inner-node update is
\[
C_{L,i}(t+\Delta t) = C_{L,i}(t) - \Delta C_i^{\text{in}} + \Delta C_i^{\text{out}} + \Delta C_{L,i},
\]
with
\[
\Delta C_{L,i} = \Delta t\,D_H \frac{C_{L,i-1}(t)-2C_{L,i}(t)+C_{L,i+1}(t)}{\Delta z^2}.
\]
The time step must satisfy the von Neumann condition
\[
\Delta t \le \frac{\Delta z^2}{2D_H},
\]
and additional local mass-balance constraints are enforced so that trapping and detrapping increments do not exceed the available lattice or trap populations [2504.17532].

Taken together, these implementations show that the McNabb–Foster equations are not tied to a single numerical method. What remains invariant is the coupled diffusion–trapping–detrapping structure and the use of the same governing system across permeation and desorption protocols.

## 5. Identifiability, parameter coupling, and non-uniqueness

The most important recent result concerning the McNabb–Foster equations is not the existence of a forward solver but the non-uniqueness of inverse parameter identification. In ferritic steels containing \((\mathrm{Ti},\mathrm{Cr})\mathrm{C}\) particles, EP and TDS curves were fitted with multiple combinations of trap density and energetic parameters, leading to the conclusion that the system was overdetermined and that it was not possible to determine the individual trapping parameters using this procedure [2507.15711].

The EP example for specimen T3 is explicit. With
\[
N=1.4\times10^{24}\ \text{m}^{-3},
\]
the fitted values were
\[
E_K = 0.4609\ \text{eV},\qquad E_a = 1.017\ \text{eV},\qquad \Delta E_T = 0.5561\ \text{eV}.
\]
With
\[
N=1.9\times10^{24}\ \text{m}^{-3},
\]
the fitted values were
\[
E_K = 0.1244\ \text{eV},\qquad E_a = 0.6351\ \text{eV},\qquad \Delta E_T = 0.5107\ \text{eV}.
\]
Across all EP measurements,
\[
E_K = 0.10\text{–}0.46\ \text{eV},
\qquad
E_a = 0.60\text{–}0.96\ \text{eV},
\qquad
\Delta E_T = 0.45\text{–}0.52\ \text{eV}.
\]

For TDS, good fits were obtained for specimen T4 over the much wider range
\[
N = 5\times10^{24} \text{ to } 5\times10^{27}\ \text{m}^{-3}.
\]
Within this range, \(E_a\) stayed relatively consistent, \(E_K\) increased approximately linearly with \(\ln N\), and \(\Delta E_T\) decreased linearly with \(\ln N\). The gradient of \(\Delta E_T\) versus \(\ln N\) was \(-0.033\) eV for the T4 example and \(-0.03\) to \(-0.045\) eV across materials [2507.15711].

The mathematical explanation is the coupling
\[
\frac{NK}{\lambda}=C,
\]
which, after substituting the Arrhenius forms, becomes
\[
-\Delta E_T = -RT\ln(N) + RT\ln(CN_L).
\]
At \(T=303\) K, \(-RT \approx -0.026\) eV; at TDS peak temperatures \(350\)–\(400\) K, \(-RT \approx -0.030\) to \(-0.034\) eV. These values are consistent with the observed slopes. A plausible implication is that EP and TDS constrain a combined trapping response more strongly than they constrain the individual quantities \(N\), \(E_K\), and \(E_a\) [2507.15711].

The comparison with Kissinger analysis reinforces this point. One study applied
\[
\frac{\partial \ln(\phi/T_M^2)}{\partial (1/T_M)} = -\frac{E_D}{R}
\]
and, after peak deconvolution, reported values roughly in the range \(0.22\) to \(0.37\) eV, with most values around \(0.27\) to \(0.31\) eV; its abstract summarizes this as a trap binding energy of about \(0.24\) eV, albeit with a high degree of uncertainty [2507.15711]. A separate sensitivity analysis also states that TDS peak temperature depends on trap binding energy, trap density, specimen thickness, and heating rate, so peak position alone is not a unique measure of binding energy [2504.17532].

## 6. Scope, applications, and terminological boundaries

Within the recent arXiv literature, the McNabb–Foster equations are applied to hydrogen transport in ferritic steels and are used specifically to interpret EP and TDS in systems containing microstructural traps. In ferritic steels containing \((\mathrm{Ti},\mathrm{Cr})\mathrm{C}\) particles, measurements showed that fine particles **< 5 nm** slowed hydrogen diffusion significantly, whereas coarser particles with average diameter **> 10 nm** had little or no effect. The McNabb–Foster framework reproduced these observations in a forward sense but did not isolate a unique trap parameter set in the inverse sense [2507.15711].

The framework also carries a clear methodological warning. One study explicitly recommends against over-interpreting a single “best-fit” McNabb–Foster parameter set from EP/TDS alone and instead suggests reporting admissible parameter ranges, correlations between parameters, and uncertainty or non-uniqueness [2507.15711]. Another reaches a related conclusion from sensitivity analysis: high-energy trap binding energy can strongly affect TDS while having little effect on permeation time lag once traps are saturated, and TDS peak temperature is influenced by geometry and trap density as well as binding energy [2504.17532].

A common terminological confusion arises from the word “Foster.” The phrase **“McNabb–Foster equations”** does not appear in the paper on the Foster–Hart measure of riskiness, where the relevant equation is
\[
E\log(1+\lambda X)=0,
\]
and it also does not appear in recent papers on non-Foster electromagnetic media, temporal metastructures, or photonic time crystals [1301.1471], [2304.03861], [2509.00795]. In materials science usage, by contrast, the term refers to the kinetic trapping–diffusion framework for hydrogen in metals.

In that restricted and technically specific sense, the McNabb–Foster equations denote a non-equilibrium, reversible, saturable trap model for hydrogen transport. Their principal value lies in unifying diffusion, trapping, detrapping, permeation, and desorption within one constitutive system. Their principal limitation, as recent arXiv studies emphasize, is that good agreement with EP and TDS does not by itself establish unique microscopic trap parameters.

Source: https://www.emergentmind.com/topics/mcnabb-foster-equations