Introduction

lindhard is a clean-room, MIT-licensed Rust engine for Monte Carlo transport of ions in matter, in the binary-collision approximation (BCA). For a beam of ions entering a layered amorphous target it computes where the ions stop, how much damage their cascades make, and what is sputtered and backscattered. A low-energy electron engine is under development; see Status and roadmap.

This book has two parts.

• Physics manual. One page per model: the equations, the assumptions, the range over which the model is meant to be used, and the papers it comes from. Each page names the input keys and the library types that select the model, so you can go from a line in an input file to its physics and back.
• User guide. Installing, a first run from start to finish (5 keV boron into silicon), the reference for the TOML input, and the output formats.

The source of truth for every equation is the code and its doc comments (cargo doc -p lindhard --open); this manual is written from them and checked against them in CI where that can be automated (every model variant must be named in the manual). The design notes, the validation results and the data provenance table stay in the repository’s docs/ directory:

• docs/validation.md: analytic checks, code-to-code comparisons and comparisons with measurement, with the known deviations.
• docs/data-provenance.md: where every number in the tree comes from, and how far it has been verified.
• docs/architecture.md: the crate layout and the milestone plan.

lindhard is early, pre-release software. The project’s rules (the clean-room license tiers in particular) are in CONTRIBUTING.md.

Conventions and the clean-room position

Units

The library works in SI internally. At the input and output boundary the units are in the key names: energies in eV (_ev), lengths in nm (_nm), angles in degrees (_deg), densities in g/cm³. Electronic stopping cross sections are reported per target atom, in the customary eV·10⁻¹⁵ cm² (to_ev_1e15_cm2 converts from J m²).

Throughout the manual \( Z_1, M_1 \) are the atomic number and mass of the moving particle, \( Z_2, M_2 \) those of the target atom, \( E \) the laboratory kinetic energy, \( N \) the atom density, \( a_0 \) the Bohr radius, \( v_0 = \alpha c \) the Bohr velocity, and \( e^2 \) stands for \( e^2 / 4\pi\varepsilon_0 \). Physical constants are CODATA 2022 (Mohr et al. 2025).

Reduced variables

The collision models use the reduced variables of Lindhard, Scharff and Schiøtt (1963). With a screening length \( a \) (see Screening lengths):

\[ x = \frac{r}{a}, \qquad \beta = \frac{b}{a}, \qquad \varepsilon = \frac{a\, E_\mathrm{cm}}{Z_1 Z_2 e^2}, \qquad E_\mathrm{cm} = E \frac{M_2}{M_1 + M_2}, \]

with \( r \) the separation and \( b \) the impact parameter. In these variables the scattering angle depends only on \( (\varepsilon, \beta) \) and on the screening function.

Determinism

Every model is a pure function of its arguments. Random numbers come from a counter-based stream (ChaCha8) keyed on the run seed and the index of the primary history, so a run gives the same bits on one thread or many. A model page that introduces randomness says which draws it makes and in which order.

The clean-room position

lindhard implements published physics from the papers. It never uses the code of closed or copyleft programs (the tiers are set out in CONTRIBUTING.md). Data is different: data is not licensed, code is (operator ruling, 2026-11-10, lindhard#629), and lindhard aims to be a complete open replacement for SRIM/TRIM. In particular:

• Tabulated stopping is admissible data, SRIM, ICRU and ZBL tables included, as long as it carries its provenance. The electronic stopping models are computed from closed-form physics by default; a table, from a user or committed to the tree, is a separately named stopping source, and a comparison against the code it came from is a consistency check, not independent validation.
• The ZBL universal screening function (eight published coefficients) and the ZBL reduced nuclear stopping fit are used. The electronic stopping coefficient sets of the same book are not used yet; they are admissible data.
• Terms whose only published coefficients are tabulated (the Barkas term, shell corrections, multi-oscillator density-effect parameters, Chu and Yang-O’Connor-Wang straggling) were declined under the former no-tables rule. They are now admissible and stay listed in Validity ranges and declined terms until each is implemented.

Every coefficient that enters the code as a number has a row in docs/data-provenance.md saying where it comes from and how far it has been verified. Where a value has only been checked against a secondary source, its page says so.

How a model page is laid out

Each page has the same sections: what the model is, the equations, the assumptions, the validity range, how to select it (the TOML key and the Rust type), its verification status, and its references. New models start from the template book/src/models/_template.md in the repository.

References

• J. Lindhard, M. Scharff, H. E. Schiøtt, Mat. Fys. Medd. Dan. Vid. Selsk. 33 (14) (1963).
• P. J. Mohr, D. B. Newell, B. N. Taylor, E. Tiesinga, Rev. Mod. Phys. 97, 025002 (2025), doi:10.1103/RevModPhys.97.025002 (CODATA 2022).
• D. J. Bernstein, ChaCha, a variant of Salsa20 (2008); J. K. Salmon, M. A. Moraes, R. O. Dror, D. E. Shaw, Proc. SC’11 (2011) (counter-based random streams).

Interatomic potentials

Code: lindhard/src/ion/potential.rs (Screening, Potential).

Model

Two atoms at separation \( r \) interact through a screened Coulomb potential

\[ V(r) = \frac{Z_1 Z_2 e^2}{r}\, \phi\!\left(\frac{r}{a}\right), \]

where \( \phi \) is a universal screening function (\( \phi(0) = 1 \), decreasing to 0) and \( a \) a screening length that depends on \( Z_1, Z_2 \) (Screening lengths). The scattering code works in \( x = r/a \), where \( \phi \) is the only input; the screening length only converts to and from SI.

Three of the four screening functions are sums of exponentials,

\[ \phi(x) = \sum_i c_i\, e^{-b_i x}, \]

and the fourth is a polynomial times an exponential.

The screening functions

The default length is the one each function was introduced or fitted with (Screening::default_length); screening_length overrides it.

ZBL universal

Four exponentials fitted by Ziegler, Biersack and Littmark (1985, ch. 2) to Hartree-Fock-Slater solid-state pair potentials:

These are the commonly printed four-digit rounding of the published values.

Kr-C

Three exponentials fitted to the Hartree-Fock Kr-Kr pair interaction (Wilson, Haggmark and Biersack 1977):

Molière

Molière’s three-exponential approximation to the Thomas-Fermi screening function (Molière 1947):

Lenz-Jensen

A polynomial-times-exponential approximation to the Thomas-Fermi-Jensen statistical model (Lenz 1932; Jensen 1932), in the form printed by Möller (2017, p. 11, eq. (28)):

\[ \phi(x) = \left(1 + q + 0.3344\, q^2 + 0.0485\, q^3 + 0.002647\, q^4\right) e^{-q}, \qquad q = \sqrt{9.67\, x}. \]

Assumptions

• Pair potentials: the interaction of two atoms does not depend on any third atom.
• The screening function is universal: \( Z_1 \) and \( Z_2 \) enter only through the screening length.
• The potentials are purely repulsive. There is no attractive well, so they are not meant for energies comparable with chemical binding (a few eV), where the cutoffs and binding energies of the BCA take over.

Validity

Screened Coulomb potentials describe the repulsive part of the interaction, from the close collisions of keV ions down to the tens of eV of cascade atoms. The ZBL universal function was fitted to pair potentials over a broad set of ion-target pairs and is the usual default. It is not exact for any particular pair: the computed ranges of B in amorphous Si run long against the measurement, and the validation attributes the offset below about 5 keV mainly to the nuclear stopping (the ZBL function being too soft for B on Si) (docs/validation.md, level 3).

Verification status

From docs/data-provenance.md: the ZBL coefficients have been verified against the SRIM book (Ziegler, Biersack and Ziegler 2015, ch. 2, Figs. 2-17 and 2-18; its Eq. 2-74 has two misprints, recorded there), and so have the Molière coefficients (Eq. 2-62; its middle prefactor is misprinted as 5.5). The Kr-C and Molière coefficients have also been cross-checked against an independent MIT-licensed implementation (ir2-lab/screened_coulomb, commit f84c3c8), a secondary source. The Lenz-Jensen polynomial agrees to within 9 % with the book’s own three-exponential fit (Eq. 2-63), which supports its scale; the set agrees with a textbook tabulation (Möller 2017) but has not been checked against the primary papers, which were not accessible.

References

• J. F. Ziegler, J. P. Biersack, U. Littmark, The Stopping and Range of Ions in Solids (Pergamon, New York, 1985), ch. 2.
• J. F. Ziegler, J. P. Biersack, M. D. Ziegler, SRIM: The Stopping and Range of Ions in Matter (SRIM Co., 2015), ch. 2.
• W. D. Wilson, L. G. Haggmark, J. P. Biersack, Phys. Rev. B 15, 2458 (1977).
• G. Molière, Z. Naturforsch. A 2, 133 (1947).
• W. Lenz, Z. Phys. 77, 713 (1932), doi:10.1007/BF01342150.
• H. Jensen, Z. Phys. 77, 722 (1932), doi:10.1007/BF01342151.
• W. Möller, Fundamentals of Ion-Solid Interaction, HZDR-073 (Helmholtz-Zentrum Dresden-Rossendorf, 2017), p. 11, eqs. (26)-(30), https://www.hzdr.de/publications/PublDoc-10091.pdf.

Screening lengths

Code: lindhard/src/ion/potential.rs (ScreeningLength).

Model

The screening length \( a(Z_1, Z_2) \) scales the separation in the screening function, \( x = r/a \). Three forms are implemented. Two of them use the Thomas-Fermi constant

\[ C_\mathrm{TF} = \left(\frac{9\pi^2}{128}\right)^{1/3} = 0.8853, \]

in the closed form of Firsov (1958, p. 535), printed as 0.8853 by Lindhard, Scharff and Schiøtt (1963, p. 8).

Without screening_length the length paired with the potential is used: universal for ZBL, Firsov for Kr-C and Molière, Lindhard for Lenz-Jensen. The run echoes the length it used in summary.json (input.physics.screening_length).

Assumptions

• The universal length is the one the ZBL screening function was fitted with; the other two come from Thomas-Fermi scaling arguments for the two-atom system.
• Mixing a screening function with a length other than its own is allowed (it is a model choice), but the fitted functions were fitted with their own length.

Validity

As for the potentials: screened Coulomb collisions from cascade energies upward.

Verification status

From docs/data-provenance.md: the Firsov and Lindhard lengths have been verified against scans of the primary papers (page references above). The universal length uses the same constant: the SRIM book (Ziegler, Biersack and Ziegler, 2015, Ch. 2) confirms the exponent 0.23 and prints the prefactor as 0.8853 in Eq. 2-73 (p. 2-30), while the labels of its Figs. 2-16 and 2-17 (p. 2-32) print 0.8854, a 10⁻⁴ relative difference. The code follows the equation and uses the Thomas-Fermi constant, 0.88534.

References

• O. B. Firsov, Sov. Phys. JETP 6, 534 (1958), pp. 535-536.
• J. Lindhard, M. Scharff, H. E. Schiøtt, Mat. Fys. Medd. Dan. Vid. Selsk. 33 (14) (1963), p. 8.
• J. F. Ziegler, J. P. Biersack, U. Littmark, The Stopping and Range of Ions in Solids (Pergamon, New York, 1985), ch. 2.
• J. F. Ziegler, J. P. Biersack, M. D. Ziegler, SRIM: The Stopping and Range of Ions in Matter (SRIM Co., 2015), ch. 2, pp. 2-29 to 2-32.

The scattering integral

Code: lindhard/src/ion/scattering.rs (theta_quadrature, theta_magic, ScatteringTable).

Model

A binary collision in a central potential is classical elastic scattering. In the reduced variables \( x = r/a \), \( \beta = b/a \), \( \varepsilon = a E_\mathrm{cm} / (Z_1 Z_2 e^2) \), the radial motion in the centre-of-mass frame obeys

\[ G(x) = 1 - \frac{\phi(x)}{\varepsilon x} - \frac{\beta^2}{x^2}, \]

with the distance of closest approach \( x_0 \) the root of \( G(x_0) = 0 \). The centre-of-mass scattering angle is (Goldstein 1980; Ziegler, Biersack and Littmark 1985, ch. 2)

\[ \theta = \pi - 2\beta \int_{x_0}^{\infty} \frac{dx}{x^2 \sqrt{G(x)}}. \]

Distance of closest approach

\( G \) is strictly increasing, so the root is unique and lies in \( [\beta,\, (1/\varepsilon + \sqrt{1/\varepsilon^2 + 4\beta^2})/2] \). Newton’s method runs inside that bracket, and a step that leaves it is replaced by bisection.

Gauss-Mehler quadrature

With \( u = x_0/x \) and then \( u = \cos t \) the integrable endpoint singularity is removed, leaving

\[ \int_0^{\pi/2} \frac{dt}{\sqrt{H(\cos t)}}, \qquad H(u) = \frac{\beta^2}{x_0^2} + \frac{\phi(x_0) - u\, \phi(x_0/u)}{\varepsilon x_0 (1 - u^2)}, \]

which the midpoint rule in \( t \) (the Gauss-Mehler, or Gauss-Chebyshev of the first kind, nodes; Mendenhall and Weller 1991) integrates with exponential convergence. The default is 64 nodes. \( H \) is evaluated in this cancellation-free form.

Precomputed angle table (the transport hot path)

The BCA does not integrate at every collision. Once per run and per \( (Z_1, Z_2, \text{potential}) \) pair it builds a table of \( y = \ln\tan(\theta/2) \) on a grid uniform in \( \ln\varepsilon \) and \( \ln\beta \), by quadrature, and interpolates bilinearly in those coordinates. \( y \) is close to linear in \( \ln\beta \) at both ends, which is what makes bilinear interpolation accurate, and an error \( \delta y \) in the table is at most \( \delta y \) radians in \( \theta \). At build time the interpolation error is measured against direct quadrature at cell centres on a deterministic subset of cells; the grid and the largest absolute and relative errors found are written to summary.json under physics.scattering_table. These are empirical bounds from sampling, not proofs. The run grid is \( 10^{-6} \le \varepsilon \le 10^{4} \), \( 10^{-5} \le \beta \le 10^{2} \), 32 points per decade (lindhard::input::TABLE_SPEC). Beyond the largest tabulated \( \beta \) the angle is taken as 0 (the constant-free-path \( p_\mathrm{max} \) is about 20 screening lengths or less at solid densities); below the smallest it is clamped; and at a reduced energy outside the table the engine falls back to direct quadrature (deterministic, only slower).

Magic formula (cross-check only)

The Biersack-Haggmark “magic formula” is kept as a fast approximate cross-check, not used in transport:

\[ \cos\frac{\theta}{2} = \frac{\beta + \rho + \Delta}{x_0 + \rho}, \qquad \Delta = \frac{A (x_0 - \beta)}{1 + G}, \]

\[ A = 2\alpha\varepsilon\beta^{b}, \quad \alpha = 1 + C_1 \varepsilon^{-1/2}, \quad b = \frac{C_2 + \varepsilon^{1/2}}{C_3 + \varepsilon^{1/2}}, \quad G = \frac{\gamma}{\sqrt{1 + A^2} - A}, \quad \gamma = \frac{C_4 + \varepsilon}{C_5 + \varepsilon}, \]

with \( \rho \) the radius of curvature of the trajectory at closest approach. Constant sets exist for the ZBL universal function (Ziegler, Biersack and Littmark 1985) and for Molière (Biersack and Haggmark 1980).

Nuclear stopping

The reduced nuclear stopping cross section follows from the angle,

\[ s_n(\varepsilon) = 2\varepsilon \int_0^\infty \sin^2\frac{\theta}{2}\, \beta\, d\beta, \]

integrated over \( \ln\beta \) by Simpson’s rule. It is used in tests and validation, not in transport (the BCA samples the collisions themselves). Two references are checked against it: the ZBL universal fit \( s_n = \ln(1 + 1.1383\varepsilon) / [2(\varepsilon + 0.01321\varepsilon^{0.21226} + 0.19593\varepsilon^{0.5})] \) for \( \varepsilon \le 30 \) (Ziegler, Biersack and Littmark 1985), and the high-energy limit for a sum of exponentials from the first-order impulse approximation of Lindhard, Nielsen and Scharff (1968, eq. (3.4)).

Assumptions

• Classical, elastic, binary collisions in a central potential; no inelastic energy loss inside the collision itself (electronic loss is handled separately, see Electronic stopping).
• The angle depends only on \( (\varepsilon, \beta) \) and the screening function; \( Z_1, Z_2 \) and the screening length enter through the reduced variables.

Validity

The classical treatment holds when the collision is well localised compared with the screening length, which is the case for heavy particles from cascade energies (a few eV) to well above the MeV range.

Verification status

The quadrature is checked against the small-angle perturbation formula of Lindhard, Nielsen and Scharff (1968) and against the ZBL nuclear stopping fit; the magic formula is checked against the quadrature for both constant sets. The ZBL nuclear stopping fit, the magic formula and its ZBL constants have been verified against the SRIM book (Ziegler, Biersack and Ziegler 2015, ch. 2, Eqs. 2-89 and 2-90, Eqs. 2-102 to 2-105 and Fig. 2-21; the figure lists C4 and C5 in the other order, which the quadrature check rejects). The ZBL constants were also cross-checked against ir2-lab/screened_coulomb (MIT). The Molière set has been verified against the book’s Table 4-1 (p. 7-6); the 1980 paper itself was not accessible. That table also lists the universal C4 and C5 in the other order, which the quadrature check rejects, as do the C4/C5 ratios of its own Molière and Kr-C columns. See docs/data-provenance.md.

References

• H. Goldstein, Classical Mechanics, 2nd ed. (Addison-Wesley, 1980).
• J. F. Ziegler, J. P. Biersack, U. Littmark, The Stopping and Range of Ions in Solids (Pergamon, New York, 1985), ch. 2.
• J. F. Ziegler, J. P. Biersack, M. D. Ziegler, SRIM: The Stopping and Range of Ions in Matter (SRIM Co., 2015), ch. 2.
• M. H. Mendenhall and R. A. Weller, Nucl. Instrum. Methods B 58, 11 (1991).
• J. P. Biersack and L. G. Haggmark, Nucl. Instrum. Methods 174, 257 (1980), doi:10.1016/0029-554X(80)90440-1.
• J. Lindhard, V. Nielsen, M. Scharff, Mat. Fys. Medd. Dan. Vid. Selsk. 36 (10) (1968), eqs. (3.3)-(3.4).

Electronic stopping

Code: lindhard/src/ion/stopping/.

A moving atom loses energy to the target electrons as well as to the nuclei. In the BCA the electronic loss is handled separately from the collisions: continuously along each free flight (nonlocal), at the collisions as a function of the distance of closest approach (local), or as a mix of the two (BCA transport, “Electronic loss”).

Conventions

• Every model returns the stopping cross section per target atom, \( S_e = -\frac{1}{N}\frac{dE}{dx} \), in J m² internally and eV·10⁻¹⁵ cm² at the boundary.
• Energies are the laboratory kinetic energy of the moving atom.
• Every model is a pure function of (ion, target \( Z_2 \), energy): no random numbers, no state.
• Each model exposes an advisory validity range (ElectronicStopping::validity). The CLI warns on stderr when the beam energy lies outside it; it does not refuse the run.
• Compounds and mixtures combine the element cross sections by Bragg additivity.

Choosing a model

A [stopping] table replaces the chosen model for the one (ion, target element) pair it declares (User stopping tables). Energy-loss straggling is described on its own page. It is a library function that the BCA engine does not call at present: the electronic loss along each flight is deterministic, \( N S_e(E)\, s \).

No electronic stopping

stopping = "none" (StoppingChoice::None, model lindhard::ion::stopping::none::NoStopping) sets \( S_e = 0 \) for every ion, target and energy, so the moving atoms lose energy only in nuclear collisions. This is the \( k = 0 \) limit of the range theory of Lindhard, Scharff and Schiott (Mat. Fys. Medd. Dan. Vid. Selsk. 33 (14), 1963), where the electronic stopping coefficient vanishes.

It is not a physical model of any target. It is there for two uses: like-for-like transport comparisons against codes run with their electronic stopping switched off (the *_nuclear level-2 problems of the validation report), and nuclear-only studies such as the level-1 range.si5k_si.* and damage checks, which use the same type. It has no validity range, so no warning is raised for it, and the electronic terms of the energy budget are exactly zero. [stopping] tables may still be given; they serve the pairs they declare.

The validity ranges side by side, and the terms that are deliberately not implemented, are on Validity ranges and declined terms.

Lindhard-Scharff

Code: lindhard/src/ion/stopping/lindhard_scharff.rs (LindhardScharff).

Model

At ion velocities below \( v_0 Z_1^{2/3} \) the electronic stopping is proportional to the velocity (Lindhard and Scharff 1961). In the reduced units of Lindhard, Scharff and Schiøtt (1963):

\[ a = 0.8853\, a_0 \left(Z_1^{2/3} + Z_2^{2/3}\right)^{-1/2}, \qquad \varepsilon = \frac{E\, a\, M_2}{Z_1 Z_2 e^2 (M_1 + M_2)}, \qquad \rho = N x\, 4\pi a^2 \frac{M_1 M_2}{(M_1 + M_2)^2}, \]

the reduced electronic stopping is

\[ \left(\frac{d\varepsilon}{d\rho}\right)_e = k_L\, \varepsilon^{1/2}, \qquad k_L = \xi_e\, \frac{0.0793\, Z_1^{1/2} Z_2^{1/2} (A_1 + A_2)^{3/2}} {\left(Z_1^{2/3} + Z_2^{2/3}\right)^{3/4} A_1^{3/2} A_2^{1/2}}, \qquad \xi_e = Z_1^{1/6}, \]

with \( A_1, A_2 \) the masses in u. The equivalent dimensional form is

\[ S_e = \xi_e\, \frac{8\pi e^2 a_0\, Z_1 Z_2}{\left(Z_1^{2/3} + Z_2^{2/3}\right)^{3/2}}\, \frac{v}{v_0}. \]

The factor \( \xi_e \approx Z_1^{1/6} \) (order 1 to 2) is part of the Lindhard-Scharff result. A test checks the reduced form, with the rounded 0.0793, against the dimensional form to 1 % for all tested ion-target pairs.

A per-element multiplicative correction \( f(Z_2) \) can be attached in the library (LindhardScharff::with_correction, default 1; it must be finite and non-negative, and 0 switches the element’s stopping off). The CLI does not expose it; use a user table for measured stopping instead.

Selecting it

stopping = "lindhard-scharff" in [physics] (the default), StoppingChoice::LindhardScharff. All of the loss is applied nonlocally, along the free flights. It is also the nonlocal half of stopping = "equipartition-ls-or" (Oen-Robinson).

Assumptions

• A free-electron-gas picture with Thomas-Fermi scaling; the stopping is smooth in \( Z_1 \) and \( Z_2 \) and has no shell structure (the measured \( Z_1 \) oscillations are not reproduced).
• No dependence on the chemical or physical state of the target.

Validity

\( v < v_0 Z_1^{2/3} \), that is \( E/A_1 \) below about \( 25\ \mathrm{keV} \cdot Z_1^{4/3} \). This is the regime of keV-range implantation and of cascade atoms. Above it the stopping peaks and turns over, which this model does not describe; see Bethe-Bloch for high velocities and Validity ranges and declined terms for why no interpolation joins the two.

The validation pages report a known offset: for B in Si the computed ranges run long against the measurement, and at 10 to 20 keV part of the offset is attributed to LS stopping being too small for that pair (docs/validation.md, level 3).

Verification status

The constants 0.8853, 0.0793 and \( \xi_e = Z_1^{1/6} \) are checked against the dimensional form in a unit test (1 %); they have not yet been checked digit by digit against the papers (docs/data-provenance.md).

References

• J. Lindhard and M. Scharff, Phys. Rev. 124, 128 (1961).
• J. Lindhard, M. Scharff, H. E. Schiøtt, Mat. Fys. Medd. Dan. Vid. Selsk. 33 (14) (1963).

Oen-Robinson and the equipartition mix

Code: lindhard/src/ion/stopping/oen_robinson.rs (OenRobinson), lindhard/src/ion/stopping/mix.rs (EquipartitionMix).

Model

Oen and Robinson (1976) make the electronic loss local: it is taken at each collision and depends on how close the two atoms came. A collision with distance of closest approach \( r_\mathrm{min} \) loses

\[ \Delta E_e(E, r_\mathrm{min}) = S_\mathrm{LS}(E)\, \frac{c^2}{2\pi a^2}\, \exp\!\left(-\frac{c\, r_\mathrm{min}}{a}\right), \qquad c = 0.3, \]

with \( S_\mathrm{LS} \) the Lindhard-Scharff cross section and \( a \) the Firsov screening length \( 0.8853\, a_0 (Z_1^{1/2} + Z_2^{1/2})^{-2/3} \). The prefactor normalises \( \int \Delta E_e\, 2\pi p\, dp = S_\mathrm{LS} \) when \( r_\mathrm{min} \) is identified with the impact parameter \( p \); that identification is the approximation of the original paper, so the impact-averaged stopping equals the LS value.

Equipartition mix

The equipartition mix splits the Lindhard-Scharff loss in two equal halves: half nonlocal, continuous along the free flight, and half local, at each collision, by the Oen-Robinson formula evaluated at the distance of closest approach. This combination is the one described in the BCA literature that builds on Oen and Robinson (1976). Averaged over impact parameter the mix equals the LS stopping (a unit test checks this to 10⁻¹⁴).

Selecting it

stopping = "equipartition-ls-or" in [physics], which is StoppingChoice::EquipartitionLsOr in the input and ElectronicLoss::EquipartitionLsOr in the engine (see BCA transport). The stopping model passed to the engine is then not used: this mode carries its own LS and Oen-Robinson losses, so it cannot be combined with [stopping] tables. The local losses are reported in the energy budget as electronic_local.

Assumptions

• \( r_\mathrm{min} \approx p \) in the normalisation (Oen and Robinson 1976).
• The local part is only sampled at the collisions the BCA makes, out to \( p_\mathrm{max} \) (to \( p_\mathrm{max}\sqrt{K + 1} \) with \( K \) weak collisions). Where that radius is not large compared with \( a / 0.3 \), the decay length of the local loss, the total electronic stopping comes out somewhat below the LS value.

Validity

The same as Lindhard-Scharff: velocities below \( v_0 Z_1^{2/3} \).

Verification status

The constants \( c = 0.3 \) and the choice of the Firsov length were entered from the contributor’s recollection of the paper and are not verified against it; spot-check them before relying on the local loss (docs/data-provenance.md). The averaging identity with LS is tested.

References

• O. S. Oen and M. T. Robinson, Nucl. Instrum. Methods 132, 647 (1976).
• O. B. Firsov, Sov. Phys. JETP 6, 534 (1958) (screening length).
• J. Lindhard and M. Scharff, Phys. Rev. 124, 128 (1961).

Bethe-Bloch and the density effect

Code: lindhard/src/ion/stopping/bethe.rs (BetheBloch, EffectiveCharge, density_effect_single_oscillator).

Model

At high velocity the electronic stopping cross section per target atom is (Bethe 1930, 1932; Fano 1963; the same equation as the Particle Data Group review of the passage of particles through matter):

\[ S = \frac{4\pi (e^2)^2 z^2 Z_2}{m c^2 \beta^2} \left[ \frac{1}{2}\ln\frac{2 m c^2 \beta^2\gamma^2 W_\mathrm{max}}{I^2}

• \beta^2 - \frac{\delta}{2} - \frac{C}{Z_2} + L_1 + L_2 \right], \]

\[ W_\mathrm{max} = \frac{2 m c^2 \beta^2\gamma^2}{1 + 2\gamma m/M + (m/M)^2}, \]

with \( m \) the electron mass, \( M \) the ion mass, \( z \) the projectile effective charge, \( I \) the mean excitation energy, \( \delta \) the density-effect correction and \( C \) the shell correction. \( \beta^2 \) is computed as \( \tau(\tau + 2)/(1 + \tau)^2 \) with \( \tau = E/Mc^2 \), which stays accurate at low energy.

• Bloch term (Bloch 1933), on by default: \( L_2 = -y^2 \sum_{n \ge 1} \frac{1}{n(n^2 + y^2)} \), \( y = z\alpha/\beta \).
• Mean excitation energy: by default the Bloch rule \( I = 10\ \mathrm{eV} \cdot Z_2 \) (Bloch 1933), a rough estimate (about 20 % in \( I \), about 2 % in \( S \)). The library takes a measured, cited value through BetheBloch::with_mean_excitation_ev.
• Shell correction \( C/Z_2 \) and density effect \( \delta \): 0 by default; the library accepts constants from a cited source (shell_over_z, density_delta).
• Barkas term \( L_1 \): not implemented (see Validity ranges and declined terms).

Effective charge

The CLI uses the bare charge; the Barkas form is available in the library and is not applied by default.

Single-oscillator density effect (opt-in)

The dielectric formulation of the density effect (Fermi 1940; Sternheimer 1952), in the form given by Fano (1963), is

\[ \delta = \sum_i f_i \ln\!\left(1 + \frac{L^2}{\omega_i^2}\right) - \frac{L^2 (1 - \beta^2)}{\beta^2\omega_p^2}, \qquad \sum_i \frac{f_i\, \omega_p^2}{\omega_i^2 + L^2} = \frac{1}{\beta^2} - 1. \]

Specialised to one oscillator (\( f = 1 \), \( \omega_0 = I/\hbar \)) the constraint solves in closed form, \( L^2 = \omega_p^2(\beta\gamma)^2 - \omega_0^2 \), and

\[ \delta = \ln\!\left[(\beta\gamma)^2 \left(\frac{\hbar\omega_p}{I}\right)^2\right] - 1

• \frac{(I/\hbar\omega_p)^2}{(\beta\gamma)^2} \quad \text{for } \beta\gamma > \frac{I}{\hbar\omega_p}, \]

and 0 below, with \( \hbar\omega_p = \hbar\sqrt{n e^2/\varepsilon_0 m} \) the free-electron plasma energy. The one-oscillator reduction is ours, not quoted from a paper. It is enabled with BetheBloch::with_density_effect (library only).

Selecting it

stopping = "bethe-bloch" in [physics], StoppingChoice::BetheBloch: bare charge, Bloch term on, Bloch-rule \( I \), no shell or density correction, all loss nonlocal.

Use it only for problems that stay above its range. Where the bracket is not positive the model returns NotApplicable, and the transport stops with that error rather than silently switching model. For any ion that slows down in the target (a semi-infinite substrate always) the run fails once the energy falls below the bracket’s zero: for example, 300 keV H into Si with stopping = "bethe-bloch" stops with “bethe-bloch is not applicable at 85306 eV”. No low-to-high energy interpolation is implemented, for the reasons given in Validity ranges and declined terms.

Assumptions

• First Born approximation plus the Bloch correction; the projectile charge is fixed by the effective-charge model.
• Without shell and Barkas corrections, expect errors of a few percent at 1 to 10 MeV/u.

Validity

\( v \ge 3 v_0 Z_1^{2/3} \) up to 1 GeV/u (the upper end because the density effect is off by default). The single-oscillator density effect switches on only at \( \beta\gamma > I/\hbar\omega_p \), about 5 to 6 for Si (a proton of about 5 GeV), so it is zero over the whole advisory range. Real materials switch the effect on much earlier (around \( \beta\gamma \approx 1.5 \) for Si), because their oscillator spectrum is spread: the one-oscillator form underestimates \( \delta \) through the transition and is only the correct asymptote for ultra-relativistic ions.

Verification status

Hand-calculated proton-in-Si values (1, 10 and 100 MeV, \( I = 140 \) eV, CODATA constants) are unit tests; they are our own arithmetic, not taken from any stopping table. The density-effect form is tested through its limits (\( \delta \to 0 \) continuously at threshold, and \( \delta \to 2\ln(\beta\gamma\, \hbar\omega_p / I) - 1 \) at large \( \beta\gamma \)). The Bloch rule, the Bloch series and the Barkas effective charge have not been checked against the original papers (docs/data-provenance.md).

References

• H. Bethe, Ann. Phys. 397, 325 (1930); H. Bethe, Z. Phys. 76, 293 (1932).
• F. Bloch, Ann. Phys. 408, 285 (1933).
• U. Fano, Ann. Rev. Nucl. Sci. 13, 1 (1963).
• Particle Data Group, “Passage of particles through matter”, in the Review of Particle Physics.
• W. H. Barkas, Nuclear Research Emulsions I (Academic Press, 1963).
• L. C. Northcliffe, Ann. Rev. Nucl. Sci. 13, 67 (1963).
• E. Fermi, Phys. Rev. 57, 485 (1940).
• R. M. Sternheimer, Phys. Rev. 88, 851 (1952).

User stopping tables

Code: lindhard/src/ion/stopping/table.rs (StoppingTable, TableOverride).

Model

A user table gives \( S_e(E) \) for one projectile (atomic number and mass) in one target element, as data the user supplies. Between the tabulated points the cross section is interpolated piecewise linearly in \( \ln S \) against \( \ln E \), which keeps the monotonicity of each segment of the data. Nothing is extrapolated or clamped: a query outside the table’s energy range is an error (OutOfTableRange) that stops the run.

provenance = "Author, Journal vol, page (year), Table N"   # required
ion_z = 5
ion_mass_amu = 11.0093            # optional; default: standard atomic weight
target_z = 14
energy_ev = [1.0e2, 1.0e3, 1.0e4]                 # strictly increasing, eV
stopping_ev_1e15_cm2 = [10.0, 30.0, 60.0]         # eV 1e-15 cm^2 per atom

(The numbers are a format illustration, not data.)

Provenance is mandatory

A table without a non-empty provenance string is rejected at load time: data without an origin is not admitted. The run records each table in summary.json (physics.stopping_tables: the path as written, the resolved path, the SHA-256 of the file’s bytes, the provenance string, the pair and the energy range), and lists it under physics.models as user-table with its provenance as the source.

Tables from any source may be loaded and committed, SRIM and ICRU tables included (data is not licensed; code is: operator ruling, 2026-11-10, lindhard#629). provenance must name the source precisely: book edition, table and page, or program, version and file.

Selecting it

Declare the files under [stopping] tables (see TOML input reference). A table replaces the [physics] stopping model for exactly the (ion_z, target_z) pair it declares, including recoils of that species when recoils are followed; every other pair uses the model. A table cannot be combined with stopping = "equipartition-ls-or", which carries its own Lindhard-Scharff and Oen-Robinson loss. With stopping = "none" the declared pairs use their tables and every other pair has no electronic stopping.

Assumptions

• The table is tied to the projectile mass it was declared for. A query for a mass that differs by more than a small relative tolerance (tight enough to separate neighbouring isotopes) is an error (TableMassMismatch): the same energy at a different mass is a different speed, and no energy-axis conversion is attempted.
• Compounds are built from element tables by Bragg additivity; a compound table is not a supported input.

Validity

Exactly the table’s energy range. The CLI checks up front that the range contains the beam energy, and for a recoil species that it starts at or below recoil_cutoff_ev; it warns if a table starts above primary_cutoff_ev, because the run fails if a projectile slows below it.

Verification status

The loader, the interpolation and the range and mass checks are unit-tested. The tests use tables generated at test time from our own Lindhard-Scharff model; no table is committed (docs/data-provenance.md).

References

The table’s own provenance string is its reference. The format and the rules are this project’s.

Bragg additivity

Code: lindhard/src/ion/stopping/bragg.rs (bragg_cross_section_per_atom, CompoundCorrection).

Model

The stopping cross section of a compound or mixture, per average atom, is the atom-fraction-weighted sum of the element cross sections (Bragg and Kleeman 1905):

\[ S_\mathrm{compound} = f \sum_j x_j S_j, \qquad -\frac{dE}{dx} = N f \sum_j x_j S_j, \]

with \( x_j \) the atom fractions, \( S_j \) the element cross sections from the chosen model (or a user table for that pair), \( N \) the total atom density and \( f \) a per-compound correction factor.

Selecting it

Always on: every layer’s electronic loss goes through the Bragg sum, which for a pure element is just that element’s cross section. summary.json lists it under physics.models as bragg-additivity. The CLI applies no correction (\( f = 1 \)). The library hook CompoundCorrection (a constant or energy-dependent factor, finite and non-negative) is there for chemical or phase corrections, such as cores-and-bonds schemes, supplied by the caller from a cited source.

Assumptions

• Each atom stops the ion independently of its chemical environment.

Validity

Wherever the element models apply. Bragg additivity ignores the chemical and physical state of the target, and no correction for it is applied unless one is supplied.

Verification status

The sum and the rejection of invalid correction factors are unit-tested.

References

• W. H. Bragg and R. Kleeman, Phil. Mag. 10, 318 (1905).

Energy-loss straggling

Code: lindhard/src/ion/stopping/straggling.rs (StragglingModel, variance_per_atom).

Model

The energy loss of an ion along a path fluctuates. For a thick target in the free-electron picture, the variance of the energy loss per unit areal density is (Bohr 1948)

\[ \Omega_B^2 = 4\pi z^2 Z_2 (e^2)^2 \quad \text{(J}^2\,\text{m}^2 \text{ per target atom)}, \]

independent of energy in the non-relativistic regime, with \( z \) the projectile charge from an effective-charge model. The relativistic variant multiplies it by \( (1 - \beta^2/2)/(1 - \beta^2) \) (Bethe and Livingston 1937; Fano 1963). For a compound the variances add by atom fraction, \( \sum_j x_j \Omega_j^2 \).

Selecting it

Library only (variance_per_atom, variance_per_atom_material). The BCA engine does not use it at present and there is no TOML key: the electronic loss along each free flight is deterministic.

Assumptions

• Free, stationary target electrons; every electron of the target contributes, whatever its binding.

Validity

High energy, \( v \gg v_0 Z_1^{2/3} \). Below about 1 MeV/u Bohr’s value is an overestimate, because the bound electrons contribute less than free ones. The corrections that reduce it at intermediate energy (Chu; Yang-O’Connor-Wang) are declined, and the Lindhard-Scharff low-velocity correction is omitted (not verified); see Validity ranges and declined terms.

Verification status

Not verified against the original papers (docs/data-provenance.md).

References

• N. Bohr, Mat. Fys. Medd. Dan. Vid. Selsk. 18 (8) (1948).
• H. A. Bethe and M. S. Livingston, Rev. Mod. Phys. 9, 245 (1937).
• U. Fano, Ann. Rev. Nucl. Sci. 13, 1 (1963).
• W. K. Chu, Phys. Rev. A 13, 2057 (1976) (declined correction).
• Q. Yang, D. J. O’Connor, Z. Wang, Nucl. Instrum. Methods B 61, 149 (1991) (declined correction).

Validity ranges and declined terms

This page is docs/stopping-models.md, included verbatim so there is one copy. It collects the validity ranges of the electronic stopping models side by side, and records the terms that are deliberately not implemented and why (the clean-room rules exclude coefficient sets that are only published as tables). The model pages are Lindhard-Scharff, Oen-Robinson, Bethe-Bloch, user tables, Bragg additivity and straggling; the references cited below are listed in full on those pages and at the end of this one.

Code: lindhard/src/ion/stopping/. All models return the stopping cross section per atom in J m² (to_ev_1e15_cm2 converts to eV·10⁻¹⁵ cm²) and expose the ranges below at run time through ElectronicStopping::validity. Ranges are advisory; the models do not refuse energies outside them (Bethe-Bloch returns NotApplicable where its bracket is not positive, and user tables return OutOfTableRange). v0 = α c is the Bohr velocity, Z1 the ion charge number.

Per-pair corrections

correction::fit_correction fits a smooth multiplicative factor g(E) of a base model to the measured points of one ion on one element target (StoppingDataset), and Correction applies it. The fit is deterministic and seeded; the held-out split, the bootstrap band (CorrectionFit::band) and the held-out residuals of the base and the corrected model side by side (HeldOutReport) come with it, and a correction that does not lower the held-out rms says so (reduces_rms_relative). The correction is applied only inside the span of its training data. The band covers the scatter of the points about the polynomial form, not measurement systematics or bias of the form, and with no stated point uncertainty it says nothing about measurement error. Example: 15 values of the stopping of Cu for 50 to 600 keV protons (Bader et al. 1956, smoothed curve values) against Lindhard-Scharff, which is far outside its range there (the base is 2 to 4 times too high); the fit is a statement about that data, not a validated model of Cu. Bethe-Bloch returns NotApplicable at 50 keV, so it cannot be fitted over this span (the fit refuses rather than dropping points).

What is deliberately not implemented

The clean-room rules (CONTRIBUTING.md) rule out tabulated coefficient sets from SRIM/ZBL, ICRU reports and similar. As a result:

• Barkas (L1) term: declined (see “Literature findings” below).
• Shell correction C/Z2: declined. Callers may set a constant from a cited source (BetheBloch::shell_over_z); default 0.
• Density effect δ: only the single-oscillator form is implemented (opt-in, BetheBloch::with_density_effect); the tabulated multi-oscillator parameter sets are declined. density_delta can still be set as a constant.
• Mean excitation energy I: defaults to the Bloch rule 10 eV · Z2 (rough); supply a cited measured value with with_mean_excitation_ev.
• Biersack-Varelas interpolation joining the low- and high-energy regimes: documented but not implemented (issue #44, closed; see “Biersack-Varelas interpolation” below). Fitted ZBL/Biersack-Varelas coefficients are not used.
• Straggling: Chu and Yang-O’Connor-Wang corrections declined (see below) and the Lindhard-Scharff low-velocity correction (not verified) is omitted; Bohr plus an optional relativistic factor only.
• Heavy-ion effective charge: the Barkas empirical form is available in bethe and is not applied by default.

Biersack-Varelas interpolation (issue #44, closed)

Status: documented, not implemented.

The joining form

The interpolation joins a low-energy and a high-energy stopping branch harmonically:

1/S = 1/S_low + 1/S_high     (equivalently S = S_low S_high / (S_low + S_high))

Two open-access sources were read (rendered pages) and show this form:

• M. V. Moro, PhD thesis, Universidade de Sao Paulo (2017), doi:10.11606/T.43.2017.tde-18092017-095345, eq. (3.7), printed p. 40. It gives s_low = A1 E^0.45 and s_high = (A2/E) ln(1 + A3/E + A4 E) (E in reduced units as defined there), attributes the form to Varelas and Biersack (1970), refined by Andersen and Ziegler (1977), and says A1 to A4 are fitted to experimental data.
• P. de Vera et al., arXiv:2608.15368 (2026), eq. (9), p. 3: the same s = s_low s_high / (s_low + s_high), citing Varelas and Biersack, and Andersen and Ziegler.

NISTIR 4999 (Berger 1992), section 3.5.2, p. 8, also mentions the “fitting formula of Varelas and Biersack (1970)” used with ICRU 49 coefficients (no equation given).

The original papers were not consulted (closed access; Unpaywall and Semantic Scholar report no open copy, ScienceDirect returned 403):

• C. Varelas and J. P. Biersack, Nucl. Instrum. Methods 79, 213 (1970), doi:10.1016/0029-554X(70)90141-2.
• J. P. Biersack and L. G. Haggmark, Nucl. Instrum. Methods 174, 257 (1980), doi:10.1016/0029-554X(80)90440-1.

Citation correction: “Biersack and Varelas, NIM 194, 93 (1982)” is not a Biersack-Varelas paper. Crossref resolves NIM 194, 93-100 (1982) to Biersack and Ziegler, “Refined universal potentials in atomic collisions” (doi:10.1016/0029-554X(82)90496-7). It is not a source for this join.

Why it is not implemented

The published s_high = (A2/E) ln(1 + A3/E + A4 E) is a fully fitted functional form (de Vera et al. eqs. (10)-(11), following ICRU 49; A1 to A4 are all fitting parameters) built so that it never crosses zero. For positive coefficients the 1 + keeps the logarithm’s argument above 1, so s_high > 0 at every E whatever A3 is. The harmonic join tends to s_low at low energy even with A3 = 0: s_high then tends to the finite value A2 A4 while s_low = A1 E^0.45 goes to 0. The fitted A3/E term only sets how fast s_high grows as E goes to 0. The coefficients come from the Andersen-Ziegler / ICRU 49 fits, which the clean-room policy (CONTRIBUTING.md) does not allow.

lindhard’s only sourced high-energy branch is BetheBloch, which has no such structure: its bracket crosses zero inside the crossover region. As S_high goes to 0+, the join S_low S_high / (S_low + S_high) goes to 0, not to S_low. Evaluated with LindhardScharff and BetheBloch defaults (Bloch term on, I = 10 eV Z2), in units of 1e-15 eV cm^2/atom (computed by a Python mirror of the formulas, not by the Rust models):

Bethe-Bloch is negative for H in Si below roughly 90 keV, for P in Si up to about 400 keV/u, and for He in Au up to about 450 keV/u. Neither workaround is acceptable. Falling back to S_low below the Bethe zero makes S jump from about 0 to the full LS value (about 27 for H in Si near 90 keV), which breaks continuity. Propagating NotApplicable leaves the model undefined over the whole low-energy regime, so ions slowing down through it cannot be transported. Restricting the join to where Bethe-Bloch is positive violates the low-energy limit and gives badly wrong stopping near the zero. A fix would need a high branch that stays positive at low E. The options are the fitted A1 to A4 (Tier C data), or a regularisation of our own, such as grafting the ln(1 + ...) structure onto Bethe by identifying A2 and A4 with Bethe quantities and setting A3 = 0. That identification is our own construction, which neither source makes. Neither option is admissible.

When it could be revisited

• A published, coefficient-free high-energy branch that stays positive below the Bethe zero (for example a Bethe form with a closed-form, non-tabulated shell correction; the review of the omitted correction terms under #41, now closed, is recorded below) is found and can be cited to an equation.
• The A1 to A4 coefficients become available from a source whose licence permits use, or are derived from our own fits to data that may be used.
• The primary papers are read and show a high branch that stays positive without fitted coefficients, or that ties the ln(1 + ...) form to Bethe quantities in a way that can be cited.

Literature findings (issue #41, closed)

Each of the four omitted items was reviewed for a closed form whose coefficients come from a paper itself. The review was from the contributor’s knowledge of the literature; no paper text was available to check digits against, so every statement below about a paper’s content is “as recalled” and the declines are on that basis. Nothing was taken from ICRU/SRIM/NIST tables.

• Barkas L1: declined. Ashley, Ritchie & Brandt (Phys. Rev. B 5, 2393 (1972)) give L1 = F(b / x^(1/2)) / (Z2^(1/2) x^(3/2)) with F a function that is tabulated (and later fitted to it) rather than closed-form. Lindhard’s 1976 treatment of the Barkas effect (NIM 132, 1) gives asymptotic forms from a harmonic-oscillator model, whose numerical coefficients we could not reproduce with confidence from memory. Jackson & McCarthy (Phys. Rev. B 6, 4131 (1972)) and later fits carry fitted coefficients that cannot be verified here. No form with verifiable coefficients exists to us; an invented one is not acceptable.
• Shell correction: declined. Walske and Bichsel shell corrections come from hydrogenic-shell calculations published as tables in η = βγ and I; the compact empirical forms (e.g. the ICRU 49 expression) are fits with coefficients tuned to data and tables, which Tier C of CONTRIBUTING.md excludes. A first-principles computation (hydrogenic shell sums) is a research task of its own, not a closed form.
• Density effect: partly implemented. The Sternheimer parameter sets (x0, x1, a, m, C̄) are tables, which we do not ingest. What is admissible is the dielectric formulation itself, computed from a model; we implement its one-oscillator specialisation (see density_effect_single_oscillator). Its high-energy limit is exact, but it switches on too late to matter in range. A multi-oscillator version needs per-material oscillator strengths, which are tabulated material data. Callers wanting a realistic δ in the range 1 to 100 GeV/u should supply one from a cited source via density_delta.
• Chu / Yang-O’Connor-Wang straggling: declined. Chu’s correction (Phys. Rev. A 13, 2057 (1976)) is built from Hartree-Fock-Slater charge densities and is published as tables; Yang, O’Connor & Wang (NIM B 61, 149 (1991)) is a fit whose coefficients are fitted to data tables. Neither has a closed form with paper-sourced coefficients that we could verify.

Deviations from the issue text

• The Lindhard-Scharff k_L quoted in the issue omits the factor ξ_e = Z1^(1/6). It is included here, because with it the reduced form k_L ε^(1/2) agrees with the dimensional LS expression to better than 1 % for all tested ion/target pairs (and without it Z1 > 1 disagrees by Z1^(1/6)).
• The Oen-Robinson local-loss constants are flagged unverified in data-provenance.md.

References

• J. Lindhard and M. Scharff, Phys. Rev. 124, 128 (1961).
• O. S. Oen and M. T. Robinson, Nucl. Instrum. Methods 132, 647 (1976).
• H. Bethe, Ann. Phys. 397, 325 (1930); H. Bethe, Z. Phys. 76, 293 (1932).
• F. Bloch, Ann. Phys. 408, 285 (1933).
• U. Fano, Ann. Rev. Nucl. Sci. 13, 1 (1963).
• W. H. Bragg and R. Kleeman, Phil. Mag. 10, 318 (1905).
• E. Fermi, Phys. Rev. 57, 485 (1940).
• R. M. Sternheimer, Phys. Rev. 88, 851 (1952).
• N. Bohr, Mat. Fys. Medd. Dan. Vid. Selsk. 18 (8) (1948).
• M. V. Moro, PhD thesis, Universidade de São Paulo (2017), doi:10.11606/T.43.2017.tde-18092017-095345.
• P. de Vera et al., arXiv:2608.15368 (2026).
• M. J. Berger, NISTIR 4999 (1992).
• C. Varelas and J. P. Biersack, Nucl. Instrum. Methods 79, 213 (1970), doi:10.1016/0029-554X(70)90141-2 (not consulted).
• J. P. Biersack and L. G. Haggmark, Nucl. Instrum. Methods 174, 257 (1980), doi:10.1016/0029-554X(80)90440-1 (not consulted for this join).
• J. P. Biersack and J. F. Ziegler, Nucl. Instrum. Methods 194, 93 (1982), doi:10.1016/0029-554X(82)90496-7.
• J. C. Ashley, R. H. Ritchie, W. Brandt, Phys. Rev. B 5, 2393 (1972).
• J. Lindhard, Nucl. Instrum. Methods 132, 1 (1976).
• J. D. Jackson and R. L. McCarthy, Phys. Rev. B 6, 4131 (1972).
• W. K. Chu, Phys. Rev. A 13, 2057 (1976).
• Q. Yang, D. J. O’Connor, Z. Wang, Nucl. Instrum. Methods B 61, 149 (1991).

Stopping residuals against measurements

Each electronic-stopping model of the stopping chapter compared with every measured stopping dataset committed under validation/data/stopping/ (format: lindhard-stopping-data/1, docs/stopping-data-format.md). Only data admitted by docs/data-provenance.md are in that directory; no NIST, IAEA or SRIM-derived value is read or shown.

Regenerate this page with

cargo run -p lindhard-plots --bin lindhard-stopping-report

and check it is current with the same command and -- --check. The output is a pure function of the committed datasets and the models, so running it twice gives the same bytes.

How to read the tables

• A residual is model - measured; positive means the model is above the measurement. Relative residuals are model / measured - 1.
• Stopping is a cross section per target atom (per average atom for a compound), in eV·10⁻¹⁵ cm². Mean and rms are the mean and root mean square of the absolute residual in that unit; rel. mean and rel. rms are those of the relative residual.
• A point the model cannot evaluate is skipped: it is counted in the skipped column, listed with the model’s reason, and is in no statistic.
• Where a system is one ion on one element, a per-pair correction of each model (a Chebyshev polynomial of degree 2 in ln E for the log of the measured-to-model ratio, lindhard::ion::stopping::correction) is fitted to the training points. 3 interior points chosen with seed 20261015 are held out of the fit. The corrected and uncorrected statistics on the held-out points are shown side by side; whether the correction helps is stated either way. The correction applies only inside the energy span of the training points.

Datasets

bader1956-p-cu

File: bader1956_p_cu.csv.

Stopping cross section of Cu for protons, 50 to 600 keV, as printed in Table I of Bader et al. (1956); values from the authors’ smooth curves through their measured points

Source: M. Bader, R. E. Pixley, F. S. Mozer and W. Whaling, Phys. Rev. 103, 32 (1956)

Note: the paper states no uncertainty for these values, so sigma is 0 (none stated)

H (1.008 u) on Cu

• Points: 15. Energy span: 50 keV to 600 keV.
• Stated uncertainty: none (sigma 0 on every point).
• Method: backscattering of protons from thin evaporated Cu layers (smoothed curve).
• Citation: M. Bader, R. E. Pixley, F. S. Mozer and W. Whaling, Phys. Rev. 103, 32 (1956), doi:10.1103/PhysRev.103.32, Table I (Cu column).

Skipped points (in no statistic above):

• bethe-bloch: line 7 (50 keV) skipped: bethe-bloch is not applicable at 50000 eV for this ion.
• bethe-bloch: line 8 (75 keV) skipped: bethe-bloch is not applicable at 75000 eV for this ion.
• bethe-bloch: line 9 (100 keV) skipped: bethe-bloch is not applicable at 100000 eV for this ion.
• bethe-bloch: line 10 (125 keV) skipped: bethe-bloch is not applicable at 125000 eV for this ion.
• bethe-bloch: line 11 (150 keV) skipped: bethe-bloch is not applicable at 150000 eV for this ion.

Held-out residuals, uncorrected and corrected

• bethe-bloch: no correction fitted (bethe-bloch is not applicable at 50000 eV for this ion).

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

srim-org-h-c8h8-1978-port

File: srim_org_h_c8h8_1978-Port.csv.

H in C8H8, 12 points of one paper (Porter, L. E. Naylor, H. Duder, J. C., Nucl. Inst. Methods, 155, 25-28 (1978)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on C8H8

• Points: 12. Energy span: 2.259 MeV to 5.598 MeV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Porter, L. E. Naylor, H. Duder, J. C., Nucl. Inst. Methods, 155, 25-28 (1978).

No per-pair correction is fitted: the target is a compound.

srim-org-h-c8h8-1980-bedn

File: srim_org_h_c8h8_1980-Bedn.csv.

H in C8H8, 7 points of one paper (Bednyakov, A. A. Bulgakov, Y. V. Nikolaev, V. S. Chernov, V. L., Sov. Phys., JETP 51, 954 (1980)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on C8H8

• Points: 7. Energy span: 120.637 keV to 480.833 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Bednyakov, A. A. Bulgakov, Y. V. Nikolaev, V. S. Chernov, V. L., Sov. Phys., JETP 51, 954 (1980).

No per-pair correction is fitted: the target is a compound.

srim-org-h-c8h8-1981-bedn

File: srim_org_h_c8h8_1981-Bedn.csv.

H in C8H8, 9 points of one paper (Bednyakov, A. A. Bulgakov, Y. V. Nikolaev, V. S. Chernov, V. L., Phys. Stat. Sol. A, 68, 187 (1981)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on C8H8

• Points: 9. Energy span: 70.145 keV to 501.695 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Bednyakov, A. A. Bulgakov, Y. V. Nikolaev, V. S. Chernov, V. L., Phys. Stat. Sol. A, 68, 187 (1981).

No per-pair correction is fitted: the target is a compound.

srim-org-h-h2o-solid-1952-wenz

File: srim_org_h_h2o-solid_1952-Wenz.csv.

H in H2O (solid), 7 points of one paper (Wenzel, W. A. Whaling, W., Phys. Rev., 87, 499-503 (1952)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on H2O

• Points: 7. Energy span: 10.32 keV to 45.171 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Wenzel, W. A. Whaling, W., Phys. Rev., 87, 499-503 (1952).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (10.32 keV) skipped: bethe-bloch is not applicable at 10320.13 eV for this ion.
• bethe-bloch: line 6 (19.421 keV) skipped: bethe-bloch is not applicable at 19420.79 eV for this ion.
• bethe-bloch: line 7 (20.489 keV) skipped: bethe-bloch is not applicable at 20489.08 eV for this ion.
• bethe-bloch: line 8 (20.489 keV) skipped: bethe-bloch is not applicable at 20489.08 eV for this ion.
• bethe-bloch: line 9 (35.868 keV) skipped: bethe-bloch is not applicable at 35868.49 eV for this ion.
• bethe-bloch: line 10 (41.049 keV) skipped: bethe-bloch is not applicable at 41048.71 eV for this ion.
• bethe-bloch: line 11 (45.171 keV) skipped: bethe-bloch is not applicable at 45170.72 eV for this ion.

No per-pair correction is fitted: the target is a compound.

srim-org-h-h2o-solid-1994-baue

File: srim_org_h_h2o-solid_1994-Baue.csv.

H in H2O (solid), 10 points of one paper (Bauer, P. Kaferbock, W. Necas, V., Nucl. Inst. Methods, B93, 132-136 (1994)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on H2O

• Points: 10. Energy span: 10.35 keV to 308.798 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Bauer, P. Kaferbock, W. Necas, V., Nucl. Inst. Methods, B93, 132-136 (1994).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (10.35 keV) skipped: bethe-bloch is not applicable at 10350.36 eV for this ion.
• bethe-bloch: line 6 (30.789 keV) skipped: bethe-bloch is not applicable at 30789.05 eV for this ion.
• bethe-bloch: line 7 (38.892 keV) skipped: bethe-bloch is not applicable at 38891.97 eV for this ion.
• bethe-bloch: line 8 (41.049 keV) skipped: bethe-bloch is not applicable at 41048.71 eV for this ion.
• bethe-bloch: line 9 (52.155 keV) skipped: bethe-bloch is not applicable at 52154.94 eV for this ion.

No per-pair correction is fitted: the target is a compound.

srim-org-h-h2o-solid-1997-baue

File: srim_org_h_h2o-solid_1997-Baue.csv.

H in H2O (solid), 9 points of one paper (Bauer, P. Golser, R. Aumayr, F. Semrad, D. Arnau, A., Nucl. Inst. Methods, B 125 102-105 (1997)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on H2O

• Points: 9. Energy span: 18.393 keV to 308.798 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Bauer, P. Golser, R. Aumayr, F. Semrad, D. Arnau, A., Nucl. Inst. Methods, B 125 102-105 (1997).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (18.393 keV) skipped: bethe-bloch is not applicable at 18392.81 eV for this ion.
• bethe-bloch: line 6 (20.489 keV) skipped: bethe-bloch is not applicable at 20489.08 eV for this ion.
• bethe-bloch: line 7 (30.789 keV) skipped: bethe-bloch is not applicable at 30789.05 eV for this ion.

No per-pair correction is fitted: the target is a compound.

srim-org-h-h2o-vapour-1990-mitt

File: srim_org_h_h2o-vapour_1990-Mitt.csv.

H in H2O (vapour), 13 points of one paper (Mitterschiffthaler, C. Bauer, P., Nucl. Inst. Methods, B48, 58 (1990)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on H2O

• Points: 13. Energy span: 28.652 keV to 307.89 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Mitterschiffthaler, C. Bauer, P., Nucl. Inst. Methods, B48, 58 (1990).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (28.652 keV) skipped: bethe-bloch is not applicable at 28652.46 eV for this ion.
• bethe-bloch: line 6 (36.856 keV) skipped: bethe-bloch is not applicable at 36856.16 eV for this ion.

No per-pair correction is fitted: the target is a compound.

srim-org-h-si-1976-grah

File: srim_org_h_si_1976-Grah.csv.

H in Si, 11 points of one paper (Grahmann, G.; Kalbitzer, S., Nucl. Inst. Methods, 132, 119-23 (1976)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on Si

• Points: 11. Energy span: 2.272 keV to 57.223 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-digitize/2.
• Citation: Digitized from the srim.org compilation: Grahmann, G.; Kalbitzer, S., Nucl. Inst. Methods, 132, 119-23 (1976).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (2.272 keV) skipped: bethe-bloch is not applicable at 2272.041 eV for this ion.
• bethe-bloch: line 6 (5.073 keV) skipped: bethe-bloch is not applicable at 5072.786 eV for this ion.
• bethe-bloch: line 7 (10.452 keV) skipped: bethe-bloch is not applicable at 10452.15 eV for this ion.
• bethe-bloch: line 8 (15.002 keV) skipped: bethe-bloch is not applicable at 15002.48 eV for this ion.
• bethe-bloch: line 9 (19.348 keV) skipped: bethe-bloch is not applicable at 19348.22 eV for this ion.
• bethe-bloch: line 10 (24.62 keV) skipped: bethe-bloch is not applicable at 24620.16 eV for this ion.
• bethe-bloch: line 11 (28.91 keV) skipped: bethe-bloch is not applicable at 28910.47 eV for this ion.
• bethe-bloch: line 12 (34.87 keV) skipped: bethe-bloch is not applicable at 34869.74 eV for this ion.
• bethe-bloch: line 13 (41.499 keV) skipped: bethe-bloch is not applicable at 41499.21 eV for this ion.
• bethe-bloch: line 14 (48.73 keV) skipped: bethe-bloch is not applicable at 48730.35 eV for this ion.
• bethe-bloch: line 15 (57.223 keV) skipped: bethe-bloch is not applicable at 57223.3 eV for this ion.

Held-out residuals, uncorrected and corrected

• bethe-bloch: no correction fitted (bethe-bloch is not applicable at 2272.041 eV for this ion).

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

srim-org-h-si-1977-cemb

File: srim_org_h_si_1977-Cemb.csv.

H in Si, 15 points of one paper (Cembali, F.; Zignani, F., Rad. Effects, 31, 169-173 (1977)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on Si

• Points: 15. Energy span: 79.969 keV to 270.389 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-digitize/2.
• Citation: Digitized from the srim.org compilation: Cembali, F.; Zignani, F., Rad. Effects, 31, 169-173 (1977).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (79.969 keV) skipped: bethe-bloch is not applicable at 79968.9 eV for this ion.

Held-out residuals, uncorrected and corrected

• bethe-bloch: no correction fitted (bethe-bloch is not applicable at 79968.9 eV for this ion).

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

srim-org-h-si-1988-mert

File: srim_org_h_si_1988-Mert.csv.

H in Si, 11 points of one paper (Mertens, P.; Bauer, P., Nucl. Inst. Methods, B33, 133 (1988)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the energy uncertainty (relative, in the compilation file) is not carried

H (1.008 u) on Si

• Points: 11. Energy span: 48.082 keV to 420.586 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-digitize/2.
• Citation: Digitized from the srim.org compilation: Mertens, P.; Bauer, P., Nucl. Inst. Methods, B33, 133 (1988).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (48.082 keV) skipped: bethe-bloch is not applicable at 48082.32 eV for this ion.
• bethe-bloch: line 6 (61.184 keV) skipped: bethe-bloch is not applicable at 61184.05 eV for this ion.
• bethe-bloch: line 7 (75.8 keV) skipped: bethe-bloch is not applicable at 75799.53 eV for this ion.

Held-out residuals, uncorrected and corrected

• bethe-bloch: no correction fitted (bethe-bloch is not applicable at 48082.32 eV for this ion).

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

srim-org-he-c8h8-1980-bedn

File: srim_org_he_c8h8_1980-Bedn.csv.

He in C8H8, 9 points of one paper (Bednyakov, A. A. Bulgakov, Y. V. Nikolaev, V. S. Chernov, V. L., Sov. Phys., JETP 51, 954 (1980)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

He (4.003 u) on C8H8

• Points: 9. Energy span: 109.351 keV to 349.267 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Bednyakov, A. A. Bulgakov, Y. V. Nikolaev, V. S. Chernov, V. L., Sov. Phys., JETP 51, 954 (1980).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (109.351 keV) skipped: bethe-bloch is not applicable at 109351.1 eV for this ion.
• bethe-bloch: line 6 (149.697 keV) skipped: bethe-bloch is not applicable at 149697.4 eV for this ion.
• bethe-bloch: line 7 (179.117 keV) skipped: bethe-bloch is not applicable at 179116.5 eV for this ion.
• bethe-bloch: line 8 (198.849 keV) skipped: bethe-bloch is not applicable at 198849.3 eV for this ion.
• bethe-bloch: line 9 (209.016 keV) skipped: bethe-bloch is not applicable at 209015.9 eV for this ion.
• bethe-bloch: line 10 (228.629 keV) skipped: bethe-bloch is not applicable at 228628.7 eV for this ion.

No per-pair correction is fitted: the target is a compound.

srim-org-he-c8h8-1981-bedn

File: srim_org_he_c8h8_1981-Bedn.csv.

He in C8H8, 8 points of one paper (Bednyakov, A. A. Bulgakov, Y. V. Nikolaev, V. S. Chernov, V. L., Phys. Stat. Sol. A, 68, 187 (1981)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the stopping is per average atom of the formula; the energy uncertainty (relative, in the compilation file) is not carried

He (4.003 u) on C8H8

• Points: 8. Energy span: 239.155 keV to 999.85 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-figures/1.
• Citation: Digitized from the srim.org compilation: Bednyakov, A. A. Bulgakov, Y. V. Nikolaev, V. S. Chernov, V. L., Phys. Stat. Sol. A, 68, 187 (1981).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (239.155 keV) skipped: bethe-bloch is not applicable at 239155.5 eV for this ion.

No per-pair correction is fitted: the target is a compound.

srim-org-he-cu-1975-nomu

File: srim_org_he_cu_1975-Nomu.csv.

He in Cu, 14 points of one paper (Nomura, A.; Kiyono, S., J. Phys. D: Appl. Phys., 8, 1551-59 (1975) (C) 2011 SRIM.org Experimental Stopping and Range Papers Page 9 of 23 Citations for Ion = H , Target = C Pub. Citation Year Authors, Title, Journal Citation and Comments Numb) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the energy uncertainty (relative, in the compilation file) is not carried

He (4.003 u) on Cu

• Points: 14. Energy span: 5.519 keV to 11.525 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-digitize/2.
• Citation: Digitized from the srim.org compilation: Nomura, A.; Kiyono, S., J. Phys. D: Appl. Phys., 8, 1551-59 (1975) (C) 2011 SRIM.org Experimental Stopping and Range Papers Page 9 of 23 Citations for Ion = H , Target = C Pub. Citation Year Authors, Title, Journal Citation and Comments Numb.

Skipped points (in no statistic above):

• bethe-bloch: line 5 (5.519 keV) skipped: bethe-bloch is not applicable at 5519.189 eV for this ion.
• bethe-bloch: line 6 (5.981 keV) skipped: bethe-bloch is not applicable at 5980.689 eV for this ion.
• bethe-bloch: line 7 (7.213 keV) skipped: bethe-bloch is not applicable at 7213.491 eV for this ion.
• bethe-bloch: line 8 (7.713 keV) skipped: bethe-bloch is not applicable at 7713.016 eV for this ion.
• bethe-bloch: line 9 (7.922 keV) skipped: bethe-bloch is not applicable at 7922.352 eV for this ion.
• bethe-bloch: line 10 (8.471 keV) skipped: bethe-bloch is not applicable at 8470.709 eV for this ion.
• bethe-bloch: line 11 (8.7 keV) skipped: bethe-bloch is not applicable at 8700.458 eV for this ion.
• bethe-bloch: line 12 (9.303 keV) skipped: bethe-bloch is not applicable at 9302.85 eV for this ion.
• bethe-bloch: line 13 (9.303 keV) skipped: bethe-bloch is not applicable at 9302.85 eV for this ion.
• bethe-bloch: line 14 (9.555 keV) skipped: bethe-bloch is not applicable at 9555.414 eV for this ion.
• bethe-bloch: line 15 (9.947 keV) skipped: bethe-bloch is not applicable at 9946.869 eV for this ion.
• bethe-bloch: line 16 (10.636 keV) skipped: bethe-bloch is not applicable at 10635.72 eV for this ion.
• bethe-bloch: line 17 (11.22 keV) skipped: bethe-bloch is not applicable at 11220.5 eV for this ion.
• bethe-bloch: line 18 (11.525 keV) skipped: bethe-bloch is not applicable at 11525.1 eV for this ion.

Held-out residuals, uncorrected and corrected

• bethe-bloch: no correction fitted (bethe-bloch is not applicable at 5519.189 eV for this ion).

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

srim-org-he-cu-1977-ande3

File: srim_org_he_cu_1977-Ande3.csv.

He in Cu, 17 points of one paper (Andersen, H. H.; Bak, J. F.; Knudsen, H.; Moller-Petersen, P.; Nielsen, B. R., Nucl. Inst. Methods, 140, 537-540 (1977)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the energy uncertainty (relative, in the compilation file) is not carried

He (4.003 u) on Cu

• Points: 17. Energy span: 6.395 MeV to 19.953 MeV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-digitize/2.
• Citation: Digitized from the srim.org compilation: Andersen, H. H.; Bak, J. F.; Knudsen, H.; Moller-Petersen, P.; Nielsen, B. R., Nucl. Inst. Methods, 140, 537-540 (1977).

Held-out residuals, uncorrected and corrected

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

Per point, for bethe-bloch: held-out points are marked held out; the others trained the fit.

srim-org-he-si-1973-lin

File: srim_org_he_si_1973-Lin.csv.

He in Si, 38 points of one paper (Lin, W. K.; Olson, H. G.; Powers, D., J. Appl. Phys., 44, 3631-34 (1973) (C) 2011 SRIM.org Experimental Stopping and Range Papers Page 1 of 5 Stopping for Ion =He , Target = Si Pub. Citation Year Authors, Title, Journal Citation and Comments Numb) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the energy uncertainty (relative, in the compilation file) is not carried

He (4.003 u) on Si

• Points: 38. Energy span: 237.462 keV to 1.969 MeV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-digitize/2.
• Citation: Digitized from the srim.org compilation: Lin, W. K.; Olson, H. G.; Powers, D., J. Appl. Phys., 44, 3631-34 (1973) (C) 2011 SRIM.org Experimental Stopping and Range Papers Page 1 of 5 Stopping for Ion =He , Target = Si Pub. Citation Year Authors, Title, Journal Citation and Comments Numb.

Skipped points (in no statistic above):

• bethe-bloch: line 5 (237.462 keV) skipped: bethe-bloch is not applicable at 237462.4 eV for this ion.
• bethe-bloch: line 6 (247.197 keV) skipped: bethe-bloch is not applicable at 247196.8 eV for this ion.
• bethe-bloch: line 7 (327.441 keV) skipped: bethe-bloch is not applicable at 327440.9 eV for this ion.
• bethe-bloch: line 8 (327.441 keV) skipped: bethe-bloch is not applicable at 327440.9 eV for this ion.
• bethe-bloch: line 9 (336.327 keV) skipped: bethe-bloch is not applicable at 336326.7 eV for this ion.
• bethe-bloch: line 10 (369.368 keV) skipped: bethe-bloch is not applicable at 369368.2 eV for this ion.
• bethe-bloch: line 11 (416.671 keV) skipped: bethe-bloch is not applicable at 416671 eV for this ion.
• bethe-bloch: line 12 (416.671 keV) skipped: bethe-bloch is not applicable at 416671 eV for this ion.
• bethe-bloch: line 13 (445.49 keV) skipped: bethe-bloch is not applicable at 445489.7 eV for this ion.
• bethe-bloch: line 14 (445.49 keV) skipped: bethe-bloch is not applicable at 445489.7 eV for this ion.
• bethe-bloch: line 15 (457.578 keV) skipped: bethe-bloch is not applicable at 457577.6 eV for this ion.

Held-out residuals, uncorrected and corrected

• bethe-bloch: no correction fitted (bethe-bloch is not applicable at 237462.4 eV for this ion).

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

srim-org-he-si-1973-zieg

File: srim_org_he_si_1973-Zieg.csv.

He in Si, 14 points of one paper (Ziegler, J. F.; Brodsky, M. H., J. Appl. Phys., 44, 188-96 (1973)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the energy uncertainty (relative, in the compilation file) is not carried

He (4.003 u) on Si

• Points: 14. Energy span: 400.26 keV to 2.751 MeV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-digitize/2.
• Citation: Digitized from the srim.org compilation: Ziegler, J. F.; Brodsky, M. H., J. Appl. Phys., 44, 188-96 (1973).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (400.26 keV) skipped: bethe-bloch is not applicable at 400260.3 eV for this ion.
• bethe-bloch: line 6 (400.26 keV) skipped: bethe-bloch is not applicable at 400260.3 eV for this ion.

Held-out residuals, uncorrected and corrected

• bethe-bloch: no correction fitted (bethe-bloch is not applicable at 400260.3 eV for this ion).

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

srim-org-he-si-1981-pear

File: srim_org_he_si_1981-Pear.csv.

He in Si, 14 points of one paper (Pearce, J. D.; Hart, R. R., J. Appl. Phys., 52, 5056 (1981)) read off the srim.org compilation plot; each has a half-pixel stopping uncertainty; the energy uncertainty (relative, in the compilation file) is not carried

He (4.003 u) on Si

• Points: 14. Energy span: 19.688 keV to 148.633 keV.
• Stated uncertainty: on every point.
• Method: digitized plot, srim_org-digitize/2.
• Citation: Digitized from the srim.org compilation: Pearce, J. D.; Hart, R. R., J. Appl. Phys., 52, 5056 (1981).

Skipped points (in no statistic above):

• bethe-bloch: line 5 (19.688 keV) skipped: bethe-bloch is not applicable at 19688 eV for this ion.
• bethe-bloch: line 6 (29.815 keV) skipped: bethe-bloch is not applicable at 29814.99 eV for this ion.
• bethe-bloch: line 7 (39.494 keV) skipped: bethe-bloch is not applicable at 39493.68 eV for this ion.
• bethe-bloch: line 8 (49.588 keV) skipped: bethe-bloch is not applicable at 49588.25 eV for this ion.
• bethe-bloch: line 9 (59.807 keV) skipped: bethe-bloch is not applicable at 59806.89 eV for this ion.
• bethe-bloch: line 10 (69.297 keV) skipped: bethe-bloch is not applicable at 69297.07 eV for this ion.
• bethe-bloch: line 11 (79.224 keV) skipped: bethe-bloch is not applicable at 79223.52 eV for this ion.
• bethe-bloch: line 12 (89.366 keV) skipped: bethe-bloch is not applicable at 89366.12 eV for this ion.
• bethe-bloch: line 13 (99.469 keV) skipped: bethe-bloch is not applicable at 99468.69 eV for this ion.
• bethe-bloch: line 14 (109.243 keV) skipped: bethe-bloch is not applicable at 109243 eV for this ion.
• bethe-bloch: line 15 (118.377 keV) skipped: bethe-bloch is not applicable at 118377 eV for this ion.
• bethe-bloch: line 16 (128.279 keV) skipped: bethe-bloch is not applicable at 128279.4 eV for this ion.
• bethe-bloch: line 17 (139.006 keV) skipped: bethe-bloch is not applicable at 139006.4 eV for this ion.
• bethe-bloch: line 18 (148.633 keV) skipped: bethe-bloch is not applicable at 148632.7 eV for this ion.

Held-out residuals, uncorrected and corrected

• bethe-bloch: no correction fitted (bethe-bloch is not applicable at 19688 eV for this ion).

Per point, for lindhard-scharff, oen-robinson, equipartition-ls-or: held-out points are marked held out; the others trained the fit.

BCA transport and free-flight conventions

Code: lindhard/src/ion/bca/ (Bca, BcaConfig, MeanFreePath, ElectronicLoss; kinematics.rs).

Model

The target is a stack of homogeneous, structureless (amorphous) layers, the first starting at depth \( x = 0 \); a semi-infinite substrate may close it. A moving atom alternates straight free flights and binary elastic collisions with target atoms, losing energy to electrons along each flight. This is the amorphous-target BCA of Biersack and Haggmark (1980), within the general BCA framework of Robinson and Torrens (1974) and the treatment of Eckstein (1991). One history:

1. The primary enters at the origin of the front face with the beam’s direction.
2. Free flight of length \( \tau\lambda \) (\( \lambda \) and the distribution of \( \tau \) depend on the free-path convention below), with the nonlocal electronic loss \( N S_e(E)\, s \) taken at the energy at the start of the segment (with Bragg additivity over the layer composition).
3. Collision with a partner drawn by stoichiometry, at impact parameter \( p = p_\mathrm{max}\sqrt{R} \) (uniform over a disc) and a uniform azimuth. The centre-of-mass angle comes from the scattering table.
4. Recoil. If the transfer \( T \) exceeds the partner’s displacement energy \( E_d \), the atom is displaced with energy \( T - E_b \) (the lattice binding \( E_b \) stays in the lattice); otherwise \( T \) stays in the lattice at the site (Biersack and Haggmark 1980; Eckstein 1991). Displaced atoms are followed in turn as full cascades.
5. A particle stops when its energy falls below its cutoff, or leaves through the front or back face if it can overcome the surface barrier; otherwise it is reflected specularly back into the target.

Kinematics

Classical elastic two-body kinematics (Goldstein 1980; Robinson and Torrens 1974), with \( \mu = M_1/M_2 \):

\[ \tan\psi = \frac{\sin\theta}{\cos\theta + \mu}, \qquad \varphi = \frac{\pi - \theta}{2}, \qquad T = \gamma E \sin^2\frac{\theta}{2}, \qquad \gamma = \frac{4 M_1 M_2}{(M_1 + M_2)^2}, \]

where \( \psi \) is the projectile’s laboratory deflection and \( \varphi \) the recoil’s, on the opposite azimuth.

Surface barrier

A planar barrier (Sigmund 1969; Eckstein 1991): a particle of energy \( E \) reaching a face with direction cosine \( c \) to the normal escapes if \( E c^2 > E_s \). Outside, its energy is \( E - E_s \), the parallel momentum is unchanged, and the normal component satisfies \( E’ c’^2 = E c^2 - E_s \) (refraction away from the normal). Target atoms use the surface binding energy \( E_s \) of the material at the face; the beam species uses primary_surface_binding_ev (default 0).

Free-path conventions

Constant

The flight length is the mean interatomic distance \( l = N^{-1/3} \) and impact parameters are uniform over a disc of radius \( p_\mathrm{max} = (\pi N^{2/3})^{-1/2} \), so \( N\pi p_\mathrm{max}^2 l = 1 \): each flight sweeps exactly one atom’s worth of target (Biersack and Haggmark 1980). So that every primary does not make its first collision at the same depth, the primary’s first flight is \( R\, l \) with \( R \) uniform in \( [0, 1) \); recoils start at an atom site and fly a full \( l \). Collisions with \( p > p_\mathrm{max} \) are dropped unless weak collisions add them.

Energy-dependent

Collisions that deflect by less than a minimum centre-of-mass angle \( \theta_\mathrm{min} \) (min_cm_angle_deg, required with this convention and only with it) are neglected. For each element \( i \), \( p_i(E) \) is the impact parameter at which the angle equals \( \theta_\mathrm{min} \), capped at the constant-convention \( p_\mathrm{max} \). The mean free path is

\[ \lambda = \frac{1}{N\pi \sum_i x_i p_i^2}, \]

flight lengths are exponential with that mean, the partner is drawn with probability \( x_i p_i^2 / \sum_j x_j p_j^2 \), and \( p = p_i\sqrt{R} \). This is the standard cross-section cut-off treatment of a Poisson collision process (Eckstein 1991). It reduces to the constant convention, with exponential flight lengths, when every \( p_i \) hits the cap. The nuclear loss of the neglected small-angle collisions is dropped, so choose \( \theta_\mathrm{min} \) small.

Layer boundaries

A flight is drawn as a dimensionless number of mean free paths \( \tau \). When it reaches an interface it is truncated there, the electronic loss for the truncated length is applied, the material is switched, and the flight continues in the new layer with the unused part \( \tau - s/\lambda \) (it is not redrawn). For exponential paths this is exact by memorylessness; for the constant path it means that splitting one layer into two of the same material changes nothing but rounding (a test checks this).

Weak collisions (optional)

With the constant free path every flight ends in one collision with \( p \le p_\mathrm{max} \), and the nuclear loss of collisions beyond \( p_\mathrm{max} \) is dropped while the electronic loss of the flight is charged in full. At cascade-tail energies that dropped part is large, so the electronic share of a cascade comes out too high. weak_collisions = K (0 to 3) adds the weak collisions of Möller and Eckstein (1988): before the hard collision, \( K \) collisions with partners at

\[ p_k = p_\mathrm{max}\sqrt{k + R}, \qquad k = 1, \dots, K, \]

each uniform over an annulus of area \( \pi p_\mathrm{max}^2 \), with its own partner and azimuth. A weak collision deflects the particle and takes its transfer, but never makes a recoil (its \( T \) stays in the lattice). A weak collision whose partner would lie in front of the front surface is skipped. Each collision is evaluated at the energy left after the previous one, which differs from the report, where all collisions of a step use the starting energy (see the ion::bca module docs for the measured effect). Weak collisions are available with the constant free path only. Default 0, because there is no single published convention and the option changes sputter yields substantially (see docs/validation.md).

Electronic loss

The local half is evaluated at the distance of closest approach of every collision, weak ones included; see Oen-Robinson.

Cutoffs and energy bookkeeping

primary_cutoff_ev and recoil_cutoff_ev are required: they have no defensible default. A recoil cutoff above the surface binding energies suppresses sputtering, so keep it below the smallest \( E_s \) when sputtering matters (Eckstein 1991). With follow_recoils = false a displaced atom stops where it was created.

Every energy change is subtracted from the particle and added to exactly one field of the energy budget (results.energy_budget_ev_per_ion), so each history conserves energy up to rounding; the largest per-history relative residual is reported. Electronic losses larger than the particle’s energy are clamped to it.

Randomness and determinism

Each primary has its own random stream, keyed on the run seed and the primary’s index; its recoils draw from the same stream, and cascades are followed from an explicit stack in a fixed order. Histories are grouped in chunks of a fixed size (64, independent of the thread count), which fixes the summation order of the tallies. The output is therefore bit-identical on any number of threads (lindhard/tests/determinism.rs).

Assumptions

• Amorphous, structureless targets: no crystal structure (no channeling).
• The target does not change with fluence (no dynamic composition).
• No refraction of the incident beam at the entrance surface (negligible at keV energies).
• Binary collisions only; many-body effects at low energy are represented only through the cutoffs, \( E_d \), \( E_b \) and \( E_s \).

Validity

The BCA is meant for energies well above the binding energies of the target. It follows cascade atoms down to a few eV, where its results depend on the cutoffs and binding energies chosen. Comparisons with other codes and with measured ranges and sputter yields, with the known deviations, are in docs/validation.md.

References

• J. P. Biersack and L. G. Haggmark, Nucl. Instrum. Methods 174, 257 (1980), doi:10.1016/0029-554X(80)90440-1.
• M. T. Robinson and I. M. Torrens, Phys. Rev. B 9, 5008 (1974).
• W. Eckstein, Computer Simulation of Ion-Solid Interactions (Springer, Berlin, 1991).
• H. Goldstein, Classical Mechanics, 2nd ed. (Addison-Wesley, 1980), sec. 3.11.
• P. Sigmund, Phys. Rev. 184, 383 (1969).
• W. Möller and W. Eckstein, TRIDYN - Binary collision simulation of atomic collisions and dynamic composition changes in solids, report IPP 9/64, Max-Planck-Institut für Plasmaphysik, Garching (1988), p. 14, eq. (26).
• O. S. Oen and M. T. Robinson, Nucl. Instrum. Methods 132, 647 (1976).

Displacement damage

Code: lindhard/src/ion/damage.rs (models), lindhard/src/tally/ (counts).

lindhard reports two different quantities side by side and never adds one to the other or substitutes one for the other:

• Model estimates of the number of displacements made by a primary knock-on atom (PKA) of a given energy: NRT and Kinchin-Pease, from the Lindhard partition of the PKA energy (results.damage.nrt).
• Event counts from the full-cascade simulation: vacancies, interstitials and replacements counted collision by collision as the BCA follows every recoil (results.damage.cascade).

Damage energy (Lindhard partition)

Of a PKA’s energy \( T \), the part eventually given to atomic motion is (Lindhard, Nielsen, Scharff and Thomsen 1963)

\[ T_\mathrm{dam} = \frac{T}{1 + k\, g(\varepsilon)}, \qquad g(\varepsilon) = 3.4008\, \varepsilon^{1/6} + 0.40244\, \varepsilon^{3/4} + \varepsilon, \]

with the analytic fit \( g \) of Robinson (1970) as adopted by the NRT standard, and

\[ \varepsilon = \frac{T a A_2}{Z_1 Z_2 e^2 (A_1 + A_2)}, \qquad a = 0.8853\, a_0 \left(Z_1^{2/3} + Z_2^{2/3}\right)^{-1/2}, \]

\[ k = \frac{0.0793\, Z_1^{2/3} Z_2^{1/2} (A_1 + A_2)^{3/2}} {\left(Z_1^{2/3} + Z_2^{2/3}\right)^{3/4} A_1^{3/2} A_2^{1/2}}, \]

the Lindhard reduced energy and the Lindhard-Scharff electronic stopping coefficient of the PKA (\( Z_1, A_1 \)) in the target (\( Z_2, A_2 \)). For a self-ion these reduce to the forms written in the NRT standard, \( \varepsilon = T / (86.931\, Z^{7/3}) \) (\( T \) in eV) and \( k = 0.1337\, Z^{1/6} (Z/A)^{1/2} \); a test checks this.

NRT and Kinchin-Pease

Norgett, Robinson and Torrens (1975):

\[ N_\mathrm{NRT} = \begin{cases} 0 & T_\mathrm{dam} < E_d, \\ 1 & E_d \le T_\mathrm{dam} < 2E_d/0.8, \\ 0.8\, T_\mathrm{dam} / (2 E_d) & T_\mathrm{dam} \ge 2E_d/0.8, \end{cases} \]

with \( E_d \) the displacement threshold and 0.8 the displacement efficiency the authors took from BCA simulations. The older Kinchin and Pease (1955) count is 0, 1, or \( T/(2E_d) \) with thresholds \( E_d \) and \( 2E_d \); the tallies evaluate it on the damage energy (the usual “modified” form).

NRT is not a cascade count. Stoller et al. (2013) discuss how BCA full-cascade vacancy counts and NRT estimates differ, and recommend that a displacement dose comparable with the NRT standard be computed from the damage energy with the NRT formula, not taken from the full-cascade vacancy count. NRT also overestimates the number of defects that survive in-cascade recombination; it is a standard exposure unit, not a prediction of surviving defects.

Cascade counts

During transport a target atom is displaced when it receives \( T > E_d \) (BCA transport, step 4). The counts are defined in lindhard::tally::CascadeDefects:

• Replacement: a moving atom comes to rest, immediately after the collision in which it displaced an atom of its own element, at that atom’s site, so it fills the site. A particle stops when its energy falls below its cutoff, so this count depends on the cutoffs (a recoil cutoff near \( E_d \) gives the most replacements); with follow_recoils = false it is always 0.
• Vacancies = displacements minus replacements.
• Interstitials: recoils that came to rest in the target and did not fill a site. Implanted beam particles are not counted here (they are in the range tally).

These are kept in total, per layer and as a depth profile (damage_profile.csv).

Selecting it

Always on: both are computed in every run. The inputs are the per-element energies: \( E_d \) (e_d_ev), \( E_b \) (e_b_ev) and \( E_s \) (e_s_ev), set per material or for an element in every layer with [physics.energies.<symbol>]. There are no model variants to choose.

Assumptions

• Compounds (this crate’s own convention). NRT and the Lindhard partition are defined for a monatomic target. For a layer with several elements the tally uses the PKA’s own \( Z_1, A_1 \), the atom-fraction-weighted mean \( Z_2, A_2 \) of the layer the PKA starts in (non-integer in general), and the PKA element’s own \( E_d \) in that layer. No published source is claimed for this averaging; it reduces exactly to standard NRT for a monatomic layer. Treat NRT numbers for compounds as an exposure index comparable within lindhard, not as a value following any compound-target standard.
• The cascade counts depend directly on \( E_d \), \( E_b \) and the recoil cutoff, and on the amorphous-target BCA itself (no recombination, no crystal structure).

Validity

NRT is a convention for damage accounting, defined for monatomic targets. The cascade counts are BCA counts of displacement events, not of surviving defects.

Verification status

The self-ion reductions of \( \varepsilon \) and \( k \) are checked in a test. The default \( E_d \) values of the element table follow the ASTM E521 convention and are not yet verified against the standard (docs/data-provenance.md); Si has no default, so an input must choose one.

Electron beams: ionisation, lattice and the knock-on threshold

An electron run reports the deposited energy in two channels (lindhard::tally::electron). Ionisation is electronic excitation, the source of electron-hole pairs: the energy left in the solid at inelastic events, plus the residual energy of electrons that fell below the cutoff (a cutoff-dependent term, placed here by convention). Lattice is LO-phonon emission and polaron traps, nonzero only with the opt-in insulator channels. Ionisation plus lattice equals the deposited energy per layer and per grid cell. Pairs per incident electron are the ionisation energy per primary divided by a per-material \( E_\mathrm{pair} \) that the input must supply with a citation; there is no default. Elastic nuclear recoil energy is not tallied by the electron engine and is in neither channel.

The knock-on report is kinematics only. For a beam electron of kinetic energy \( E \) and a free atom of mass \( M \) at rest, head on,

\[ T_\mathrm{max} = \frac{2E(E + 2mc^2)}{Mc^2}, \]

and the atom can be displaced when \( T_\mathrm{max} \ge E_d \), which sets the threshold beam energy \( E_\mathrm{th} = -mc^2 + \sqrt{(mc^2)^2 + E_d Mc^2/2} \). It ignores lattice binding, thermal motion and crystal direction, and gives no displacement cross section or dpa. An element without \( E_d \) is reported as “not set”.

References

• J. Lindhard, V. Nielsen, M. Scharff, P. V. Thomsen, “Integral equations governing radiation effects”, Mat. Fys. Medd. Dan. Vid. Selsk. 33 (10) (1963).
• M. T. Robinson, in Nuclear Fusion Reactors (British Nuclear Energy Society, London, 1970), p. 364.
• M. J. Norgett, M. T. Robinson, I. M. Torrens, Nucl. Eng. Des. 33, 50 (1975).
• G. H. Kinchin and R. S. Pease, Rep. Prog. Phys. 18, 1 (1955).
• R. E. Stoller, M. B. Toloczko, G. S. Was, A. G. Certain, S. Dwaraknath, F. A. Garner, Nucl. Instrum. Methods B 310, 75 (2013).

Dynamic composition (bookkeeping only)

Code: lindhard/src/ion/dynamic.rs (CompositionGrid, Relaxation).

Status. This is the start of milestone M3 (targets that change with fluence). Only the bookkeeping exists: a grid of slabs whose composition callers change with explicit inventory deltas. The fluence loop, and the adapter that turns transport tallies into deltas, are not written yet, so no CLI run uses it and there is no TOML key.

Model

The finite layers of a target are held as slabs. Each slab stores an areal inventory \( A_i \) (atoms/m²) per element, kept sorted by \( Z \) so every derived quantity is deterministic; a slab built from a layer of atom density \( n_i \) and thickness \( t \) starts with \( A_i = n_i t \). An optional semi-infinite substrate backs the grid and never changes.

Inventories change only through CompositionGrid::apply. A retained primary adds one atom of its species where it came to rest; a recoil subtracts an atom where it was created and adds it where it stops (an escaping recoil is a loss). The update either succeeds completely or leaves the grid untouched.

Volume relaxation

After every update each slab’s thickness is recomputed from its inventory by one of two stated conventions:

For an element in its own solid, atomic_volume_from_density gives \( v = M / (N_A \rho) \). After relaxation the slab boundaries are rebuilt from the front surface at \( x = 0 \). Swelling moves the interior interfaces and the back face. With sputter erosion on, the sputtered atoms of each element are removed from the front slabs (the loss they would otherwise cause in the slab where they were displaced is cancelled, so nothing is removed twice), the thickness they occupied under the chosen convention is the recession of that step, and the grid is re-anchored so the surface is again \( x = 0 \); a depth plus the cumulative recession is the depth in the original frame. Erosion is off by default and then changes nothing.

Retention limit

Optionally, each implanted species Z gets a saturation limit: the largest atomic fraction \( f_Z \) of a slab it may make up. Dynamic-composition codes (TRIDYN, SDTrimSP; see References) describe a saturation concentration of the implanted species above which further atoms are not retained and are released again. The update here is this crate’s own statement of that idea; no numeric limit is built in, it is always given by the caller. Per accepted step, each slab is looked at with the composition it would have after the step. With \( A_Z \) of the limited species and \( O \) atoms/m² of everything else, the slab may hold at most

\[ A_\mathrm{max} = \frac{f_Z}{1 - f_Z}, O \]

(that is \( A_Z / (A_Z + O) \le f_Z \)), and the excess \( A_Z - A_\mathrm{max} \) is re-emitted: removed from the slab through the same inventory deltas as any other change, and it leaves the target. A slab holding only the limited species loses all of it. With several limited species in one slab they are met in ascending \( Z \), each seeing the amounts left by the earlier ones, so the limit of a later species is met exactly and an earlier one can end marginally below its own after a later removal.

The atom balance of each element over a step is \( \Delta A_Z^\mathrm{total} = \text{transport net} - \text{re-emitted} \), where the transport net is the retained primaries and the moved, sputtered and transmitted recoils under the conventions above. The re-emission is a pure function of the merged integer tally and the grid and draws no random numbers, so results stay bit-identical at any thread count. The model is off by default and then changes nothing (the CLI adds no columns either).

Assumptions

• The retention limit is a bookkeeping cap on a slab’s composition. It has no diffusion, no temperature or time dependence, no trapping chemistry and no record of where re-emitted atoms go; the saturation fraction is a free parameter, not a prediction. A limit reached within one step is applied to the step’s result, so steps small enough that no slab jumps far past the limit give a better approximation (the adaptive policy sees the removal).

• Ideal mixing is additivity of atomic volumes (a Vegard-type rule applied to atoms); it ignores chemistry, voids and amorphisation swelling, and a compound’s real density generally differs from its prediction.

• A fixed number density reproduces a known compound density but cannot respond to composition, and one value serves every slab.

• Neither is an equation of state with a pressure or phase model from the literature; both are this crate’s stated conventions, and the module claims no published source for them.

Validity

Bookkeeping for compositions that stay close to the phase whose volumes or density the caller supplies. The caller must give the volume of every species that occurs (none is inferred, in particular not for gases).

Verification status

Conservation, transactional updates and the relaxation conventions are covered by lindhard/tests/dynamic.rs. The fluence stepping loop (DynamicRun, with the tally-to-delta adapter) is covered by lindhard/tests/dynamic_run.rs: the low-fluence limit equals the static engine, refining the step size converges, steps continue one global random stream, results are bit-identical on 1, 2 and 8 threads, adaptive rejection consumes no indices, and inventory follows the event conventions. The retention limit is checked for the analytic plateau (a saturated slab holds exactly \( f_Z \) and none exceeds it), for the atom balance (implanted − sputtered − re-emitted = retained), for bit-identity at 1, 2 and 8 threads, and for leaving results unchanged when off or out of reach. There is no comparison with a measured dynamic profile yet.

The fluence loop

DynamicRun delivers the beam’s primaries in steps. A step runs n primaries on the current target, scales the integer atom counts of the events by the fluence one primary stands for, applies them to the grid and relaxes it. The target is held fixed within a step, so the step must be small enough that the composition does not change much inside it: the adaptive policy bounds the largest relative change of a slab per step and retries a too-large step with fewer ions from the same first primary. The front surface stays at x = 0; the time series reports the interface depths and total thickness measured from it, not a receding surface. The step-size policy is this crate’s own design, not a published scheme. In the CLI it is the [dynamic] table (docs/cli.md).

References

• W. Möller and W. Eckstein, Nucl. Instrum. Methods B 2, 814 (1984), and report IPP 9/64 (1988): TRIDYN, dynamic composition changes. Read for the idea of a saturation concentration only; nothing was ported.
• A. Mutzke et al., SDTrimSP Version 6.00, report IPP 2019-02 (2019): the same, for the saturation of implanted species.

No published source is claimed for the two relaxation conventions or for the exact form of the retention update (see above). The atomic volume uses the Avogadro constant of CODATA 2022 (P. J. Mohr, D. B. Newell, B. N. Taylor, E. Tiesinga, Rev. Mod. Phys. 97, 025002 (2025)) and the element densities of lindhard/src/elements.rs, whose sources are in docs/data-provenance.md.

History and literature of ion stopping

This chapter is a reading guide to how the theory of ion stopping and range was built, from the first alpha-particle experiments to the codes of the 1980s. It names the papers that matter, says in a sentence what each contributed, and links every model that lindhard implements to the manual page where it lives. It is written for a reader who wants to know where a model came from before trusting it. It is not a substitute for the papers.

Where this chapter comes from

The shape of the story, and many of the references, follow the historical review that James F. Ziegler wrote for The Stopping and Range of Ions in Solids [Z85] and later put on the SRIM web site, srim.org [S03]. That review, and the bibliography of measured stopping data that Ziegler compiled over decades, are the reason a field this large has a shared literature at all. The text below is our own summary, with every statement tied to a primary paper. The srim.org site is no longer maintained, so each page we cite is linked twice, to the live site and to a copy in the Internet Archive:

Prefer the archived copies: pages on the live site have since changed, and some content on them is unrelated to the subject. lindhard builds on this published literature as an independent implementation. It shares none of SRIM’s code, uses SRIM’s published data with citations where useful, and is not endorsed by its authors (see Conventions and the clean-room position). This chapter contains no stopping data, curves or figures from srim.org.

The Monte Carlo side of the lineage, from the 1947 neutron histories to the binary-collision codes, is told separately in docs/history.md. It cites some of the same papers; this chapter is about the stopping physics.

Vocabulary

A fast ion in matter loses energy two ways. Electronic stopping is the energy given to the target’s electrons. Nuclear stopping is the energy given to target atoms as a whole in elastic collisions. The rate of loss per unit path, \( dE/dx \), is the stopping power (older papers) or stopping force (newer ones); divided by the atomic density it is the stopping cross section. Stopping depends on the ion’s velocity much more than on its energy, which is why the history below is organised by velocity.

From radioactivity to the nuclear atom (1900–1915)

Soon after radioactivity was discovered, Marie Curie proposed that alpha rays are material projectiles that lose speed in matter [C00]. Bragg and Kleeman showed in 1905 that the range of alpha particles in a compound is roughly the sum of contributions from its atoms [BK05]. That is the origin of what lindhard calls Bragg additivity.

Geiger and Marsden found that a small fraction of alpha particles bounce back from a thin foil [GM09]. Rutherford showed that this needs a single close collision with a compact, heavy, charged core [R11], which is the nuclear atom and the scattering problem that the scattering integral still solves.

Thomson [T12] and Bohr [B13, B15] then wrote the first theories of the energy loss itself. Bohr treated target electrons as bound oscillators with the frequencies of optical spectra, separated nuclear from electronic losses, and saw that the unknown charge state of the moving ion limited his accounting. The problems he listed (energy transfer to a quantised medium, polarisation of the medium, the ion’s charge state, screened nuclear scattering, relativity) organise most of what follows.

Quantum mechanics and fast particles (1930–1934)

Bethe [Be30, Be32] and Bloch [Bl33a, Bl33b] redid the problem with quantum mechanics in the Born approximation. The result, the Bethe-Bloch formula, is the standard description of a fast, fully stripped ion, and the base of all high-energy stopping tables. Bethe and Heitler [BH34] extended the treatment to radiative losses of fast electrons. Fermi’s 1940 paper [F40] added the polarisation of the medium, the origin of the density effect. lindhard implements the formula, the effective charge options and the density effect in Bethe-Bloch and the density effect. For a survey of the relativistic regime see Ahlen [A80].

Fission fragments and the effective charge (1939–1941)

The discovery of fission [HS39] produced heavy, partially stripped ions with several MeV of energy, and a new question: how to apply a stopping theory built for bare charges to an ion carrying its own electrons. Bohr [B40, B41] proposed that an ion keeps the electrons that orbit more slowly than it moves, and used the Thomas-Fermi atom [Th27] to turn this into an effective charge that scales with velocity. Lamb [L40] argued from binding energy instead, and Knipp and Teller [KT41] used the idea to scale proton stopping to other ions. lindhard’s effective-charge options are on the Bethe-Bloch page.

Screened nuclear collisions

Nuclear stopping needs the potential between two atoms. Thomas-Fermi statistics [Th27] gave a universal screening function, and analytic approximations to it followed: Lenz and Jensen in 1932 [Le32, J32] and Molière in 1947 [M47]. Bohr [B48] argued for a screening length that combines the two nuclear charges, and Firsov [Fi58] obtained a different combination by fitting numerically computed potentials of colliding Thomas-Fermi atoms. Later potentials replaced the statistical atom by Hartree-Fock calculations: the Kr-C potential of Wilson, Haggmark and Biersack [WHB77], and the universal potential of Ziegler, Biersack and Littmark, fitted to many such calculations [BZ82, Z85]. lindhard offers each of these in Interatomic potentials, with the screening lengths in Screening lengths.

The slow-ion regime and the free electron gas (1947–1960)

For ions slower than the target’s electrons the loss is proportional to the ion’s velocity. Fermi and Teller derived this for a degenerate electron gas [FT47]. Lindhard’s 1954 dielectric treatment [Li54], with the independent work of Neufeld and Ritchie [NR55], put the electron gas on a full quantum footing; Lindhard and Scharff used it in 1961 to give the velocity-proportional stopping of keV ions [LS61]. That model is Lindhard-Scharff. Northcliffe [N60] collected measurements for many ions and targets into a ratio to proton stopping at the same velocity, from which an empirical effective charge follows. His review [N63] is a good entry point to the experimental side.

LSS theory (1961–1963)

Lindhard, Scharff and Schiøtt [LSS63] joined the nuclear and electronic parts into one theory with a universal reduced energy and reduced range, so that range distributions could be computed for any ion-target pair from a single curve. The reduced units, and the velocity-proportional electronic loss, are the starting point of lindhard’s Lindhard-Scharff model and its screening-length conventions. LSS theory rests on the statistical atom, which is why it shows no shell structure. Later work kept its framework and replaced its ingredients [Z85].

Realistic atoms and computers (1970–1985)

Computers made it possible to replace statistical atoms by Hartree-Fock ones. Rousseau, Chu and Powers did this for electronic stopping [RCP71], and Wilson, Haggmark and Biersack for nuclear stopping [WHB77]. Brandt and Kitagawa [BK82] revised Bohr’s criterion for how many electrons an ion keeps: they compared the ion’s electron velocities with the Fermi velocity of the solid, not with the ion velocity alone. In the same period Oen and Robinson [OR76] gave a local, impact-parameter-dependent electronic loss, which lindhard offers, mixed with the Lindhard-Scharff loss, in Oen-Robinson and the equipartition mix.

On the transport side, Biersack and Haggmark’s TRIM [BH80] made binary-collision Monte Carlo practical for amorphous targets and introduced the analytic scattering-angle formula that the scattering integral keeps only as a cross-check. BCA transport and free-flight conventions describes the transport loop. Ziegler, Biersack and Littmark’s 1985 book [Z85] and the SRIM program [ZZB10] that grew from it are the best-known products of this period. For the earlier computer studies of ions in crystals, and the theory of channeling, see docs/history.md.

Straggling and compounds

The fluctuation of energy loss about its mean, straggling, was already treated by Bohr [B48] for free electrons; Fano’s review [Fa63] covers the corrections. lindhard’s version is on the Energy-loss straggling page. Stopping in compounds and mixtures starts from Bragg’s additivity rule [BK05] and corrects it for chemical binding; the srim.org page on compounds (linked above) discusses the corrections, and lindhard’s treatment is Bragg additivity. Ziegler [Z85] and srim.org describe a body of range and stopping data for many ion-target pairs. Data and their provenance are tracked separately, in docs/data-provenance.md.

Where to read next

• Bohr’s 1948 review [B48] is long but readable and covers most of what was known after 50 years.
• Northcliffe [N63] and Fano [Fa63] are the standard reviews of heavy-ion and light-ion stopping around the time of LSS.
• Whaling [W58] reviews the experiments of the 1950s. Ahlen [A80] covers relativistic ions.
• Sigmund’s monograph Particle Penetration and Radiation Effects [Si06] is a modern textbook treatment.
• The International Commission on Radiation Units and Measurements publishes evaluated tables and their rationale: Report 49 for protons and alpha particles (1993) [ICRU49], and Report 37 for electrons and positrons (1984) [ICRU37]. lindhard does not use their tables; the reports are listed here for the physics in their text.
• A tribute to the man who compiled much of this literature: James F. Ziegler (1935–2024).

Corrections to the srim.org reference list

The list at the end of the srim.org review has inconsistencies. We checked each reference below against the journal or publisher record (Crossref unless noted) and used the corrected form in this chapter.

References

Entries with a DOI were checked against Crossref. Entries marked (no DOI) are not in Crossref; their citation is given as in the sources named, and where we did not see the original we say so.

• [A80] S. P. Ahlen, “Theoretical and experimental aspects of the energy loss of relativistic heavily ionizing particles,” Rev. Mod. Phys. 52, 121 (1980), doi:10.1103/RevModPhys.52.121.
• [B13] N. Bohr, “On the theory of the decrease of velocity of moving electrified particles on passing through matter,” Phil. Mag. 25, 10 (1913), doi:10.1080/14786440108634305.
• [B15] N. Bohr, “On the decrease of velocity of swiftly moving electrified particles in passing through matter,” Phil. Mag. 30, 581 (1915), doi:10.1080/14786441008635432.
• [B40] N. Bohr, “Scattering and stopping of fission fragments,” Phys. Rev. 58, 654 (1940), doi:10.1103/PhysRev.58.654.
• [B41] N. Bohr, “Velocity-range relation for fission fragments,” Phys. Rev. 59, 270 (1941), doi:10.1103/PhysRev.59.270.
• [B48] N. Bohr, “The penetration of atomic particles through matter,” Mat. Fys. Medd. Dan. Vid. Selsk. 18, no. 8 (1948). (no DOI; reprinted in Niels Bohr Collected Works vol. 8, doi:10.1016/S1876-0503(08)70172-5.)
• [Be30] H. Bethe, “Zur Theorie des Durchgangs schneller Korpuskularstrahlen durch Materie,” Ann. Phys. 397, 325 (1930), doi:10.1002/andp.19303970303.
• [Be32] H. Bethe, “Bremsformel für Elektronen relativistischer Geschwindigkeit,” Z. Phys. 76, 293 (1932), doi:10.1007/BF01342532.
• [BH34] H. Bethe and W. Heitler, “On the stopping of fast particles and on the creation of positive electrons,” Proc. R. Soc. A 146, 83 (1934), doi:10.1098/rspa.1934.0140.
• [BH80] J. P. Biersack and L. G. Haggmark, “A Monte Carlo computer program for the transport of energetic ions in amorphous targets,” Nucl. Instrum. Methods 174, 257 (1980), doi:10.1016/0029-554X(80)90440-1.
• [Bl33a] F. Bloch, “Zur Bremsung rasch bewegter Teilchen beim Durchgang durch Materie,” Ann. Phys. 408, 285 (1933), doi:10.1002/andp.19334080303.
• [Bl33b] F. Bloch, “Bremsvermögen von Atomen mit mehreren Elektronen,” Z. Phys. 81, 363 (1933), doi:10.1007/BF01344553.
• [BK05] W. H. Bragg and R. Kleeman, “On the α particles of radium, and their loss of range in passing through various atoms and molecules,” Phil. Mag. 10, 318 (1905), doi:10.1080/14786440509463378.
• [BK82] W. Brandt and M. Kitagawa, “Effective stopping-power charges of swift ions in condensed matter,” Phys. Rev. B 25, 5631 (1982), doi:10.1103/PhysRevB.25.5631.
• [BZ82] J. P. Biersack and J. F. Ziegler, “Refined universal potentials in atomic collisions,” Nucl. Instrum. Methods 194, 93 (1982), doi:10.1016/0029-554X(82)90496-7.
• [C00] M. Curie, Comptes Rendus 130, 76 (1900). (no DOI; citation as given on srim.org, original not seen.)
• [F40] E. Fermi, “The ionization loss of energy in gases and in condensed materials,” Phys. Rev. 57, 485 (1940), doi:10.1103/PhysRev.57.485.
• [Fa63] U. Fano, “Penetration of protons, alpha particles, and mesons,” Ann. Rev. Nucl. Sci. 13, 1 (1963), doi:10.1146/annurev.ns.13.120163.000245.
• [Fi58] O. B. Firsov, Sov. Phys. JETP 6, 534 (1958). (no DOI; as cited on the screening-lengths page.)
• [FT47] E. Fermi and E. Teller, “The capture of negative mesotrons in matter,” Phys. Rev. 72, 399 (1947), doi:10.1103/PhysRev.72.399.
• [GM09] H. Geiger and E. Marsden, “On a diffuse reflection of the α-particles,” Proc. R. Soc. A 82, 495 (1909), doi:10.1098/rspa.1909.0054.
• [HS39] O. Hahn and F. Strassmann, “Über den Nachweis und das Verhalten der bei der Bestrahlung des Urans mittels Neutronen entstehenden Erdalkalimetalle,” Naturwissenschaften 27, 11 (1939), doi:10.1007/BF01488241.
• [ICRU37] ICRU, Stopping Powers for Electrons and Positrons, ICRU Report 37 (1984). (no DOI; title confirmed from the review by D. K. Brice, Nucl. Instrum. Methods B 12, 187 (1985), doi:10.1016/0168-583X(85)90718-9.)
• [ICRU49] ICRU, Stopping Powers and Ranges for Protons and Alpha Particles, ICRU Report 49 (1993). (no DOI; title confirmed from the review by J. O. Deasy, Med. Phys. 21, 709 (1994), doi:10.1118/1.597176.)
• [J32] H. Jensen, “Die Ladungsverteilung in Ionen und die Gitterkonstante des Rubidiumbromids nach der statistischen Methode,” Z. Phys. 77, 722 (1932), doi:10.1007/BF01342151.
• [KT41] J. Knipp and E. Teller, “On the energy loss of heavy ions,” Phys. Rev. 59, 659 (1941), doi:10.1103/PhysRev.59.659.
• [L40] W. E. Lamb, “Passage of uranium fission fragments through matter,” Phys. Rev. 58, 696 (1940), doi:10.1103/PhysRev.58.696.
• [Le32] W. Lenz, “Über die Anwendbarkeit der statistischen Methode auf Ionengitter,” Z. Phys. 77, 713 (1932), doi:10.1007/BF01342150.
• [Li54] J. Lindhard, “On the properties of a gas of charged particles,” Mat. Fys. Medd. Dan. Vid. Selsk. 28, no. 8 (1954). (no DOI.)
• [LS61] J. Lindhard and M. Scharff, “Energy dissipation by ions in the keV region,” Phys. Rev. 124, 128 (1961), doi:10.1103/PhysRev.124.128.
• [LSS63] J. Lindhard, M. Scharff, H. E. Schiøtt, “Range concepts and heavy ion ranges (Notes on atomic collisions, II),” Mat. Fys. Medd. Dan. Vid. Selsk. 33, no. 14 (1963). (no DOI.)
• [M47] G. Molière, “Theorie der Streuung schneller geladener Teilchen I. Einzelstreuung am abgeschirmten Coulomb-Feld,” Z. Naturforsch. A 2, 133 (1947), doi:10.1515/zna-1947-0302.
• [N60] L. C. Northcliffe, “Energy loss and effective charge of heavy ions in aluminum,” Phys. Rev. 120, 1744 (1960), doi:10.1103/PhysRev.120.1744.
• [N63] L. C. Northcliffe, “Passage of heavy ions through matter,” Ann. Rev. Nucl. Sci. 13, 67 (1963), doi:10.1146/annurev.ns.13.120163.000435.
• [NR55] J. Neufeld and R. H. Ritchie, “Passage of charged particles through plasma,” Phys. Rev. 98, 1632 (1955), doi:10.1103/PhysRev.98.1632.
• [OR76] O. S. Oen and M. T. Robinson, “Computer studies of the reflection of light ions from solids,” Nucl. Instrum. Methods 132, 647 (1976), doi:10.1016/0029-554X(76)90806-5.
• [R11] E. Rutherford, “The scattering of α and β particles by matter and the structure of the atom,” Phil. Mag. 21, 669 (1911), doi:10.1080/14786440508637080.
• [RCP71] C. C. Rousseau, W. K. Chu, D. Powers, “Calculations of stopping cross sections for 0.8- to 2.0-MeV alpha particles,” Phys. Rev. A 4, 1066 (1971), doi:10.1103/PhysRevA.4.1066.
• [S03] J. F. Ziegler, “Historical review,” srim.org, 2003 (the date of the chapter): live page, Internet Archive copy. The page says it is taken from [Z85].
• [Si06] P. Sigmund, Particle Penetration and Radiation Effects, Springer Series in Solid-State Sciences 151 (Springer, 2006), doi:10.1007/3-540-31718-X.
• [T12] J. J. Thomson, “Ionization by moving electrified particles,” Phil. Mag. 23, 449 (1912), doi:10.1080/14786440408637241.
• [Th27] L. H. Thomas, “The calculation of atomic fields,” Proc. Camb. Phil. Soc. 23, 542 (1927), doi:10.1017/S0305004100011683.
• [W58] W. Whaling, “The energy loss of charged particles in matter,” in Handbuch der Physik 34, 193 (Springer, 1958), doi:10.1007/978-3-642-45898-9_5.
• [WHB77] W. D. Wilson, L. G. Haggmark, J. P. Biersack, “Calculations of nuclear stopping, ranges, and straggling in the low-energy region,” Phys. Rev. B 15, 2458 (1977), doi:10.1103/PhysRevB.15.2458.
• [Z85] J. F. Ziegler, J. P. Biersack, U. Littmark, The Stopping and Range of Ions in Solids, vol. 1 of the series The Stopping and Ranges of Ions in Matter (Pergamon, New York, 1985). (no DOI; record checked in Open Library, which gives the title and publisher. The series title is written both with “Range” and “Ranges” in the sources we saw.)
• [ZZB10] J. F. Ziegler, M. D. Ziegler, J. P. Biersack, “SRIM – The stopping and range of ions in matter (2010),” Nucl. Instrum. Methods B 268, 1818 (2010), doi:10.1016/j.nimb.2010.02.091.

James F. Ziegler (1935–2024)

James F. Ziegler died on 14 May 2024, aged 89 [O24]. A great deal of what the field of ion stopping shares today, its compiled measurements, its bibliography, its reference books and its best-known simulation program, exists because of his work. lindhard builds on that legacy, as an independent implementation written from the published physics. This page records his professional life and work. It is not written or endorsed by Ziegler or his estate, and nothing here implies that they endorse lindhard.

Life and career

Ziegler was born on 20 April 1935 in Gary, Indiana [O24]. He studied at Yale University, where he completed a Ph.D. in physics in 1967 [B-JFZ]. The obituary adds that he began in electrical engineering, served in the US Navy, and studied physics at the University of Cincinnati before returning to Yale [O24].

He joined IBM Research in 1967 and stayed until 2000, mostly at Yorktown Heights, New York. He worked on ion channeling, ion implantation for making semiconductor devices, and radiation damage, and from 1975 to 2000 managed IBM’s radiation science group [B-JFZ]. He then taught as a professor at the US Naval Academy in Annapolis, Maryland [B-JFZ, O24].

He was elected a Fellow of the American Physical Society (1980) and a Fellow of the Institute of Electrical and Electronics Engineers (1996), and received a von Humboldt Senior Scientist prize (1987) for work on the theory of ion implantation [B-JFZ, O24].

Work

The Stopping and Ranges of Ions in Matter series. Between 1977 and 1985 Ziegler wrote or co-wrote a run of volumes that collected stopping powers and range distributions for ions in all elements, with Andersen (hydrogen), himself (helium), and Littmark (range distributions) [AZ77, Z78, LZ80, Z80]. The first volume of the series, with Biersack and Littmark, is The Stopping and Range of Ions in Solids (Pergamon, 1985) [ZBL85]. It contains the universal interatomic potential known as ZBL and the computational methods behind SRIM; lindhard’s ZBL potential is taken from the published description of it.

SRIM. With Biersack and other co-authors he developed the program, first TRIM [BH80] and then SRIM, that made Monte Carlo calculation of ion ranges routine in laboratories and fabs [ZZB10]. SRIM has been distributed since the mid-1980s. The SRIM site, srim.org, also carries his historical review of the field [S03]; see the history chapter.

Compiled measurements. He assembled, over decades, a bibliography of published measurements of stopping and range, and used the measurements it indexes in his stopping-power volumes [Z77]. Preserving that literature is part of why this chapter exists.

Cosmic rays and electronics. In 1979 Ziegler and Lanford showed that cosmic-ray secondaries can flip bits in computer memories at sea level [ZL79]. Ziegler returned to the subject for years afterwards [Z96], and the 1979 paper became a starting point for the study of soft errors in terrestrial electronics [B-JFZ].

Sources and links

• Obituary, The Block Island Times, 7 June 2024, on Legacy.com: O24. “New Shoreham, RI” there is the newspaper’s location, not his. This page does not use the family or guestbook sections of the obituary, which are for the reader to visit there.
• Professional biography on srim.org: live URL (not found at the time of writing), Internet Archive copy, 2024-08-10 [B-JFZ].
• srim.org: home page and historical review.

References

• [AZ77] H. H. Andersen and J. F. Ziegler, Hydrogen Stopping Powers and Ranges in All Elements (Pergamon, New York, 1977). (No DOI; title and series as listed in the srim.org review and Ziegler’s biography; no library record checked.)
• [BH80] J. P. Biersack and L. G. Haggmark, Nucl. Instrum. Methods 174, 257 (1980), doi:10.1016/0029-554X(80)90440-1.
• [B-JFZ] “James F. Ziegler, professional biography,” srim.org (Internet Archive copy of 2024-08-10, linked above).
• [LZ80] U. Littmark and J. F. Ziegler, Handbook of Range Distributions for Energetic Ions in All Elements (Pergamon, New York, 1980). (No DOI; as for [AZ77].)
• [O24] Obituary of James Ziegler, The Block Island Times, 7 June 2024. (The page returned an access-denied status to an automated fetch, so the facts taken from it were supplied by the project operator from a copy of the obituary, and those that also appear in [B-JFZ] were cross-checked there.)
• [S03] J. F. Ziegler, “Historical review,” srim.org; see the history chapter.
• [Z77] J. F. Ziegler, Bibliography of Experimental Range and Stopping Power Data (Pergamon, 1977), as listed in [B-JFZ]. (No DOI; no library record checked.)
• [Z78] J. F. Ziegler, Helium Stopping Powers and Ranges in All Elements (Pergamon, New York, 1978). (No DOI; as for [AZ77].)
• [Z80] J. F. Ziegler, Handbook of Stopping Cross-Sections for Energetic Ions in All Elements (Pergamon, New York, 1980). (No DOI; as for [AZ77].)
• [Z96] J. F. Ziegler, “Terrestrial cosmic rays,” IBM J. Res. Dev. 40, 19 (1996), doi:10.1147/rd.401.0019.
• [ZBL85] J. F. Ziegler, J. P. Biersack, U. Littmark, The Stopping and Range of Ions in Solids, vol. 1 of The Stopping and Ranges of Ions in Matter (Pergamon, New York, 1985); see the history chapter [Z85].
• [ZL79] J. F. Ziegler and W. A. Lanford, “Effect of cosmic rays on computer memories,” Science 206, 776 (1979), doi:10.1126/science.206.4420.776.
• [ZZB10] J. F. Ziegler, M. D. Ziegler, J. P. Biersack, “SRIM – The stopping and range of ions in matter (2010),” Nucl. Instrum. Methods B 268, 1818 (2010), doi:10.1016/j.nimb.2010.02.091.

Validation

What lindhard has been compared with, against which published reference, and how far off it is. The full record, with every row, tolerance and investigation, is docs/validation.md; this page summarises it. Where lindhard disagrees with a reference the disagreement is reported, not hidden. Provenance of every number in the tree: docs/data-provenance.md.

The numbers below are the level 1 results (full statistics, lindhard 0.0.1) as of 2026-11-01. Level 1 runs in CI on every pull request.

Level 1: analytic limits and internal consistency

Each row compares the engine with an equation or limit that can be evaluated independently. “Deviation” is the measured difference; “tolerance” is the bound the check asserts (three standard errors plus solver uncertainty plus a stated systematic allowance for Monte Carlo rows).

Results are byte-identical on 1, 2 and 8 threads (engine.determinism).

Known deviations

These are measured and not asserted away:

• ZBL universal fit, 0.3 < ε < 1000: the fit runs up to 5.5 % low against the quadrature, because its high-energy branch drops the screening constant. The quadrature agrees with the exact asymptote to 0.1 %.
• Range of light ions in a semi-infinite target: 10 keV B into Si comes out 3.97 % above the infinite-medium equation, because backscattered ions leave the target. A difference of definition, not of transport.
• Range against measurement: for B in amorphous Si the default stopping inputs (Lindhard–Scharff electronic, ZBL nuclear) give a projected range 24 to 36 % too long. Because the transport reproduces the LSS equation for the same inputs, the offset comes from those inputs.
• Sputter yield against displacement energy: Ar at 1 keV on Si gives 0.30 with E_d = 15 eV and 0.49 with E_d = E_s. Set E_d ≤ E_s when low-energy yields matter.
• Dynamic composition: Cu–Ni fluence-dependent targets are compared with measurement, and the engine does not reproduce the measured profile. They are compared, not validated.

Levels 2 and 3

Level 2 compares with other codes through oracle adapters (summaries are committed; the third-party outputs are not). Level 3 compares with measurement: stopping data, sputter yields, backscatter coefficients, secondary-electron yields and inelastic mean free paths. Each comparison has its data source, residuals and verdict in docs/validation.md, and the page also says where a comparison is missing. The residuals against stopping measurements are in the manual.

Status and roadmap

lindhard is early, pre-release software (version 0.0.1; no release is published yet, so it is built from source). State as of 2026-11-01:

The plan of record is WORK_PLAN.md and the design is docs/architecture.md.

About and the clean-room rule

lindhard is a Monte Carlo engine for the transport of ions and electrons in matter. Every model comes from published papers and cites its source where it is implemented. The code is clean-room: nothing is ported from code whose license does not allow it, and no tabulated output of other programs is used as data. What may be ported, what may only be read about in papers, and what may never enter the tree are set out as license tiers in CONTRIBUTING.md, together with the disclosure rules.

lindhard is released under the MIT license (LICENSE). The source is at https://github.com/2AMLogic/lindhard.

Installing

lindhard is pre-release and is not yet published as a package; build it from source. It is pure Rust with no C or Fortran dependencies in the default build.

Prerequisites

• A stable Rust toolchain. The repository’s rust-toolchain.toml selects the stable channel, so rustup installs the right one on first use.
• git.

Build

git clone https://github.com/2AMLogic/lindhard.git
cd lindhard
cargo build --release -p lindhard-cli

The binary is target/release/lindhard. Build in release mode for real runs; a debug build is many times slower. You can run it in place, through cargo, or install it into cargo’s binary directory:

target/release/lindhard --version
cargo run --release -p lindhard-cli -- --version
cargo install --path lindhard-cli      # puts `lindhard` on your PATH

--version prints the crate version and the git describe of the source it was built from, for example lindhard 0.0.1 (c2da637). The same two values are written into every summary.json, so a result can be traced to the binary that produced it.

Commands

lindhard check input.toml                  # parse and validate, no transport
lindhard run input.toml --out dir/         # run, write dir/summary.json and CSVs
lindhard run input.toml --out dir/ --ions 200 --seed 7 --threads 4
lindhard --version

--ions and --seed override run.ions and run.seed, and the output records the override. --threads overrides run.threads and never changes the results. The first run walks through both commands.

The library

The physics is in the lindhard library crate; the command line is a thin front end on it. To browse the library API:

cargo doc -p lindhard --open

First run: 5 keV B into Si

This walkthrough runs the example examples/b_5keV_si.toml: 10 000 boron ions of 5 keV into amorphous silicon, 7° off the surface normal. Every block of command output below is taken from an actual run of that file: the blocks are generated by validation/book_walkthrough.py, and CI reruns the example and fails if they no longer match. Run the commands from the root of the repository, with lindhard built as in Installing (or replace lindhard with cargo run --release -p lindhard-cli --).

1. The input

# 5 keV boron into amorphous silicon, 7 degrees off normal.
#
# Run:   lindhard run examples/b_5keV_si.toml --out out/b_5keV_si
# Check: lindhard check examples/b_5keV_si.toml
#
# Schema: docs/cli.md. Energies are in eV, lengths in nm, angles in degrees.

[beam]
ion = "B"
energy_ev = 5000.0
tilt_deg = 7.0

[target]
# An element symbol names the pure element at its tabulated density
# (lindhard/src/elements.rs).
substrate = "Si"

[physics]
potential = "zbl"
stopping = "lindhard-scharff"
free_path = "constant"
primary_cutoff_ev = 5.0
recoil_cutoff_ev = 2.0

# Si has no default displacement energy in the element table, so the input
# must choose one. 15 eV is an illustrative model parameter, not a
# recommendation: choose E_d for your problem.
[physics.energies.Si]
e_d_ev = 15.0

[run]
ions = 10000
seed = 1

[tally]
depth_bin_nm = 0.5
depth_bins = 200

# Fit a dual-Pearson profile to the depth histogram (off by default).
dual_pearson = true

Section by section:

• [beam]: the projectile by element symbol (its mass defaults to the standard atomic weight), its energy in eV and its tilt from the surface normal in degrees.
• [target]: substrate = "Si" names pure silicon at its tabulated density as a semi-infinite substrate. There are no finite layers, so nothing can be transmitted.
• [physics]: the model choices, spelled out even where they are the defaults: the ZBL universal potential, Lindhard-Scharff electronic stopping, and the constant free path. The two cutoffs are required: the beam ion is dropped below 5 eV, recoils below 2 eV (kept below the surface binding energy of Si, so that sputtering is not cut off).
• [physics.energies.Si]: the displacement energy \( E_d \) of Si. The element table has no default for Si, so the input must choose one; 15 eV here is an illustrative value, not a recommendation. Si’s surface binding energy \( E_s \) comes from the element table and its lattice binding \( E_b \) defaults to 0.
• [run]: the number of histories and the seed. Together with the input and the binary, the seed fixes every bit of the output.
• [tally]: a depth grid of 200 bins of 0.5 nm, and a dual-Pearson fit of the depth profile, which is off by default.

The full list of keys is in the TOML input reference.

2. Check the input

lindhard check examples/b_5keV_si.toml

examples/b_5keV_si.toml: OK
  beam: B at 5000 eV, tilt 7 deg, azimuth 0 deg; 10000 ions, seed 1
  layer 0: Si (Si 1.0000), semi-infinite, 4.9940e22 atoms/cm^3
  transport: amorphous-bca
  screening function: zbl-universal
  screening length: universal
  scattering angle: gauss-mehler-quadrature-table
  electronic stopping: lindhard-scharff
  compound stopping: bragg-additivity
  free path: constant
  displacement criterion: e_d-e_b
  surface barrier: planar
  random numbers: chacha8-per-history-stream

check parses and validates the input without transporting anything. It prints the resolved beam and target (the atom density comes from the tabulated mass density of Si) and every model the run will use. An invalid input exits non-zero with a message naming the offending key, for example target.layers[0].thickness_nm: -5 nm must be finite and positive.

3. Run it

lindhard run examples/b_5keV_si.toml --out out/b_5keV_si

10000 ions: 9440 stopped, 560 backscattered, 0 transmitted, 2549 sputtered atoms; wrote out/b_5keV_si

The one-line summary counts the fates of the 10 000 primaries and the target atoms that left through the front face. The run takes a few seconds in a release build. Add --threads N to choose the number of worker threads; the results are bit-identical whatever N is.

The output directory holds:

damage_profile.csv
depth_profile.csv
escape_spectra.csv
ions.csv
lateral_profile.csv
summary.json

4. Read the summary

summary.json starts with the format and software version, then the input exactly as it was run (defaults filled in, so the run can be reproduced from its own header), then the physics and the results. A few parts of it follow; the Output files page describes every key.

The models used. physics.models names every model with its citation (the citations are omitted here):

{
  "models": [
    {
      "role": "transport",
      "name": "amorphous-bca"
    },
    {
      "role": "screening function",
      "name": "zbl-universal"
    },
    {
      "role": "screening length",
      "name": "universal"
    },
    {
      "role": "scattering angle",
      "name": "gauss-mehler-quadrature-table"
    },
    {
      "role": "electronic stopping",
      "name": "lindhard-scharff"
    },
    {
      "role": "compound stopping",
      "name": "bragg-additivity"
    },
    {
      "role": "free path",
      "name": "constant"
    },
    {
      "role": "displacement criterion",
      "name": "e_d-e_b"
    },
    {
      "role": "surface barrier",
      "name": "planar"
    },
    {
      "role": "random numbers",
      "name": "chacha8-per-history-stream"
    }
  ]
}

Where the primaries went. results.primaries:

{
  "primaries": {
    "stopped": 9440,
    "backscattered": 560,
    "transmitted": 0,
    "stopped_depth_mean_nm": 24.39118364690772,
    "stopped_depth_std_nm": 12.729758786642684
  }
}

A few percent of the boron ions backscatter; the rest stop in the silicon.

The range distribution. results.range has the moments of the depth at which the stopped primaries came to rest, with their standard errors. mean_nm is the projected range \( R_p \) and std_dev_nm the straggle \( \Delta R_p \); the kurtosis is 3 for a Gaussian.

{
  "depth": {
    "n": 9440,
    "mean_nm": 24.39118364690772,
    "std_dev_nm": 12.729758786642684,
    "skewness": 0.37117262070685103,
    "kurtosis": 2.7137519074674534,
    "mean_std_err_nm": 0.13101896888723405,
    "std_dev_std_err_nm": 0.08575869371084056,
    "skewness_std_err": 0.02029360121408467,
    "kurtosis_std_err": 0.04691467531104645
  },
  "pearson_iv": null,
  "pearson_iv_error": "moments outside the Pearson IV region: skewness = 0.37117262070685103, kurtosis = 2.7137519074674534 (needs kurtosis > Some(3.260100391096362))"
}

When the moments lie outside the region of the Pearson IV family, as they do here, pearson_iv is null and pearson_iv_error gives the reason instead of a silent fallback. The dual-Pearson fit requested in [tally] (results.range.dual_pearson, shortened here) splits the profile into a near-surface head and a tail:

{
  "head_fraction": 0.19726158947438133,
  "head": {
    "mean_nm": 14.329373762517326,
    "std_dev_nm": 10.43582470528939
  },
  "tail": {
    "mean_nm": 26.250774127979774,
    "std_dev_nm": 12.442005256086365
  },
  "chi_square": 148.43612744959518
}

docs/validation.md (level 3) compares the projected range of this problem with a SIMS measurement of B in amorphous Si. With these default models the computed range runs long, and the page discusses why.

Damage. results.damage.per_ion puts the two kinds of damage figure side by side (Displacement damage): the NRT and Kinchin-Pease estimates from the damage energy of the primary knock-on atoms, and the vacancies, interstitials and replacements counted in the simulated cascades. They are different quantities and are not expected to agree.

{
  "per_ion": {
    "pka": 21.3829,
    "pka_energy_ev": 3375.856719281543,
    "damage_energy_ev": 2727.0751458675286,
    "nrt_displacements": 75.17911873954918,
    "kinchin_pease_displacements": 92.13933231739158,
    "displacements": 124.9056,
    "replacements": 8.3356,
    "vacancies": 116.57,
    "interstitials": 116.3151
  }
}

Sputtering and backscatter.

{
  "sputtering": {
    "yield_per_ion": 0.2549
  },
  "escapes": {
    "backscatter_coefficient": 0.056,
    "transmission_coefficient": 0.0,
    "energy_reflection_coefficient": 0.00983649666019863
  }
}

Where the energy went. Every history’s energy is accounted for, per ion: electronic loss, energy left in the lattice (sub-threshold transfers; \( E_b = 0 \) here), work against the surface barrier, energy carried out by backscattered ions and sputtered atoms, and the kinetic energy of particles when they fell below their cutoff (rest). The largest relative bookkeeping residual of any history is at the level of rounding.

{
  "energy_budget_ev_per_ion": {
    "incident": 5000.0,
    "electronic_nonlocal": 2172.1786901529567,
    "electronic_local": 0.0,
    "lattice": 2628.2107420032125,
    "surface_barrier": 1.1801870000000008,
    "backscattered": 49.18248330099315,
    "sputtered": 15.101149480441952,
    "transmitted": 0.0,
    "rest": 134.14674806239364,
    "max_relative_residual": 3.4560798667371273e-15
  }
}

5. The profiles

The CSV files hold the profiles on the grids set in [tally]. The first rows of depth_profile.csv, the stopped primaries per 0.5 nm bin and that count per incident ion per nm:

depth_lo_nm,depth_hi_nm,stopped_primaries,fraction_per_nm
0.0,0.5,13,0.0026
0.5,1.0,31,0.0062
1.0,1.5,38,0.0076
1.5,2.0,38,0.0076
2.0,2.5,50,0.01
2.5,3.0,45,0.009
3.0,3.5,68,0.0136

The last row has an upper edge of inf and collects everything deeper than the grid, so nothing is dropped. damage_profile.csv gives the cascade defects on the same grid, lateral_profile.csv the lateral and radial spread of the stopped primaries, escape_spectra.csv the energy and angle spectra of everything that left the target, and ions.csv the final state of every primary. Their columns are described in Output files.

6. Change something

Some edits to try, each a one-line change to the input:

• potential = "kr-c" and screening_length = "lindhard": another potential and screening length (Interatomic potentials).
• stopping = "equipartition-ls-or": half of the electronic loss taken locally at the collisions (Oen-Robinson); electronic_local in the energy budget becomes non-zero.
• weak_collisions = 3 under [physics]: weak collisions beyond \( p_\mathrm{max} \) (BCA transport).
• A finite layer in front of the substrate, as in examples/as_50keV_si_sio2.toml.

lindhard check shows the effect of each edit on the model list before you spend time on a run.

TOML input reference

An input file is TOML with these tables: [beam], [materials] (optional), [target], [physics], [stopping] (optional), [run] and [tally] (optional). An electron run has an [electron] table instead of [beam], [physics] and [tally] (section “Electron runs” below). The physics behind each [physics] choice is in the physics manual: potentials and screening lengths, electronic stopping, the free-path conventions and weak collisions, and user stopping tables. The first run walks through a complete file, and the repository’s examples/ directory has more.

This page is the input section of docs/cli.md, included here so there is one copy.

Units are in the key names: energies in eV (_ev), lengths in nm (_nm), angles in degrees (_deg), densities in g/cm³. Every table rejects unknown keys. The schema types are lindhard::input (shared with future front ends).

[beam]

[beam.divergence] (optional)

A finite angular spread of the ion beam about the nominal tilt_deg / azimuth_deg direction. Absent: every primary has the nominal direction and the run is bit-for-bit what it was before the table existed. Static ion runs only; with a [dynamic] target it is an error (beam.divergence), and an [electron] run has no such key (unknown key error).

Widths are finite and in [0, 10] degrees: the Gaussian law is a small-angle (plane-angle) reading, and the inward conditioning below is documented for small spreads only. The sampling laws are those of lindhard::ion::crystal::Divergence (docs/crystal-orientation.md).

How it is applied. Each primary’s direction is drawn once, at the entry point, from the nominal direction and the history’s own random stream, so it applies to the primary from the start, including while it crosses an amorphous screen layer before a crystal substrate; recoils keep their collision-generated directions. The crystal orientation is built from the nominal tilt and azimuth and does not move. The result is the inward- conditioned distribution: a drawn direction that does not point into the target (depth component not positive) is rejected and redrawn, at most 1000 attempts per primary, after which the run fails with an error instead of hanging (for the accepted widths one draw is inward with probability about one half or more, so this does not happen in practice). Near-grazing nominal tilts and wide spreads therefore give a distribution truncated at the surface, not the plain Gaussian or cone. The draws come from a copy of the history’s stream at word 2^65, disjoint from the transport draws, the thermal-vibration segment (2^66) and the lattice-shift segment (2^67), so turning divergence on changes no other random draw, and results do not depend on the thread count.

The summary’s physics.beam_divergence records the resolved model, width (width_rad and width_deg, with width_kind sigma_per_plane or cone_half_angle), the incidence policy inward-conditioned, the attempt bound and the stream segment; the echoed input.beam.divergence is the input as written. Both are absent without the table, and format.version stays 1 (see “Compatibility and extension”). In Python, Beam(..., divergence_model="gaussian" | "uniform-cone", divergence_deg=...) is the same table. Example: b_5keV_si_crystal_divergence.toml.

This is an input capability: it makes a finite-spread run possible and reproducible. It is not an experimentally validated channeling prediction, and it does not resolve the known deviations of the crystal model (#225).

[materials.<name>]

Named materials, in the form of lindhard::material::MaterialSpec:

[materials.SiO2]
density_g_cm3 = 2.2                    # required for compounds
elements = [
  { symbol = "Si", atom_fraction = 1.0 },
  { symbol = "O", atom_fraction = 2.0, e_d_ev = 20.0, e_s_ev = 2.0 },
]

Each element takes symbol or z, atom_fraction or mass_fraction (relative weights, normalised), and optional e_d_ev, e_b_ev, e_s_ev. Defaults for the energies come from the element table where one is tabulated; see lindhard/src/material.rs.

[target]

[target]
substrate = "Si"            # optional semi-infinite substrate

[[target.layers]]           # finite layers, front to back
material = "SiO2"
thickness_nm = 10.0

A material (in a layer or as the substrate) is either a key of [materials], an element symbol (the pure element at its tabulated density), or an inline material table with the same keys as [materials.<name>]. A [materials] key wins over an element symbol of the same name. Without a substrate the target has a back face and particles can be transmitted.

[physics]

stopping = "none" switches electronic stopping off: particles lose energy only in nuclear collisions, and the energy budget’s electronic terms are zero. It is not a realistic setting for any target. It exists for like-for-like transport comparisons against codes run with their electronic stopping off (the *_nuclear level-2 problems, validation.md) and for nuclear-only studies. It has no energy range, so no validity warning is ever raised for it; physics.models lists it as none, with nothing for Bragg additivity to sum unless [stopping] tables are given.

[physics.energies.<symbol>] sets e_d_ev, e_b_ev and/or e_s_ev for that element in every layer that contains it, after (so overriding) the material’s own values. An element with no tabulated default and no value set is an error naming the layer and the key to set.

tuning (phenomenological calibration, not a published model choice). "none" or omission leaves the physics and the echoed input exactly as without the key. A named set scales the surface binding energy E_s by a per-element factor fitted to measured data. Rules: the factor multiplies the resolved E_s (an explicit [physics.energies] or material value if given, else the elemental default), once per layer, after overrides; the global element table and the collision algorithm are untouched. The pilot supports static ion runs on single-element layers only, with a beam species the set was fitted for and target elements the set lists; compounds, [dynamic] targets, other beams, unlisted elements and unknown set names are rejected with a physics.tuning error. A beam energy outside the set’s fitted range, or a tilted beam, runs with a warning (an extrapolation). summary.json then has physics.tuning with the set, its version and provenance, and per layer the original E_s, the factor and the effective E_s. Tuned results must be reported next to, never in place of, untuned ones.

Shipped sets (fit record and held-out scores: docs/data-provenance.md, “Tuning factor sets”; a new version ships under a new name):

The factors are a calibration of yields under those settings; other settings (potential, E_d, cutoffs, weak collisions) were not part of the fit, and the set does not claim better accuracy for them.

[stopping]

Optional. Supplies user stopping tables for the electronic stopping of particular (ion, target element) pairs.

[stopping]
tables = ["tables/b_in_si.toml", "tables/p_in_si.toml"]

Absent, the input means what it always meant (format.version is unchanged and nothing is echoed).

Table file. The lindhard::ion::stopping::table::StoppingTable format:

provenance = "Author, Journal vol, page (year), Table N"   # required
ion_z = 5
ion_mass_amu = 11.0093            # optional; default: standard atomic weight
target_z = 14
energy_ev = [1.0e2, 1.0e3, 1.0e4]                 # strictly increasing, eV
stopping_ev_1e15_cm2 = [10.0, 30.0, 60.0]         # eV 1e-15 cm^2 per atom

Interpolation is piecewise linear in ln S versus ln E. provenance is mandatory (data without an origin is not admitted). Tables from any source, SRIM and ICRU included, are admissible; the provenance names the source (book edition, table and page, or program, version and file).

Paths. A relative path resolves against the directory of the input file (not the current directory). The echoed input keeps the path as written.

Composition with [physics] stopping. A table replaces the [physics] stopping model for exactly the pair it declares (ion_z, target_z), including recoils of that species when follow_recoils is on. Every other pair uses the [physics] stopping model. A pair with a table is never silently served by the model: a query outside the table’s energy range, or for a different ion mass, is an error that stops the run. Nothing is extrapolated. Declare each pair once. A table for a pair that cannot occur in the run (including a table for a target element’s recoils when follow_recoils = false) is accepted with a warning that it is unused.

Recoil species. Tables are keyed by (ion_z, target_z) only, so with follow_recoils = true a table whose ion_z is a target element also serves every recoil of that element. Recoils carry the standard atomic weight and are followed down to physics.recoil_cutoff_ev, so such a table is checked up front against both: its ion_mass_amu must be the standard weight, and it must start at or below recoil_cutoff_ev; either failure is an error. A consequence is that an isotopic self-ion beam (e.g. beam.mass_amu = 27.9769 for Si into Si) cannot take a table for its own pair while recoils are followed: drop the table, use the standard weight, or set follow_recoils = false. A recoil-species table that ends below the largest energy the beam can transfer to that element warns. Tables cannot be combined with stopping = "equipartition-ls-or" (that mode carries its own Lindhard-Scharff/Oen-Robinson loss and would ignore them). They can be combined with stopping = "none": the declared pairs use their tables and every other pair has no electronic stopping. That is coherent but easy to misread, so check physics.models (it lists both none and each user-table).

Errors name the field (stopping.tables[0]): an unknown key in [stopping], a missing or unreadable file, invalid table contents (including a missing provenance), a duplicate pair, a table whose ion mass differs from the beam ion’s, a table whose range does not contain the beam energy, and the recoil-species checks above. A table that starts above physics.primary_cutoff_ev warns, because the run fails if a projectile slows below it.

Provenance in the output. Tables are user data, so the run records them. summary.json has physics.stopping_tables, one entry per table: path (as written), resolved_path (absolute where possible), sha256 of the file’s bytes, the table’s provenance string, ion_z, ion_mass_amu, target_z and the energy range. physics.models lists each as user-table with the path and provenance as its source. The key is absent without [stopping].

[[crystal]] (optional)

Optional. Runs the named stack layers on an explicit cubic lattice (the engine’s crystal flight model, Bca::with_crystal) instead of the amorphous random-target model. Every layer not named stays amorphous; absent, the input means what it always meant and nothing is echoed. One entry per crystal; several entries give different lattices or orientations in different layers.

[beam]
ion = "B"
energy_ev = 5000.0
tilt_deg = 7.0        # the crystal's tilt
azimuth_deg = 22.0    # the crystal's twist, from `reference` towards n x r

[[crystal]]
layers = [0]            # stack layer indices; the substrate is the last one
preset = "Si"           # "Si", "Ge", "GaAs" or "3C-SiC"
normal = [0, 0, 1]      # (hkl) of the surface; the inward normal
reference = [0, 1, 0]   # [uvw] in the surface plane: the zero of the azimuth
wafer_rotation_deg = 0.0

[crystal.thermal]       # optional; absent: a static lattice
temperature_k = 300.0
# debye_temperature_k = 640.0

Angles. beam.tilt_deg and beam.azimuth_deg are the engine’s tilt and twist; there is no second beam direction. Tilt is the polar angle from the inward surface normal, azimuth is measured in the surface plane from reference towards n x reference; docs/crystal-orientation.md states the full convention and a worked example.

Checks (lindhard check and run, before transport). Errors name the key (crystal[i].normal, .reference, .layers, .preset, .thermal.temperature_k, …): a zero or non-integer index triple, a reference not in the surface plane, a layer that does not exist or is already assigned, a preset whose elements are not in the layer’s material or whose atom density differs from the lattice’s by more than 5 %, and a non-positive search length. Unsupported combinations are errors, never a silent amorphous fallback: a [dynamic] target, an [electron] run, physics.free_path other than "constant", physics.weak_collisions > 0, and physics.tuning (the only set was fitted for amorphous runs). physics.stopping = "equipartition-ls-or" is allowed with a warning: its local Oen-Robinson loss has constants that are not verified against the paper.

Output. summary.json gains physics.crystal, one object per crystal in input order, exactly as the engine reports it (Bca::crystal_metadata): regions, the lattice constant and its temperature, normal_hkl, reference_uvw, tilt_rad, twist_rad, wafer_rotation_rad, the resolved search parameters p_max_m, q_max_m, search_length_m (metres), thermal (its input, the per-species RMS displacement, the search margin and the sampling rule) and, when it applies, electronic_constants_unverified: true. physics.models gains a crystal transport entry. An amorphous input gets neither. format.version stays 1: the key is an addition (see “Compatibility and extension”).

Scope and caveats. Static, layered ion targets with the cubic presets only; hexagonal lattices, custom cells and dose-dependent crystal damage are not exposed (beam divergence is: see [beam.divergence]). This is an interface to the existing engine, not a validation: the known deviations of its channeled ranges are tracked separately (issue #225) and nothing here claims they are resolved. Examples: b_5keV_si_crystal.toml and as_50keV_sio2_on_si_crystal.toml (an amorphous oxide over a crystal substrate); neither is a validated prediction.

[dynamic] (optional)

Makes the run fluence-dependent: the target composition is updated as the fluence builds up (sputter erosion, build-up of implanted atoms). Without the table nothing changes: the run, its output files and its bytes are those of a static run. The model, its conventions and its limits are in lindhard::ion::dynamic and the book chapter on dynamic composition.

run.ions is the number of histories of the whole run and fluence_cm2 the fluence they represent, so each ion stands for fluence_cm2 / run.ions ions/cm². The ions are delivered in steps; after each step the grid is updated from that step’s events (a recoil is subtracted where it is created and added where it stops, a stopped beam ion is added, an escaped atom is a loss; the substrate is an immutable reservoir) and relaxed. Primary i of the run always uses the random stream (seed, i), whatever the step sizes and thread count.

With erosion = false the front surface stays at x = 0: swelling moves the interior interfaces and the back face of the slabs, not the front surface (the surface_nm column is that fixed frame, always 0). With erosion = true the lost thickness is removed from the front and the grid is re-anchored so the current surface is again x = 0; surface_nm is then the cumulative recession R in nm and a depth x in the output is x + R in the original frame. The recession of a step is the volume of the removed atoms per area under the chosen relaxation (sum Z removed_Z v_Z, or sum removed / n), and removal equals the sputtered counts per element, so no atom is created or lost. If a step sputters more of an element than the slabs hold, the excess is not removed and the element is counted in the clamped column. With a substrate, atoms sputtered from the substrate remove nothing from the slabs, so the recession falls short of Y F / n once the film is thin. With a retention_limit the excess over the limit is removed after each step (model and assumptions: the book’s dynamic-composition page); the steps file gains cum_re_emitted_atoms_per_cm2_<Sym> columns, the atoms/cm² re-emitted so far per element. Without one the output is unchanged. dynamic_summary.json totals gain recession_nm only with erosion on. Depths in the output are measured from the front surface of that step. A dynamic run needs E_d for every element that can occur, including the beam species. The Python bindings run static inputs only.

[run]

[tally]

[tally.voxel] takes origin_nm (corner of the first voxel, [x, y, z]), spacing_nm (voxel edge lengths) and dims (voxel counts), all three-element arrays; x is the depth and the beam starts at y = z = 0. It is accepted for static ion runs only. Events outside the grid are summed in the file header, never dropped.

The depth grid (depth_bin_nm, depth_bins, from the front face) is shared by the range histogram, the dual-Pearson fit and the defect profiles. Particles outside any grid are counted in explicit underflow and overflow entries, never dropped. Every count and per-ion value is for the same incident ions.

Electron runs ([electron])

An input with an [electron] table runs the low-energy electron engine (lindhard::electron::transport) instead of the ion BCA, with the full electron tally (lindhard::tally::FullElectronTally). It exposes what the library does and adds no physics; every choice and every data provenance is written to the output, so a result can be reproduced from its own header. The schema types are lindhard::input::electron. [materials] and [target] are the ion run’s tables (above); [beam], [physics], [stopping], [tally] and [dynamic] are not accepted. Example: ../examples/electron/e_10keV_si.toml.

[electron.beam]
energy_ev = 10000.0

[electron.transport]
cutoff_ev = 1.0
cutoff_reference = "vacuum-level"
secondaries = "kieft-bosch"
boundary = "step-barrier"

[electron.elastic]
potential = "thomas-fermi-yukawa"

[electron.inelastic]
model = "penn-single-pole"

[electron.materials.Si]
optical_elf = "si_elf.toml"
band = { kind = "insulator", valence_band_width_ev = 10.0, band_gap_ev = 2.0, affinity_ev = 3.0, provenance = "..." }

[target]
substrate = "Si"

[run]
histories = 1000
seed = 1

[electron.beam]

Primaries start on the front face (just outside it with the step barrier) at the landing point, which defaults to y = z = 0. [electron.beam.profile] selects the lateral profile by model: point (the default: every primary lands at the landing point), gaussian (sigma_nm, the standard deviation of each of y and z), disc (radius_nm, uniform in area) or rectangle (width_y_nm, width_z_nm, full widths, uniform). Sizes must be positive. Each history draws its start position, then its direction, from its own random stream, so results do not depend on the thread count. Divergence is conditioned on directions into the target, as for ions. The default (a point beam, no divergence) draws no extra random numbers, so existing inputs give identical results for a fixed seed. Through the library, use Transport::run_beam with a BeamSpread; Transport::run is the point beam. The slab has no lateral edge, so a spot is never partly off the face here.

[electron.transport] (lindhard::electron::transport::TransportConfig)

[electron.elastic]

The corrections are solved per grid energy with the chosen potential’s own Poisson density: the stand-in’s, or the DHFS density of Salvat et al. (1987) Eq. (12) (AtomicElastic::compute_corrected); the elastic table’s model and provenance strings name them and every polarizability with its source.

[electron.inelastic]

Energy axis of the inelastic table. The transport reads a table at the electron’s kinetic energy E above the band bottom of its layer. For every material with a band the inelastic table is built on that axis, with the model’s Fermi energy taken from the band: the row at E is the model’s at E - E_F, so the kinematics use E and no loss exceeds E - E_F. E_F is the band’s minimum excitation energy: the Fermi energy of a metal (fermi_ev, or the free-electron value), the conduction-band bottom valence_band_width_ev + band_gap_ev of an insulator. Table energies at or below it have no losses. This is the convention of the table compiler of Nebula (compile_full_imfp_icdf, cstool commit 0c739eb; see electron::inelastic::table, “Energy axis”, and electron::transport, “Energy reference of the inelastic table”). fermi_energy_ev plays no part for such a material, which is why it must be 0: adding it would count the Fermi energy twice. A material without a band gets the table of the model’s own axis with fermi_energy_ev, read at E as it is.

Migration. Before this rule, fermi_energy_ev applied to every material and tables ignored the band; losses beyond E - E_F were clamped by the secondary model. An input with a band and a nonzero fermi_energy_ev now fails, naming the key: remove the key. Results of runs with a band change (see CHANGELOG.md).

[electron.tables]: one log-spaced energy grid shared by the elastic and inelastic tables of every material. The inelastic table of a material with a band holds extra rows besides (#339; “Grid for a material with a band”, below).

The transport holds the rates of the first and last rows beyond the grid; a grid that starts above the lowest stopping threshold or ends below the largest possible energy warns.

Grid for a material with a band (#291). On the band-bottom axis the inelastic rate is zero at E_F and climbs steeply above it (within the first few eV, and for the Al single pole again at 15 to 17 eV), and the transport interpolates it linearly between grid rows. None of 10, 20 or 40 points per decade resolves that region to 1 %; at 80 the Mermin tables of Al and Si pass and the single poles still miss it. The largest relative error of the interpolated 1/λ over E - E_F = 2 to 20 eV (Al and Si, 400 eV inputs) is 4 to 75 % at the default 20 points per decade and 2 to 27 % at 40, depending on material and model; the single-pole model is worse than Mermin. Full table in validation.md, “Table grid above the Fermi level”; reproduce with lindhard-cli/examples/inelastic_low_energy.rs, mode rate. The library can add rows to the inelastic table from E_F upward (InelasticTableOptions::rate_refinement, #339), which brings every one of those tables to about 1 % at both the default grid and the grid of validation/experiments/se_yield.py. lindhard run uses it for the inelastic table of every material with a band (#339): every cell of the grid overlapping E - E_F = 2 to 20 eV is bisected while the linear interpolation of the rate misses the exact rate by more than 1 % at the cell’s quartiles, midpoint or estimated error peak. The table then has the grid’s rows and the added ones (a few tens: Al and Si, 5 to 25 rows at these grids); physics.materials.inelastic_table.energies counts them. The elastic table, and the inelastic table of a material without a band, keep exactly the grid’s rows. A denser grid is not needed for this region.

[electron.materials.<name>]: the electron data of each material the target uses, keyed by the name the target gives it (a [materials] key or an element symbol; an inline target material is an error in an electron run). Every name the target uses needs an entry; an unused entry warns.

[electron.materials.<name>.mermin_fit] (lindhard::electron::inelastic::MerminFitOptions): the mermin-melf model fits Drude-Lorentz oscillators to each material’s optical ELF (fit_mermin_oscillators) and extends them in momentum with the Mermin dielectric function (Mermin, Phys. Rev. B 1, 2362 (1970); de Vera et al., Int. J. Mol. Sci. 23, 6121 (2022)). This optional table sets the fit for one material; without it the fit uses MerminFitOptions::default() and results are unchanged. With any other electron.inelastic.model the table is an error, since nothing would use it.

The starting point and the iteration cap of the fit are not exposed. The fit is recorded, oscillator by oscillator, in the inelastic table’s model string; the fit also reports its f-sum and P_eff sum rules against the table’s (MerminFit), which are worth checking when choosing the options. Three oscillators do not describe every ELF: the default fit of the Cu ELF in validation/data/optical recovers 57 % of its f-sum and 67 % of its P_eff (docs/validation.md, “Cu: the Mermin IMFP and the default oscillator fit”), so such a material should use more. More oscillators can also give a very wide one, far wider than the table’s energy range, whose f-sum lies above the table (#311). Measured on that Cu ELF (same section, “Weighting options and the wide oscillator”): at 10 oscillators the default weighting gives one (f-sum 188 times the table’s, unconverged), relative_floor = 0.1 and 1 do not (f-sum 0.79 and 0.87 of the table’s, P_eff 0.98 and 1.03, converged), and uniform is worse (4e4 times). At 16 oscillators none of relative_floor 1e-3, 1e-2, 0.1, 1 or uniform avoids it (31 to 3e33 times the table’s f-sum), including fits that report convergence. So, whatever the options, check the fit’s f-sum and P_eff against the table’s in the fit report (MerminFit); convergence alone and a plausible IMFP do not show a sound fit. For example:

[electron.inelastic]
model = "mermin-melf"

[electron.materials.Cu]
optical_elf = "cu_elf.toml"

[electron.materials.Cu.mermin_fit]
oscillators = 6

The echoed input repeats the table, with its defaults filled in, only when it is given.

No optical or band data of any real material is committed (data-provenance.md); the data files are the user’s, and their terms are the user’s concern. Subshell binding-energy tables (SubshellBindingTable) are read through the optional binding key of inner_shells, below.

[electron.materials.<name>.inner_shells] adds inner-shell ionisation channels to a material (lindhard::electron::transport::Transport::with_inner_shells, attached to every layer of that material):

[electron.inelastic]
model = "penn-single-pole"

[electron.materials.Si.inner_shells]
subshells = ["K", "L1", "L2", "L3"]
# binding = "my_binding.toml"   # optional override of the EADL2017 table

Each channel is the shell’s hydrogenic optical ELF (electron::inelastic::hydrogenic_shell_elfs, atom density from the material) tabulated with single-pole Penn on the energy axis of the material’s inelastic table (build_shell_channel_tables). Limits, as in the library: single-pole Penn per channel, no exchange correction, and the material’s valence ELF must not already contain the shells (no partition rule is applied: an ELF that covers the shell edges counts them twice). A material with no inner_shells table runs exactly as before. The summary records the subshells and the table under physics.materials[].inner_shells (binding_source is "file", with the file’s path and sha256, or "embedded EADL2017", plus binding_provenance), the channel in physics.models, and the channels per layer in physics.transport. The channel tables are built on every run and are not part of the table cache.

[electron.tally] (lindhard::tally::ElectronTallyConfig)

[electron.tally.psf] (lindhard::tally::PsfConfig, fit_psf): the energy deposited in the slab depth_lo_nm <= x < depth_hi_nm, binned in the distance from the beam axis on a central disc plus log-spaced bins, with per-history errors, and double- and triple-Gaussian fits to it.

[electron.tally.dose] (lindhard::tally::DoseConfig, module lindhard::tally::dose): the energy and the dose in each layer per primary, a depth-dose profile through the stack and the electrons crossing each interface, all with per-history standard errors. Absent, none of it is tallied and nothing else changes.

Conventions. The target is laterally infinite, so the energy deposited in a layer per primary divided by the layer’s areal mass density ρ t is the absorbed dose per unit incident fluence, in Gy m² (multiply by a fluence in m⁻²). A semi-infinite last layer has no finite areal mass and is reported as energy only. Deposits use the tally’s point convention (an inelastic loss at the collision point, the sub-threshold remainder where the electron stops) and are measured from the band bottom, so with a secondary model they include the Fermi-sea energy of liberated conduction electrons and can sum to more than the energy the beam imparted. A deposit is binned by the layer it was made in, so the bins of a layer sum to its layer total. A crossing is counted each time an electron crosses an interface, at its kinetic energy on the new side. The standard error of a quantity is that of the mean of the per-history totals (the primary and its secondaries together), NaN for fewer than two histories.

[electron.tally.voxel]: the same fields as the ion [tally.voxel] (origin_nm, spacing_nm, dims; same validation). Tallies the deposited energy per voxel, by channel (electronic_ev, lattice_ev, rest_ev), and writes it as damage_export.csv with engine = electron (see Damage export format). The voxel energies plus outside sum to results.budget.deposited_ev. The species, defect and displacement columns are empty or zero: the electron tally has no species axis. Deposits are point events, so a flight is not split across voxel faces. Without the table no file is written, and one left in a reused output directory is removed.

[electron.tally.exposure] (lindhard::tally::exposure): the energy X that one primary history (the primary and its secondaries) deposits in each target volume, binned into a histogram with an exact bin at X = 0, and the distribution of S = X_1 + ... + X_N for N ~ Poisson(lambda) electrons. E[S] = lambda E[X] and Var[S] = lambda E[X^2] use the exact moment sums; quantiles and threshold probabilities come from a Poisson-weighted n-fold convolution of the binned X and are reported as low / central / high: the three values with every bin at its lower edge, centre and upper edge (the true quantile and P(S > b) lie in [low, high], the true P(S < a) in [high, low]; histories above x_max_ev are clamped into the last bin and counted in overflow_fraction). Valid only for a laterally uniform target: the landing offset of the beam is drawn per history and the volume shifted by it. A volume extending past the film collects only the part inside it. A lambda too large for the convolution gives the exact moments and a note instead of the distribution.

[run] of an electron run: histories (required), seed (required) and threads (all cores; not echoed, never changes results).

Output files

lindhard run input.toml --out dir/ writes summary.json and up to five CSV profiles into dir/. The first run shows excerpts of each from a real run. The models behind the damage figures are on Displacement damage.

This page is the output section of docs/cli.md, included here so there is one copy.

summary.json

software.git_describe freshness. The CLI’s build script stamps git describe --always --dirty, and reruns it when a build input that the description reads changes: tracked sources and manifests of lindhard and lindhard-cli (plus the workspace Cargo.toml, Cargo.lock and rust-toolchain.toml), the index, HEAD, the branch ref, packed-refs and refs/heads/refs/tags of the common git directory (so linked worktrees follow shared tags and branches). A non-empty LINDHARD_GIT_DESCRIBE in the build environment is used verbatim (trimmed) and git is not consulted; outside a git checkout the value is unknown. Limitations: -dirty is a flag, not a fingerprint, so two different dirty builds report the same string; untracked files never affect the -dirty flag, while git describe --dirty covers the whole tracked checkout, so an edit to any tracked file (docs, tests, other crates included) marks the stamp -dirty whenever the script next reruns, but a tracked file outside the watched set does not by itself trigger a rerun, so the stamp can stay stale until a watched input or the index changes (for example after git add); and it describes the checkout at build time, not the commit of a binary built elsewhere.

depth_profile.csv

depth_lo_nm,depth_hi_nm,stopped_primaries,fraction_per_nm: stopped primaries per depth bin, and that count per incident ion per nm. The last row’s upper edge is inf (overflow) and its density is empty.

lateral_profile.csv

quantity,lo_nm,hi_nm,count,per_ion_per_nm: stopped primaries by y, z and radial distance (quantity is y, z or radial), for the grid set by lateral_bin_nm and lateral_bins. Each quantity ends with an underflow row (lo_nm = -inf) and an overflow row (hi_nm = inf) with empty density.

damage_profile.csv

depth_lo_nm,depth_hi_nm,vacancies,interstitials,replacements, then the same three per incident ion per nm: the cascade defect counts by depth, on the depth grid. vacancies is displacements minus replacements. The last row, with depth_hi_nm = inf, holds everything deeper than the grid (empty densities). The totals equal results.damage.cascade.

escape_spectra.csv

species_z,symbol,beam,face,spectrum,lo,hi,count,per_ion_per_unit: the energy (spectrum = energy_ev, lo/hi in eV) and polar-angle (polar_deg, degrees from the outward normal) spectra of every species leaving each face (front or back). Each spectrum ends with -inf and inf rows for entries outside the grid, with empty densities. The density is per incident ion per eV or per degree.

voxel_tally.csv and voxel_tally.vtk

Written only with a [tally.voxel] table (and removed from a reused output directory when the table is absent); summary.json lists them under files. They are the library’s VoxelReport exports (format version 2, see Voxel tallies): per occupied voxel, the beam atoms implanted, the aggregate displacements, the energy deposited by channel, and for every species index of the run (0 is the beam, then the target elements in order of first appearance) vacancies_<i>_<symbol>, interstitials_<i>_<symbol> and replacements_<i>_<symbol>. The header # species = 0:B,1:Si,2:O names the species. Vacancies are placed at the displacement site and interstitials where the displaced atom stopped; the definitions and the cutoff dependence of replacements are those of Displacement damage. Summed over all voxels, the header’s outside sums and all species they equal the results.damage.cascade totals exactly.

damage_export.csv

Written with a [tally.voxel] table, removed from a reused output directory without it, and listed under files in summary.json. The same voxel map as voxel_tally.csv under the versioned exchange header (engine, version, seed, SHA-256 of the input echo, units); see Damage export format.

ions.csv

index,fate,x_nm,y_nm,z_nm,energy_ev,dir_x,dir_y,dir_z,layer: the final state of every primary, by history index. fate is stopped, backscattered or transmitted; for escaped primaries the position is on the face and the energy and direction are outside the target. x is depth.

Dynamic runs ([dynamic])

A run with a [dynamic] table writes dynamic_summary.json, dynamic_steps.csv and dynamic_composition.csv instead of the static files.

dynamic_summary.json (format lindhard-dynamic-summary): the echoed input (with dynamic), physics-model list models, the species, totals (ions, steps, rejected attempts, fluence, slab count and thickness before and after, yields per ion), files and the timing run object.

dynamic_steps.csv: one row per accepted step; step 0 is the initial target. step, first_index (global index of the step’s first ion), ions, ions_done, fluence_cm2 (delivered so far), attempts (more than 1 if the adaptive bound rejected the step), max_change, clamped, removed_slabs (slabs that emptied), n_slabs, surface_nm (cumulative surface recession in nm; always 0 with erosion = false), thickness_nm (total of the finite slabs), then cumulative counts since the start: cum_backscattered, cum_transmitted, cum_stopped_in_target, cum_stopped_in_substrate, cum_sputtered, sputter_yield (atoms per ion so far), cum_recoils_transmitted and cum_sputtered_<Sym> per element.

dynamic_composition.csv: the slab profile after every step (step 0 is the initial target), one row per step and slab: step, slab (0 is the front), front_nm, back_nm, thickness_nm, then per element atoms_per_cm2_<Sym> and fraction_<Sym> (atom fraction). Slabs that emptied are gone from later steps.

Electron runs: electron_summary.json and electron_*.csv

An electron run writes these instead of the ion files.

electron_summary.json (format {"name": "lindhard-electron-summary", "version": 1}):

results.damage (electron runs; an ion run’s results.damage is the NRT block above) splits the deposited energy into ionisation (electronic excitation, the source of electron-hole pairs) and lattice heating, and tests the beam against the displacement threshold. beam_energy_ev is the beam energy of the input and residual_channel is "ionisation". Per entry of layers: layer, material, ionisation_ev and lattice_ev (summed over histories) and the same _per_primary; with an e_pair for the layer’s material also e_pair_ev, e_pair_source and pairs_per_primary = ionisation_ev_per_primary / e_pair_ev (the keys are absent without e_pair; the top-level pairs_per_primary, the sum over layers, is present only when every layer that has ionisation energy has an e_pair); and knock_on, one object per element: element, z, t_max_ev = 2E(E + 2mc²)/(Mc²) (relativistic, free atom at rest, head-on, M the atomic weight times the atomic mass unit), e_d_ev, can_displace (T_max >= E_d) and threshold_beam_ev (the beam energy at which T_max = E_d). An element without E_d shows "not set" in those three fields. E_d is the material’s value ([materials] e_d_ev, or the element default where the element table has one). It is a threshold test only: no displacement cross section or dpa. For a step-barrier run E is the vacuum energy of the beam.

Conventions of the split (lindhard::tally::electron, “The ionisation and lattice channels”). Ionisation is the energy left in the solid at inelastic events plus the residual energy of electrons that fell below the cutoff (budget.residual_ev, the cutoff-dependent term; placed here by convention). Lattice is LO-phonon emission plus polaron traps (budget.phonon_emitted_ev + budget.polaron_ev; both are zero without the insulator channels). Ionisation plus lattice equals deposited_ev in the budget, per layer (deposition.per_layer_ionisation_ev, per_layer_lattice_ev) and per cell, to floating-point rounding, and like all deposits includes Fermi-sea energy with secondaries. Elastic nuclear recoil energy is not tallied by the electron engine and is in neither channel. budget.ionisation_ev and budget.lattice_ev are the totals.

results.table_coverage is a numerical diagnostic: one entry per layer (layer), and for its elastic and inelastic table the grid bounds energy_min_ev and energy_max_ev and the counts below (E < energy_min_ev: the first row’s rate and distribution were used), within (both bounds included: interpolated, or a row read exactly) and above (E > energy_max_ev: the last row’s were used), over primaries and secondaries. They count rate evaluations, not collisions: the transport evaluates both tables of the electron’s layer once before every free flight, including flights cut short at a layer face, the flight after a face reflection and flights with zero total rate, so the totals exceed the number of elastic and inelastic events, and a layer’s elastic and inelastic totals are equal. They are not fractions of path length or of deposited energy either. Nonzero below or above counts say that part of the transport used the constant continuation of a table beyond its grid ([electron.tables]); they do not say how much that changed the result. Summaries written before the key existed lack it.

results.event_caps is another numerical diagnostic, of the collision cap ([electron.transport] max_events): secondary_tracks is the number of secondary electrons the cap cut off (one per capped track) and affected_histories the number of primary histories in which the primary or at least one secondary was cut off (each history counted once, however many of its electrons were capped). Both are counts, not energies; the energy the capped electrons still carried is budget.event_cap_ev. The primaries’ own caps stay in fates.event_capped, so fates.event_capped = 0 alone does not show that the secondary cascades ran to completion: check event_caps.affected_histories. A capped track is a truncated one, not a physical fate. Summaries written before the key existed read back with both counts zero.

electron_escape_spectra.csv: face,spectrum,class,lo,hi,count,per_primary_per_unit. For each face (front, back): the energy spectrum of all escaping electrons (spectrum = energy_ev, class = all, eV) and the polar-angle spectrum of each class (polar_deg, slow or fast, degrees from the outward normal). Each spectrum ends with -inf and inf rows for entries outside the grid, with empty densities. The density is per primary per eV or per degree.

electron_deposition_cylindrical.csv (with tally.cylindrical): ir,ix,r_lo_nm,r_hi_nm,depth_lo_nm,depth_hi_nm,energy_ev,ev_per_primary_per_nm3,ionisation_ev,lattice_ev, one row per cell. electron_deposition_cartesian.csv (with tally.cartesian): ix,iy,iz,x_lo_nm,x_hi_nm,y_lo_nm,y_hi_nm,z_lo_nm,z_hi_nm,energy_ev,ev_per_primary_per_nm3,ionisation_ev,lattice_ev. ionisation_ev + lattice_ev = energy_ev per cell (the two columns are new and last, in eV summed over histories). Energy deposited outside a grid is outside_ev in the summary (with outside_ionisation_ev and outside_lattice_ev). Like budget.deposited_ev, the cell energies are measured from the band bottom and, with secondaries, include Fermi-sea energy the beam did not supply.

Charge (results.charge, and with a grid electron_charge_cylindrical.csv / electron_charge_cartesian.csv). Electron number, not energy, counted as integers per incident electron:

• every electron that ends inside the target counts +1 (stopped) where it ends, whatever the fate except escape: below the stopping threshold, trapped without an interaction, absorbed at the back face, polaron-trapped, or cut off by the collision cap;
• every liberated electron counts -1 (created, the hole it leaves) at the collision, followed or not. One below the stopping threshold is not followed: it is dropped where it was made and counted +1 there, so the two cancel. The stopping threshold follows electron.transport.cutoff_ev and cutoff_reference, so the cutoff decides how much of the creation and stopping is attributed to followed secondaries and where;
• an electron leaving through a face carries its charge out (escaped_front, escaped_back, by generation).

primary is the incident electron, secondary any electron the secondary model created. A created_by_* count is attributed to the generation of the electron whose collision made the hole. The balance is exact: incident + created = stopped + escaped, results.charge.imbalance is 0 on every run, and the net electron number left in the target, stopped - created, equals incident - escaped. results.charge.per_layer has the counts and the squared per-history net charge for each layer; with a grid, results.charge.cylindrical and .cartesian keep the grid and the totals inside_stopped, inside_created, and cells.outside for positions outside it.

The CSV columns after the cell edges (same cells and order as the deposition files) are stopped_primary,stopped_secondary,created_by_primary,created_by_secondary,net,net_per_primary,net_per_primary_se. net_per_primary is net per incident electron (an excess of electrons is positive); its standard error uses the net charge each history left in the cell, so correlated stops in one cell are not treated as independent. Summing the cylindrical rows over ir gives the depth profile. Counts are integers, so the files are identical on any thread count. Summaries written before the key existed read back with an empty charge report.

electron_psf_profile.csv (with tally.psf): one row per radial bin (the central disc first), r_lo_m,r_hi_m,r_center_m,area_m2,energy_ev,std_err_ev,density_ev_per_m2,density_std_err_ev_per_m2, then a model_<model>_ev column (the fitted energy of the bin) for each fit that succeeded. electron_psf_parameters.csv: model,parameter,value,std_error, the fitted parameters (lengths in m, energies in eV), then chi2, dof and reduced_chi2 rows per fit, then converged (1 if the iterations stopped on a converged step, 0 if the iteration budget ran out) and iterations (the integer count used). These metadata rows, like the chi2, dof and reduced_chi2 rows, have an empty std_error cell and are not fitted parameters, so a consumer that reads every row as a parameter must skip them. Parameters from a fit with converged 0 need inspection of the fit diagnostics before use.

electron_depth_dose.csv (with tally.dose.depth_bins_per_layer): one row per depth bin, bin,layer,x_lo_nm,x_hi_nm,energy_ev,energy_ev_std_err,ev_per_nm,ev_per_nm_std_err,dose_gy_m2,dose_gy_m2_std_err, per primary; the dose cells are empty in a semi-infinite layer. electron_interface_transmission.csv (with tally.dose.transmission and at least two layers): one row per interface, direction, class (primary or secondary) and spectrum bin, interface,from_layer,to_layer,class,count_per_primary,count_std_err,energy_ev_per_primary,energy_std_err,energy_fraction,bin,e_lo_ev,e_hi_ev,crossings, the totals repeated on each row of a spectrum.

electron_tables.csv: material,energy_ev,elastic_inverse_mfp_per_nm,inelastic_inverse_mfp_per_nm,inelastic_mean_loss_ev,inelastic_stopping_ev_per_nm, the tables the run used, per material and inelastic-table energy (the stopping power is λ⁻¹ ⟨W⟩ of the stored loss distribution). The elastic column is empty at the rows the rate refinement adds to the inelastic table of a material with a band (#339), which the elastic table does not have.

Cross-section table cache (--table-cache)

Building the elastic and inelastic tables is the slow part of a short electron run (tens of seconds with penn-single-pole, far longer with penn-full; see docs/validation.md, “Electron oracles”). With --table-cache DIR, lindhard run looks each table up in DIR (created if missing) and builds and stores only the ones it does not find, so a series of runs that differ only in seed, history count, tallies or transport settings builds its tables once.

Each entry is three files, named by the SHA-256 of a key document: <kind>-<sha256>.toml, the table in the versioned cache form of lindhard::electron::data::CrossSectionTable, read back with that type’s loader; <kind>-<sha256>.sha256, the hex SHA-256 of the table file’s bytes; and <kind>-<sha256>.key.json, the key itself. The key spells out every input the table depends on:

• the key schema version and the table cache format_version;
• the build: crate version, git describe --always --dirty, and the SHA-256 of the running lindhard executable, so any rebuild that changes the code (an uncommitted edit included) misses, and a rebuilt binary never reuses an older binary’s table;
• the table kind and the exact energy grid (electron.tables, after the defaults are filled in);
• the material’s name, composition and density;
• elastic: the model (electron.elastic.model), the potential, the exchange and correlation-polarization corrections with all their inputs, the starting probability grid and the refinement tolerance;
• inelastic: the model (electron.inelastic.model), its Fermi energy, the SHA-256 and provenance of the optical ELF file, the material’s band parameters (they set the table’s axis and Fermi energy for a material with a band) and, with mermin-melf, the material’s resolved Mermin fit options (mermin_fit; no table and a table of the defaults give the same key), and the rate refinement of a material with a band (rate_refinement: its window, tolerance and row bound; null without a band).

Every f64 is written in shortest round-trip form, so a change in the last bit of any number is a different key. A lookup must find the stored key equal, byte for byte, to the run’s own (a mismatch under the same hash means the file was edited, and is an error naming the differing field). The table file’s bytes must then hash to the SHA-256 stored beside it: a table edited or corrupted anywhere, even in a single cross-section or probability value that still parses, is refused with an error naming the file, before it is parsed, and a missing .sha256 file is an error too. The table is then loaded and validated by the library, and its axis and energy grid are checked against the run: an elastic table must have exactly the run’s grid, an inelastic table must contain it in order (it may hold the refinement’s rows as well). A file found under the run’s key that fails any of these checks is an error, never a silent rebuild; remove the entry’s files to rebuild it. A table of another cache format_version can never be found, since the version is part of the key. Writes go to a temporary file renamed into place, the table first, then its hash, then its key, so a key on disk always has its table and hash beside it.

Concurrent runs. Several runs may share one cache directory at once (threads of one process or separate lindhard processes, such as a seed sweep started in parallel on a cold cache). Each entry has a lock file, <kind>-<sha256>.lock, and a run holds the entry’s lock while it looks the entry up, builds the table and stores it. So for one key:

• one run builds the table; every other run waits for it, then finds the stored entry (the lookup is made only once the lock is held) and reads it. A run that has to wait for another process prints note: waiting for another process to finish the table cache entry <lock file> on standard error. A run whose build fails stores nothing, and the next waiting run builds in turn;
• runs that need different tables (another material, kind or key) take different locks and do not wait for each other;
• there is no timeout: a run waits as long as the holder is alive, since a build can take minutes. A holder that exits or is killed releases its lock (the operating system drops it with the process), so it never blocks the next run. A holder that is stopped but alive (for example suspended with Ctrl-Z) does block it. A waiting run can be interrupted safely;
• a holder killed while storing can leave temporary files (.<file>.tmp.<pid>.<n>) and a table without its key. The key is the commit marker, so that entry was never stored: the next run to take the lock removes the temporary files and rebuilds it. A stored entry that fails any check above is still an error, never a rebuild.

Temporary files are created under names unique to the process and the write, and removed if the write fails. The lock files are empty and stay in the directory; removing them while no run uses the cache is harmless.

Read-only caches. A cache the run can read but not write (a read-only directory or file system, or one whose lock files belong to another user) is still used for what it holds. When an entry’s lock file cannot be created or opened for writing because permission is denied or the file system is read-only, the run looks the entry up without the lock. An entry that is there is read after every check above, the SHA-256 of the table’s bytes included; the key is the commit marker, so a key that can be read always has its complete table beside it. An entry that is not there is an error naming the directory (... has no entry for this run's <kind> table of <material> and is not writable ...): the run does not build a table it could not store. Make the directory writable, use another directory, or run without --table-cache. If only the directory is read-only and the entry’s lock file exists and is writable by the run, the lock is taken as usual; the outcome is the same (the run checks that it can create a file in the directory before it builds).

The lock is an operating-system advisory lock on the lock file (flock on Unix, LockFileEx on Windows, through Rust’s std::fs::File::lock), and inside one process the run also excludes its own threads without relying on it. This is supported on a local file system. On a network file system (NFS, SMB) whether locks are honoured between machines depends on the server and mount options, and is not guaranteed; share a cache between machines at your own risk, or give each machine its own. If the file system refuses the lock altogether the run fails with an error naming the lock file (use a local directory, or run without --table-cache). Older lindhard builds take no lock, so do not point one at a directory a current run is writing.

A cached table is the built one bit for bit (the TOML cache form round-trips every f64), so outputs do not depend on whether the tables were built or read, nor on the thread count. Only physics.materials.*_table.source and .cache, and the timings in run, differ. lindhard-cli/tests/examples.rs checks a building run, a storing run and reading runs on 1 and 8 threads against each other; lindhard-cli/tests/table_cache_concurrency.rs checks the concurrency contract above (threads and processes on one cold entry build once, different entries build at the same time, a killed holder does not block, a read-only cache serves its entries and refuses a miss). Remove the directory to reclaim space; nothing else prunes it.

Exposure: lindhard exposure

Turns a fitted point-spread function into the dose of a list of shots (a pattern exposure). It does not transport anything: the dose is the analytic sum of Gaussians described in the module documentation of lindhard::exposure, so it is fast and exact for point shots. --direct adds an independent transport run of the same shots (see “Direct mode” below).

lindhard exposure --shots shots.csv \
    --psf-parameters out/electron_psf_parameters.csv \
    --blur-m 1e-8 --x0-m -5e-7 --y0-m -5e-7 --dx-m 1e-8 --dy-m 1e-8 \
    --nx 100 --ny 100 --threshold-fraction 0.5 -o exposure/

PSF. Either a parameters file written by an electron run with tally.psf (--psf-parameters; --model double_gaussian|triple_gaussian picks a fit, default the triple fit if present), or typed in: --total-ev, --alpha-m, --beta-m, --eta and, for the triple form, --gamma-m and --nu. The dose is only as good as the fit; check its reduced chi-squared in electron_summary.json first. The film is assumed laterally uniform.

Shot list (--shots): CSV with a header; required columns x_m, y_m, weight (number of primaries, or a relative dwell); optional blur_m (extra beam width of this shot), shape (point, disc, rect), size_m (disc radius, or rectangle width) and size2_m (rectangle height, default size_m). Blank lines and # lines are skipped. Lengths are in metres.

Beam blur (--blur-m, b in exp(-r^2/b^2)) adds in quadrature to every width of the PSF (w^2 -> w^2 + b^2, exact for Gaussians; derivation in the module docs), together with each shot’s blur_m. Shaped spots are a fixed --subsamples x --subsamples midpoint sub-sampling of the spot (default 16); the error is stated in the module docs.

Grid: --x0-m --y0-m --dx-m --dy-m --nx --ny; the dose is evaluated at the grid points. --threads never changes the results (each grid point sums the shots in list order).

Contours. With --threshold (eV/m^2) or --threshold-fraction (of the largest grid dose), iso-dose contours are found by marching squares. A grid point is inside when its dose is at or above the threshold; each contour has the inside on its left, so outer boundaries have positive area and holes negative. A contour that leaves the grid is open and has no area or centroid (enlarge the grid). The area is the polygon area, so it differs from the smooth contour’s by O(h^2) in the grid spacing.

Files written to the output directory:

The large-area dose of a uniform exposure of w primaries per m^2 is w E (the PSF integrates to E), of which the forward term gives the share 1 / (1 + eta [+ nu]).

Direct mode: --direct

The convolution above assumes a laterally uniform film. --direct target.toml also evaluates the same shot list by electron transport, as an independent answer at the same grid points. The model, frame and error are in the module documentation of lindhard::exposure (“Direct Monte Carlo exposure”); in short:

• target.toml is an electron input ([electron], see “Electron runs”). Its layers, cross-section tables and beam energy are used; where the primaries land (and the spot) comes from the shot list, so the input’s own beam profile, landing point and divergence are not used. The target is that input’s laterally uniform layer stack: a laterally varying (voxel) target is available through the library and not yet through TOML.
• Each shot runs --direct-histories primaries (split over --direct-replicas independent replicas, at least 2) from the shot position with its spot shape; a point shot with blur b is a Gaussian beam exp(-r^2/b^2). A shaped spot with a nonzero blur is refused. Exposure coordinates are offsets from the target’s entry point, x and y being the geometry’s y and z.
• The energy deposited in the film slab --slab-lo-m (default 0) to --slab-hi-m (required) is tallied per grid cell and divided by the cell area; the point is the cell centre. The result is the mean over the replicas with its standard error (Monte Carlo noise only), reported next to the convolution averaged over the same cell.
• A point is flagged in exposure_comparison.csv when the doses differ by more than --tolerance-sigmas (default 5) combined standard errors, sqrt(se_direct^2 + se_conv^2), where se_conv is --conv-rel-err times the convolution dose (default 0: the PSF’s own uncertainty is not propagated). Points where both doses are below --floor-fraction (default 0.01) of the largest dose are not judged, as the run cannot resolve them.
• Seeds follow the run’s thread-count-independent scheme (one SplitMix64 mix per replica and shot), so the files are identical for any --threads.

The PSF flags (--psf-parameters or the typed-in values) are still needed for the comparison.

Reusing an output directory

--out may name an existing directory; the run overwrites the files it writes and creates the directory if needed. The CLI also owns the reserved optional file names of the run’s mode. After a successful run, an optional file the run did not produce is removed if present: ions.csv (without tally.per_ion), and electron_deposition_cartesian.csv or electron_deposition_cylindrical.csv (without the matching deposition grid), and electron_psf_profile.csv and electron_psf_parameters.csv (without tally.psf), electron_exposure_histogram.csv (without tally.exposure), damage_export.csv (without electron.tally.voxel), and electron_depth_dose.csv and electron_interface_transmission.csv (without the matching tally.dose part). A missing file is not an error; a failed removal is, and names the path. The summary is written last (also for dynamic runs, after both CSVs, so a failed CSV write leaves no new summary) and lists only files that exist.

Once output starts, the run first removes the previous summary of its own mode (summary.json, electron_summary.json or dynamic_summary.json) and only then replaces CSVs. The new summary is written to a freshly created temporary sibling (created exclusively, so an existing file or symlink of the same name is never touched; another name is used) and renamed into place after every other write succeeds. So, after a failed run, the absence of the mode’s summary means the output set is incomplete (it may hold a mix of old and new CSVs); a summary that is present describes a completed run. If the previous summary cannot be removed, the run fails naming its path before any CSV changes. A failure before output starts (simulation or serialization) leaves existing outputs as they were. This is not a transaction: the CSVs are replaced in place, and concurrent runs into one directory are not supported. Other files in the directory are never touched, and no cleanup happens between ion, electron and dynamic runs. Do not keep your own data under a reserved name.

Reproducibility

Everything in summary.json except the trailing run object, and every CSV file, is a function of the input and the binary only: byte-identical for the same input and seed at any thread count (for a dynamic run: the same three files, with run the only thread-dependent part of dynamic_summary.json). For an electron run the same holds for electron_summary.json (apart from run) and every electron_*.csv: the tables are built bit-identically on any thread count, and histories run in chunks of a fixed size (16, recorded as physics.transport.chunk_size) merged in chunk order. lindhard-cli/tests/examples.rs checks this on 1 and 4 threads (1, 2 and 8 for the dynamic example; 1 and 4 for the electron example, and with --table-cache 1 and 8). Floats are written in shortest round-trip form.

Compatibility and extension

Consumers must ignore keys they do not know. New results are added as new keys, never by changing existing ones:

• New tallies become new objects under results, with their settings as new keys under [tally] and their profiles as new CSV files listed under files. results.range, results.damage, results.sputtering and results.escapes were added this way, without a version bump.
• physics.crystal (the [[crystal]] metadata) was added this way, without a version bump; an amorphous run does not carry it. physics.beam_divergence (the [beam.divergence] metadata) likewise.
• New model choices become new values of the existing [physics] keys, or new keys with defaults, so existing inputs keep their meaning.
• format.version is bumped only when an existing key is removed or changes meaning.

The electron schema and its output follow the same rules, with lindhard-electron-summary counting its versions separately:

• A new library model becomes a new value of an existing key (electron.inelastic.model, electron.elastic.potential, electron.transport.secondaries, electron.transport.boundary, a band kind, a phonon preset) and a new physics.models entry; existing values keep their meaning and defaults never change.
• New per-material data (for example a precomputed cross-section cache) becomes a new optional key of [electron.materials.<name>] that names a file. It must be read with the lindhard::electron::data loader of its type, so data without a provenance is refused, and recorded under physics.materials with its path, SHA-256 and provenance.
• The cross-section table cache is the one exception to the rule above, by decision (#168): it is a command-line flag (--table-cache DIR), not an input key, because it never changes a result (as with --threads, the input and its echo stay the same whether or not tables are reused), and a content-addressed directory keyed on every input of the build cannot name a stale file the way a hand-written path can. Its tables are still read with the lindhard::electron::data loader and echoed under physics.materials with their path, SHA-256 and provenance.
• A new tally becomes a key under [electron.tally], an object under results and, for profiles, a new CSV file listed under files.
• Every table keeps deny_unknown_fields, and every default is echoed, so an input written today still means the same thing.

Range calculator in the browser

The range calculator runs the ion transport of this book in your browser, as WebAssembly on one thread. Give it an ion, an energy and an elemental target; it returns the projected range $R_p$, the straggle $\Delta R_p$ and a depth profile. It is a static page: there is no server, and nothing you enter leaves your machine.

A run of $10^4$ 5 keV B ions into Si takes a few seconds in a current browser. The calculator uses the choices of examples/b_5keV_si.toml (amorphous target at its tabulated density, ZBL screening, Lindhard–Scharff electronic stopping, 7° tilt), and a 15 eV displacement energy for a target element the element table gives none. It is not a replacement for the command line: it has no layers, compounds, crystals or electron mode.

Build it yourself

rustup target add wasm32-unknown-unknown
# wasm-bindgen command-line tool, release 0.2.129 (the version pinned in
# lindhard-wasm/Cargo.toml)
lindhard-wasm/build.sh web out/pkg
cp lindhard-wasm/web/index.html out/index.html
python3 -m http.server -d out        # then open http://localhost:8000

Agreement with the native build

The same input and seed give the same $R_p$ and $\Delta R_p$ as the native build, to within one standard error (in practice to about $10^{-15}$ relative): cargo test -p lindhard-wasm checks the native build and lindhard-wasm/tests/node/parity.test.mjs the WebAssembly build against one recorded result. The tolerance and its reason are in lindhard-wasm/tests/parity.rs: wasm32-unknown-unknown uses Rust’s bundled libm where a native build uses the platform’s, the two differ in the last bit, and a collision cascade amplifies that until individual histories diverge into statistically equivalent ones.

Simulator in the browser

The simulator is the layered-target companion of the range calculator. It runs the same engine as WebAssembly in a Web Worker, so the page stays responsive during a run, and it is static: there is no server and nothing you enter leaves your machine.

The page has a Setup view and a Run view. The Run view shows live counters, a trajectory view, the energy partition, live distribution plots and downloads.

Run

• Counters: ions done, R_p and ΔR_p, lateral (y, z) and radial spread, ions stopped, backscattered and transmitted per ion, sputtered atoms per ion, vacancies per ion (the followed cascade in full-cascade mode, the NRT estimate in quick mode), and the share of the incident energy that went to electrons, to the lattice, out of the target, and elsewhere. Brackets give one standard error.
• Trajectories: the paths of the first 50 ion histories, in a side view (depth across, lateral y down, layers shaded and named) or an end-on y × z view. Primary paths are always drawn; recoil cascades can be switched off. Species are told apart by colour and by the shape of the marker at the track’s end, and the generation by the dash pattern (solid primary, dashed first recoils, dotted later ones). A history cut at the per-track point limit is flagged.
• Distributions, updated as the run proceeds: the ion range distribution (total and, with more than one layer, per layer, with the Pearson IV density of the four depth moments drawn over it when it exists; Hofker 1975, Heinrich 2004), the lateral and radial distributions, damage versus depth (vacancies, interstitials, replacements; zero in quick mode), and escapes (counts per species and face, with energy and polar-angle spectra). Series are told apart by colour and by dash pattern, with a legend; counts that fall outside the plotted range are stated under each plot.
• Downloads: each plot’s data as CSV (one row per bin per series, with the under- and overflow counts as rows), the full report as JSON, and the trajectory view as a PNG. Results are shown with their statistical error, the models in use, and a pointer to docs/validation.md.
• Limit: the page runs up to 100 000 ions so that it stays responsive; use lindhard run for larger runs.

Setup

• Ion: element, mass (the standard atomic weight unless set), energy in eV, keV or MeV, and the tilt from the surface normal (7° by default).
• Target layers, surface first: a name, a thickness (Å, nm or µm) or “semi-infinite” for the last layer, a composition by atom fraction, a density and the displacement, lattice-binding and surface-binding energies of each element. A blank energy uses the element table’s value, shown in the field; a pure element takes its tabulated density; a compound needs a density. Layers can be added, deleted and reordered, and a preset fills one layer (Si, Cu, Au, SiO₂, GaAs). The preset densities are the ones recorded in docs/data-provenance.md; the O energies of the SiO₂ preset are the illustrative model parameters of examples/as_50keV_si_sio2.toml, not recommendations.
• Calculation: “ion distribution, quick damage” (follow_recoils = false) or “full cascades” (follow_recoils = true), the ion count and the seed.
• Physics: the potential and the electronic-stopping model, each with its citation, and the two energy cutoffs.

The input file

The page shows the TOML input its setup produces, with Download and Load buttons. The file runs unchanged under lindhard run, and gives the same report for the same seed and ion count (checked natively in lindhard-wasm/tests/cli_equivalence.rs). Loading a file fills the form. The form exposes a subset of the input keys; a key it cannot represent (a crystal, a dynamic target, a beam azimuth, mass_fraction, and so on) is reported by name and the file is not loaded, rather than silently dropped.

Limits

Amorphous, static targets only. Results carry statistical error (one standard error is shown) and the models listed on the page; the validation status is in docs/validation.md.

Voxel tallies and VTK export

The 1D tallies of Output files report results against depth. For a target that is not laterally uniform (a voxel grid, a triangle mesh or a CSG solid), the library’s lindhard::tally::VoxelTally bins the same events in 3D. It is a library tally: attach it to a run through Bca::run(|| VoxelTally::new(grid)) and read it with VoxelTally::report(). The command-line driver writes it when the input has a [tally.voxel] table (CLI reference): voxel_tally.csv and voxel_tally.vtk in the output directory.

What is tallied

The tally uses its own regular grid (VoxelTallyGrid: origin, voxel edge lengths and voxel counts), independent of the target geometry. For a voxel grid target, VoxelTallyGrid::from_voxel_grid gives the target’s own voxels. Each voxel owns the half-open box [o + i h, o + (i + 1) h) on every axis, and x is the depth.

• Per-element defects. A tally built with VoxelTally::with_species(&bca.species_z()) (the driver does this) also counts, for each species index of the run, the vacancies_<i>_<symbol> at the displacement site, the interstitials_<i>_<symbol> where the displaced atom came to rest without filling a site, and the replacements_<i>_<symbol>, with the definitions of Displacement damage and the same replacement rule as the 1D tally. A replacement is credited to the species of the atom it replaced, so vacancies + replacements = displacements per species and voxel. Beam atoms are implanted, not interstitials, and a displaced atom that leaves the target leaves a vacancy and no interstitial. Summed over voxels, outside and species the counts equal the run’s cascade totals exactly (integer equality, lindhard/tests/voxel_tally.rs). With follow_recoils = false there are no replacements. The aggregate displacements column is unchanged.
• Damage is the simulated cascade displacement count, the voxel version of the cascade counts in Displacement damage. The NRT / Kinchin-Pease model estimate is not tallied per voxel, because it needs the Lindhard partition with the mean Z and A of each material; the two are never mixed.
• Deposited energy is split by channel so that each can be checked against the engine’s energy budget. The work done against the surface barrier and the energy carried out of the target are not deposited in any voxel and are not tallied.
• The nonlocal electronic loss of a free flight is spread uniformly along the straight segment, as the engine defines it, so a segment crossing voxel faces is split in proportion to its length in each voxel. All other events go to the voxel that contains them. Events outside the grid are summed separately (outside in the report and the CSV header).

Sparse storage and determinism

Only voxels that received an event are stored (a BTreeMap on the flat voxel index i + nx (j + ny k)), so a grid of 10¹² voxels costs only its occupied voxels, and iteration, merging and export always follow the same order. Counts are integers and energies are f64 sums merged in chunk order, so for a fixed chunk size the report is bit-identical at any thread count. The conservation tests check that the totals over all voxels equal the engine’s own 1D tallies of the same run for a laterally uniform target (lindhard/tests/voxel_tally.rs).

CSV export

VoxelReport::write_csv writes a # key = value header that carries the grid, then one row per occupied voxel in flat-index order, with the voxel centre in metres and the energies in eV. Floats are written in the shortest form that reads back to the same bits, and VoxelReport::read_csv returns the report unchanged.

The format is version 1 (the banner # lindhard voxel tally v1) for a tally without species and version 2 with them: v2 adds the header lines # species = 0:B,1:Si,2:O and # outside_defects = ... (vacancies, interstitials, replacements per species) and the per-species columns after deposited_ev. The VTK title line carries species= and the same columns become cell scalars. Both readers accept both versions; v1 files and the checked-in examples below are unchanged. The file examples/voxel-tally/example.csv is a small synthetic report (made-up numbers, not a run):

# lindhard voxel tally v1
# origin_m = 0e0,0e0,0e0
# spacing_m = 5e-9,5e-9,5e-9
# dims = 4,2,2
# histories = 100
# outside = 0,1,3.5e0,0e0,0e0
i,j,k,x_m,y_m,z_m,implanted,displacements,electronic_ev,lattice_ev,rest_ev,deposited_ev
0,0,0,2.5e-9,2.5e-9,2.5e-9,0,12,1.8005e3,9.625e1,0e0,1.89675e3
1,0,0,7.500000000000001e-9,2.5e-9,2.5e-9,3,41,5.210125e3,3.105e2,4.75e0,5.525375e3
...

VTK export

VoxelReport::write_vtk writes a legacy ASCII VTK 3.0 file (Kitware, The VTK User’s Guide, “VTK File Formats”), which needs no library: an UNSTRUCTURED_GRID with one VTK_VOXEL cell (type 11) per occupied voxel, so the file is as sparse as the tally, and the columns above as cell scalars. ParaView and VisIt open it directly. The grid and the history count are in the title line. The file does not carry the outside sums; use the CSV for those. VoxelReport::read_vtk reads back what write_vtk writes, and the checked-in examples/voxel-tally/example.vtk is the same synthetic report as the CSV example. Both example files are compared with the exporters’ output in the test suite, and read back.

Layout cross-section import (finding)

Status: blocked upstream. No import code exists. This page records what the pinned layout tool offers, so the next step is a precise upstream request rather than a second layout reader in this tree.

The pin

Layout reading is not reimplemented in lindhard; it is meant to come from the klayout-tools package at one pinned version. The pin lives in one file, tools/klayout-tools.pin (currently 0.7.0, the latest PyPI release when this was written). Bump it deliberately, in its own change, and re-run the check below.

Run a pinned release without touching any installed tool:

uvx --from "klayout-tools==$(cat tools/klayout-tools.pin)" klt --help

What 0.7.0 offers

Checked by running each verb on a two-rectangle layout written with klt draw.

The Python modules of the package (klayout_tools._layout and friends) are underscore-private and not documented as an interface, and klayout.db is a transitive dependency, not part of the package’s contract. None of them is a supported polygon export.

Finding

The pinned release has no machine-readable path from a layout region to polygons. The only way to get geometry out of a region is klt clip, which writes a stream, and reading that stream back would be a second GDS/OASIS reader in this tree. That is ruled out, so the cross-section import (process stack, 2D/3D target build, area cross-check against the layout tool) is not started.

What is needed upstream

A read-side verb or documented Python function in klayout-tools that, for a layout, an optional --region and a layer list, returns per-layer merged polygons as JSON:

• units of micrometres and the database unit, as klt layers already reports;
• one entry per (layer, datatype) with a list of polygons, each a hull plus holes as vertex lists, hierarchy flattened and merged (same rule as merged_area_um2);
• a schema_version, so the seam can fail on an unsupported version.

With that, the import needs only a process-stack input format (layer to material, thickness, mask/etch rule), its validation, the target builders on the Phase 1 geometry kernels, and a test comparing imported areas with klt stats for a fixture layout.

Examples gallery

Four scenarios, each with a runnable input under examples/, the command that runs it, and a figure and table of the result:

1. Implant through a screen oxide: 50 keV As into Si through 10 nm of SiO₂.
2. Sputtering yield against incidence angle: 1 keV Ar on Cu, swept over the tilt.
3. Electron deposition in a resist on Si: 10 keV electrons into a 100 nm layer on Si, with synthetic electron data.
4. Dose in a buried oxide: electrons of 1 to 30 keV on a metal and dielectric stack over a thin oxide, with synthetic electron data.

Every figure and table is generated from a run of the lindhard binary by one committed script, validation/gallery.py; nothing is drawn by hand. From the root of the repository:

cargo build --release -p lindhard-cli
python3 validation/gallery.py --lindhard target/release/lindhard

The script takes about 20 seconds and rewrites the files under book/src/gallery/generated/. Each page states its seeds and particle counts. Results do not depend on the thread count; the trailing digits of a table can differ between platforms.

Every page says what is validated and what is only an illustration. The inputs are examples of how to use the engine. Their energies and the E_d values are illustrative model parameters, not recommendations, and no figure here is a validated prediction unless its page says which comparison supports it (the comparisons are in docs/validation.md).

Implant into Si through a screen oxide

50 keV arsenic into silicon at 7° off the surface normal, through a 10 nm SiO₂ screen layer. The input is examples/as_50keV_si_sio2.toml:

# 50 keV arsenic into silicon through a 10 nm SiO2 screen layer.
#
# Run: lindhard run examples/as_50keV_si_sio2.toml --out out/as_50keV_si_sio2
#
# Schema: docs/cli.md. Energies are in eV, lengths in nm, angles in degrees.

[beam]
ion = "As"
energy_ev = 50000.0
tilt_deg = 7.0

# A compound needs an explicit density. 2.2 g/cm^3 is the usual handbook
# density of vitreous (amorphous) silica. The O energies have no defaults in
# the element table; the values below are illustrative model parameters,
# not recommendations.
[materials.SiO2]
density_g_cm3 = 2.2
elements = [
  { symbol = "Si", atom_fraction = 1.0 },
  { symbol = "O", atom_fraction = 2.0, e_d_ev = 20.0, e_s_ev = 2.0 },
]

[target]
substrate = "Si"

[[target.layers]]
material = "SiO2"
thickness_nm = 10.0

[physics]
potential = "zbl"
stopping = "lindhard-scharff"
free_path = "constant"
primary_cutoff_ev = 5.0
recoil_cutoff_ev = 2.0

# Applies to Si in both the screen and the substrate (illustrative value).
[physics.energies.Si]
e_d_ev = 15.0

[run]
ions = 2000
seed = 1

[tally]
depth_bin_nm = 1.0
depth_bins = 150

It shows a compound material with an explicit density, a finite layer over a substrate, and per-element energy overrides ([materials.SiO2] and [physics.energies.Si]). Run and check it with

lindhard check examples/as_50keV_si_sio2.toml
lindhard run examples/as_50keV_si_sio2.toml --out out/as_50keV_si_sio2

(or replace lindhard with cargo run --release -p lindhard-cli --). The run takes a fraction of a second.

Result

2000 histories, seed 1, 1 nm depth bins. The figure is depth_profile.csv of that run: the stopped As per incident ion per nm, against depth from the front face of the oxide. The shaded band is the screen layer.

Nearly all ions pass the oxide and stop in the Si, so the profile is the bare-Si profile shifted by the oxide’s own stopping, with the interface visible as the edge of the shaded band. The statistical error of the moments is the standard error in the table; 2000 histories is a coarse profile, chosen to keep the page quick.

What is validated, and what is an illustration

• Validated: the projected range of 50 keV As in bare amorphous Si against the Lindhard-Scharff-Schiøtt first-moment equation, on the ZBL potential and Lindhard-Scharff stopping used here (docs/validation.md, level 1, check range.as50k_si.rp). That is an internal consistency check against the theory the model implements, not a comparison with a measurement.
• Illustration only: the oxide layer, the E_d and E_s values (illustrative; the O energies have no defaults in the element table), the 7° tilt and the sputtered-atom count. The amorphous target has no channeling tail; the crystalline variant (examples/as_50keV_sio2_on_si_crystal.toml) is likewise an input capability, not a validated channeling prediction. This page is not a validated prediction of a measured profile.

Sputtering yield against incidence angle

The yield of Cu under 1 keV Ar as the beam is tilted away from the surface normal. The base input is examples/sputtering_ar_cu_angle.toml, the physics of ar_1keV_cu.toml (KrC potential, Lindhard-Scharff stopping, E_s and E_d from the element table) with a short depth grid and no escape or lateral tallies:

# Sputtering yield of Cu under 1 keV Ar at one incidence angle: the base input
# of the angle sweep in the book's examples gallery.
#
# Run:   lindhard run examples/sputtering_ar_cu_angle.toml --out out/sputter_angle_0
# Sweep: validation/gallery.py runs this input at a series of `beam.tilt_deg`
#        values (the CLI has no flag for the tilt, so it writes one copy of
#        the input per angle) and tabulates results.sputtering.yield_per_ion.
#
# Same physics as ar_1keV_cu.toml (KrC potential, Lindhard-Scharff stopping,
# E_s and E_d from the element table) with a short depth grid and no escape or
# lateral tallies, so that a sweep over angles and seeds stays small. Illustration
# of how the BCA yield depends on the angle, not a validated prediction: see
# docs/validation.md for the sputtering comparisons that are validated.
#
# Schema: docs/cli.md. Energies are in eV, lengths in nm, angles in degrees.

[beam]
ion = "Ar"
energy_ev = 1000.0
tilt_deg = 0.0

[target]
substrate = "Cu"

[physics]
potential = "kr-c"
stopping = "lindhard-scharff"
free_path = "constant"
primary_cutoff_ev = 2.0
# Keep the recoil cutoff below E_s (3.49 eV for Cu), or sputtered atoms are
# lost.
recoil_cutoff_ev = 1.0

[run]
ions = 10000
seed = 1

[tally]
depth_bin_nm = 0.5
depth_bins = 40

The sweep

The CLI has no flag for the tilt, so the sweep needs one copy of the input per angle that differs only in beam.tilt_deg. The generator writes them to a scratch directory and runs each with --seed:

lindhard run examples/sputtering_ar_cu_angle.toml --out out/sputter_angle_0 --seed 1

is the run at 0°; validation/gallery.py does the same for 0°, 15°, 30°, 45°, 60°, 70° and 80°, and for seeds 1 to 4 at each (28 runs of 10 000 ions, a fraction of a second each), then reads results.sputtering.yield_per_ion from each summary.json. The mean over the four seeds is plotted; the error bar is the standard error of that mean from the spread of the four seeds, so it includes the correlation of atoms sputtered by one cascade.