Module table 

Inelastic sampling tables from the single-pole Penn DIIMFP (issue #94): an electron::data::CrossSectionTable with SamplingAxis::InelasticEnergyLoss (the inverse mean free path and the inverse CDF of the energy loss W at every energy of a caller-chosen grid), and a sampler of the momentum transfer q given (E, W), from which the deflection angle follows.

Both are built from SinglePolePenn only; no data and no equation beyond the model of super::penn enter here (the elementary kinematics of a binary collision aside, see below). build_inelastic_table_for_model runs the same energy-loss table procedure on any of the three models of super::model::PennInelastic (single pole, full Penn, Mermin-ELF), using that model’s own DIIMFP and inverse IMFP; nothing else changes.

Energy-loss table

At incident energy E the loss distribution is the DIIMFP at fixed E, p(E, W) (SinglePolePenn::diimfp_per_m_ev), normalised: F(W) = ∫_0^W p dW' / ∫_0^E p dW'. The table stores

• λ⁻¹(E) from SinglePolePenn::imfp_and_stopping (the closed-form reduction of #18, not the integral of p computed here; the two agree to the model’s tolerance, 1e-5 in the #18 tests), and
• the quantiles W(u_j) with F(W(u_j)) = u_j on a probability grid shared by all rows.

Density model. W spans decades, so the integrand is handled in s = ln W, where p dW = p W ds. g(s) = p W is sampled on adaptively chosen nodes between the lowest loss of the model’s rate (the lowest tabulated ELF energy, below which p = 0, for the single-pole model; see “Rows below the ELF table” for the Mermin and full Penn models) and E, and taken to be linear between nodes. Starting from INITIAL_LOSS_PANELS equal panels, a panel is split at its midpoint while the Simpson-minus-trapezoid estimate of its error in either ∫ g ds or ∫ W g ds exceeds tol · I · h/H (I the whole integral, h the panel width, H the whole width; tol = InelasticTableOptions::density_tolerance). The quantiles are then the exact inverse of this piecewise-quadratic CDF (a quadratic per panel, solved without cancellation), so the density model is the only approximation of the loss distribution and its error is controlled by that tolerance.

Full Penn rows. For the full Penn model the density is built on the same nodes from the DIIMFP of a DiimfpGrid tabulated once for the table’s energies (FullPenn::diimfp_grid, #256), not from the direct nested integrals: its interpolation error, held to the model tolerance (1e-4, relative or relative to the row’s mean density), is a second approximation of the loss distribution, next to the density tolerance. It is also smooth in W, where the direct DIIMFP carries quadrature noise at the model tolerance, which is above the default density tolerance and would drive the panel refinement to MAX_NODES.

Unresolved grid cells. Where the grid could not hold its interpolation to the model tolerance even at its narrowest cell (1e-4 in ln ω), it leaves the cell to the direct model (DiimfpGrid::unresolved_cells). The row integrand then switches model at each edge of such a cell, a step of the order of the model tolerance, and inside it carries the direct model’s quadrature noise. The split test above does not depend on the panel width for a step (both sides scale with h), so a step above about 3 · tol of the row’s mean density would be bisected down to floating-point resolution and then padded with duplicate nodes up to MAX_NODES, each a direct evaluation. So a panel within one grid-cell width (1e-4) of an unresolved cell is not split once it is no wider than that width: that is the finest scale at which the grid tested the DIIMFP, and below it the split test measures the noise and the step, not the shape. What such a panel can still miss is bounded by its width times the step, 1e-4 · (model tolerance) · g, per cell edge, far below the model tolerance of the row. Rows that do not come near an unresolved cell, and every model without a grid, are refined exactly as before. loss_density_node_count reports the node count of a row.

Floating-point resolution (#345). Every model is a numerical integral at its own tolerance, so its DIIMFP carries small jumps where the adaptive quadrature changes its subdivision between two adjacent W. The Mermin DIIMFP (super::mermin, integration tolerance 1e-5, the same as the default density tolerance) has them: the Si row at T = 1000 eV of the committed Yang et al. (2019) ELF (model Fermi energy 13.46 eV) has 18, of 0.4 to 7 times 1e-5 of the local density, between W = 20 and 976 eV. A jump is not a shape the density can resolve: as for the model switch at an unresolved cell, the split test of a panel across it does not fall with the width, so the panel is bisected until its midpoint rounds onto one of its ends. Such a panel (one ulp of ln W wide) is not split: its midpoint is not a new node, and splitting it would only add zero-width panels, up to MAX_NODES of them. This is this project’s choice of floor, with no tuned width: it is the point where bisection can no longer place a node. A jump bisected to it costs about 50 nodes (the levels from the starting panel to one ulp), and leaves nothing of the density unresolved: the panel holding it has no interior. Before #345 the zero-width panels were added until the row reached MAX_NODES; they carry no integral (h = 0), so a row whose other panels had converged is stored the same, with fewer evaluations (the Si row above: 20 099 nodes before, 5 792 now).

Node bound. A row whose refinement still has a failing panel when it reaches MAX_NODES nodes is refused, InelasticTableError::DensityNodeLimit, rather than stored with that panel unresolved (the convention of the rate refinement’s row bound, below).

Probability grid. The stored inverse CDF is interpolated linearly in u (CrossSectionTable::inverse_cdf). Its first moment, which is what sets the stopping power S = λ⁻¹ ⟨W⟩ with ⟨W⟩ = ∫_0^1 W(u) du, is dominated by the upper tail of W, so the starting grid (logistic_probability_grid) is uniform in t = ln(u/(1-u)), and it is refined adaptively, shared by all rows: a segment is split at its logistic midpoint when δ Δu > tol μ max(Δu, 1e-4) at any energy, with δ the difference between the exact quantile there and the linear interpolation and μ the row’s mean loss (the criterion of the elastic tables, crate::electron::elastic::table). The test stopping_power_from_the_loss_cdf_matches_the_model holds stopping_power_ev_per_m to 1e-3 of S(E) of #18 at every grid energy.

Rows where λ⁻¹ = 0 (no loss allowed) are stored empty, and so are rows with T at or below the lowest loss of the model’s rate (module docs, “Energy axis”). The table covers the direct single-pole model only: no exchange and no surface or band-gap structure (see super::penn). Inner shells get one table each, beside the valence table, from build_shell_channel_tables. A model with the exchange correction enabled (SinglePolePenn::with_exchange) is rejected by build_inelastic_table and MomentumTransferSampler.

Energy axis

The models take T, the kinetic energy above the model’s Fermi level, with kinematics on T' = T + E_F and losses up to T (S2017 eqs. (2)-(3), super::penn). By default (EnergyAxis::ModelFermiLevel) a table is stored on that axis: the row at E is the model’s at T = E.

EnergyAxis::BandBottom stores it on the band-bottom axis instead: the row at E is the model’s at T = E - E_F, so its kinematics are on T' = E and its losses stop at E - E_F; rows with T <= 0 are stored empty. With E_F the band’s, this is the convention of the table compiler of Nebula, compile_full_imfp_icdf in cstool/dielectric_function/compile.py (Nebula-simulator/cstool commit 0c739eb3fcc3fe5297e74c601ac4a9546db596cf, BSD-3-Clause): it evaluates every row at the electron’s kinetic energy K (the momenta of the kinematic limits are those of K and K - ω) and keeps only the losses ω < K - F (its “Fermi correction”), F the energy its caller compile_full_penn (apps/cstool.py, same commit) passes, the band’s get_min_excitation() (cstool/input_data/band_structure.py): the Fermi energy of a metal, the conduction-band bottom W_v + E_g of an insulator or semiconductor. The model’s own E_F is the F of the axis, so there is one Fermi energy, not two. The table builder knows no band; the caller sets the model’s Fermi energy (lindhard run does so per material, docs/cli.md).

Rows below the ELF table. The sampled loss window of a row is (W_min, T], W_min the lowest loss the model’s rate λ⁻¹(T) counts, so that the stored rate and the sampled losses cover the same window:

• single-pole Penn: W_min = elf_min, the lowest tabulated ELF energy. Its DIIMFP is zero below it (super::penn), so its rate counts no loss there;
• full Penn: W_min = 1e-6 T (super::full_penn, its LOWEST_LOSS_FRACTION). The DIIMFP of the full model is not zero below elf_min: each Lindhard gas of the expansion has an electron-hole continuum down to ω -> 0 at q > 0 (S2017 eq. (4)), and the model’s rate and stopping power integrate every loss from 0. A window opened at elf_min (until #361) left the losses below it out of the sampled distribution while the stored rate counted them: on the Si ELF tabulated from 0.5 eV the share of λ⁻¹ below 0.5 eV is 4.7e-3 at T = 10 eV, and the table’s S = λ⁻¹ ⟨W⟩ was biased high by the same amount. The lower limit is a numerical floor, not a model cutoff (a literal zero cannot bound the ln W axis): p ∝ ω at small ω, so the rate below W falls as W² and the stopping power as W³, and the omitted share of both is below the model tolerance by orders of magnitude (tests/penn_full.rs bounds it against the direct DIIMFP). The row’s rate is the model’s own; the sampled distribution and ⟨W⟩ cover (1e-6 T, T], which the model’s 0 < ω <= T equals up to that omitted share. The DIIMFP grid reaches down to the floor of the lowest row energy, so these losses are read from it, not computed directly;
• Mermin-ELF: W_min = 1e-8 T, the lower limit of the ω integral of its rate and stopping power (super::mermin, the same constant). Its fitted Drude-Lorentz ELF has no cutoff at elf_min (the omission of the threshold step is a documented choice of that model), so the DIIMFP below elf_min is the fit’s extension below the tabulated range and is part of both the rate and the sampled loss distribution (#342). Before #342 the window opened at elf_min for every model, and a Mermin row drew the losses its rate counts below elf_min as losses above it: for the Si ELF tabulated from 0.5 eV (default fit, model Fermi energy 13.46 eV) the window (0.5 eV, T] held 28 % of the rate at T = 0.6 eV, 72 % at 1 eV, 92 % at 2 eV and 98.4 % at 5 eV.

A row with 0 < T <= W_min has an empty window, so no loss can be drawn there, and it is stored empty (inverse mean free path 0, no quantiles), as a row with T <= 0 is, on either axis (#340): without this the ln W window [ln W_min, ln T] is reversed (or of zero width at T = W_min) and the density integrates to a negative value (or 0). For the Penn single-pole model these are the rows with 0 < T <= elf_min. A Mermin or full Penn row has W_min < T at every T > 0, so it is never empty for this reason: its rows with 0 < T <= elf_min, stored empty from #340 to #342 (Mermin) or #361 (full Penn), carry the model’s rate and losses in (W_min, T].

Rate refinement

On the band-bottom axis the rate λ⁻¹(E) is zero at E_F and climbs by two orders of magnitude within about 20 eV above it, and the transport interpolates it linearly between rows (crate::electron::transport, “Free flight”). A log grid fine enough there is far finer than needed everywhere else (#291, #339). InelasticTableOptions::rate_refinement (RateRefinement, off by default) adds rows to this one table instead: every cell of the grid that overlaps E - E_F in RateRefinement::window_ev (RATE_WINDOW_EV by default) is bisected at its midpoint in E while the linear interpolation of the exact row rate between the cell’s ends differs from the exact rate at any of the cell’s test points by more than RateRefinement::tolerance (RATE_TOLERANCE) of it. The test points are the cell’s quartiles and midpoint (rate_test_points) and, where those pass, the estimated position of the largest error: the vertex of the parabola through the largest of the sampled errors and its two neighbours (rate_error_peak). Children that no longer overlap the window are not tested, so the bisection never approaches the zero of the rate at E_F, and cells outside the window are never split. Only the scalar rate is evaluated while refining (the row rule of the table: zero for an empty row); each level’s test points are evaluated in parallel and merged in sorted order, so the grid is the same on any thread count. The rows are then built as any other, from the merged grid, which contains the given grid. A refinement that would add more than RateRefinement::max_rows rows is an error naming the cell that still fails (a bound against a rate the bisection cannot resolve, not a tuned value), and so is a cell still to test that can no longer be bisected because its midpoint rounds to an end. A step in the rate splits one cell per level, so it reaches rounding (about 50 levels) before any useful row bound: that second error, not the row bound, is what stops it.

Rate onset (#357). The single-pole rate is identically zero for E - E_F <= elf_min (“Rows below the ELF table”, W_min = elf_min) and positive above it: a step at the onset E_F + elf_min, which the transport’s linear interpolation from the zero row below it cannot follow at any row density. Where the onset lies at or above the window start (an ELF tabulated from 2 eV or more, at the default window), a cell is tested only if it overlaps the window and its lower end is at or above the onset: the cell straddling the onset and the cells below it are left as given, the interpolation across the straddling cell being the model’s own step. Testing them would bisect towards the step until its cell is exhausted to rounding, the error above, and so refuse an ELF that builds without the option. Where the onset lies below the window start nothing changes: the bisection’s children that leave the window are not tested, so it never reaches the step, as it never reaches the zero at E_F. Every committed ELF starts at or below 0.5 eV (validation/data/optical/), so for none of them does the onset change the refined grid. The Mermin rate is positive for every T > 0 (it counts losses down to 1e-8 T), and so is the full Penn rate (down to 1e-6 T, #361: it was grouped with the single-pole onset before), so those models have no onset and keep the full window.

The test is a criterion, not a bound. A test of the midpoint alone (#351) left cells whose error peaked off the midpoint above the tolerance, and so did one of the quartiles and midpoint (the Al single pole from 5 eV at 10 per decade, 1.0008 % at 6.1 eV, between the first quartile and the midpoint of its cell); the parabolic vertex locates such a peak without a tuned number of samples. The acceptance check is the dense evaluation of the rate mode of lindhard-cli/examples/inelastic_low_energy.rs (every 0.1 eV over the window, docs/validation.md). The table’s other grid (the elastic table of a run, or a table without the option) is unchanged: the transport reads each table on its own CrossSectionTable::energy_ev.

Momentum-transfer sampler

Given (E, W), the DIIMFP integrand in q is, in u = ln q (S2017 eq. (2), as in super::penn), Im[-1/ε(q, W)] on q- <= q <= q+, q± = k ± k', k = sqrt(2T'), k' = sqrt(2(T' - W)) (Hartree units, T' = E + E_F), cut to the q where the plasma frequency of the pole lies in the tabulated ELF range. MomentumTransferSampler is a direct sampler built per (E, W): the integrand is sampled on adaptive nodes (breaks at the kinematic limits and at every ELF knot, 8 panels per smooth segment to start, the same refinement rule as above with the zeroth moment only and tolerance DEFAULT_Q_TOLERANCE), linear in ln q between nodes, and inverted exactly. The density used is therefore the model’s to within that tolerance, and the histogram test compares it with an independent quadrature of SinglePolePenn::loss_function at 1e6 samples. It is not a table: building one costs the integrand evaluations of one DIIMFP slice, so a transport loop should keep the sampler for repeated draws at one (E, W) and expect that cost for each new pair. A tabulated (E, W, u) form is left for later work.

The deflection angle follows from momentum conservation of the primary, q² = k² + k'² - 2 k k' cos θ, i.e. cos θ = (k² + k'² - q²)/(2 k k') (MomentumTransferSampler::cos_theta), clamped to [-1, 1]; q = q- is θ = 0 and q = q+ is θ = π. This is the elementary kinematics of the loss W (the final wavenumber k' fixed by energy conservation) with the model’s own q±.

Determinism

Rows (and the profiles of a full-Penn DiimfpGrid) are built in parallel with rayon and collected in grid order; each row, and the probability refinement (a parallel map over rows followed by a serial split decision), is a fixed sequence of floating-point operations, so tables are bit-identical on any thread count (tested). Samplers consume uniforms u supplied by the caller (for example from crate::rng::stream) and keep no state of their own.

Structs

InelasticTableOptions

Grids and tolerances of a table build.

MomentumTransferSampler

Samples the momentum transfer q of an inelastic collision with incident kinetic energy E and loss W (module docs). Built once per (E, W).

RateRefinement

Extra rows that resolve the inelastic rate λ⁻¹(E) near the Fermi level on the band-bottom axis (module docs, “Rate refinement”; #339).

ShellResolvedTables

The tables of a ShellResolvedChannels: the valence channel’s and one per inner shell, in the order of ShellResolvedChannels::inner_shells (build_shell_channel_tables).

Enums

EnergyAxis

The energy axis a table is stored on (module docs, “Energy axis”).

InelasticTableError

Errors from the table builder and the samplers.

Constants

DEFAULT_DENSITY_TOLERANCE

Default tolerance of the loss-density model (see the module docs).

DEFAULT_LOGISTIC_STEP

Default logistic step of the starting probability grid.

DEFAULT_PROBABILITY_TAIL

Default probability tail of the starting grid (see logistic_probability_grid).

DEFAULT_Q_TOLERANCE

Default tolerance of the momentum-transfer density model.

DEFAULT_RATE_REFINEMENT_MAX_ROWS

The default bound on the rows the rate refinement may add (a safety bound against a rate that fails the test in many cells at once; not a tuned value): 1000 rows, so cells of 18 meV on average across the 18 eV of RATE_WINDOW_EV, far below the 0.1 eV step of the rate gate. A step in the rate splits one cell per level and is stopped earlier, when its cell is exhausted to rounding (module docs, “Rate refinement”).

DEFAULT_REFINE_TOLERANCE

Default target of the probability-grid refinement, relative to the mean loss of each row.

INITIAL_LOSS_PANELS

Number of equal panels (in ln W) the loss density starts from.

MAX_NODES

Upper bound on the number of nodes of one density model and on the size of the refined probability grid.

RATE_TOLERANCE

The tolerance of the rate refinement (#339) and of the rate gate of lindhard-cli/examples/inelastic_low_energy.rs (#291), fixed there before any value was measured: the largest relative error of the interpolated 1/λ over RATE_WINDOW_EV above the Fermi level must not exceed 1 %. It is the statistical error of the most precise committed default δ entries, rounded up: in validation/experiments/se_yield_results.json the Poisson floor of Al default is 0.0585 on δ = 6.85 at 600 eV (0.85 %) and 0.0571 on δ = 6.52 at 400 eV (0.88 %). Other configurations are more precise: Al barrier-off reaches 0.60 % (600 eV), and at 400 eV the Al barrier-off, cutoff-band-bottom, phi-low and phi-high entries are 0.63 to 0.90 %; the Au, Cu and Si default entries at 400 eV are 1.4, 1.7 and 2.8 %. A grid error below the tolerance cannot be told from the noise of those tables if δ responds at most in proportion to the rate in that window; that response is not measured here.

RATE_WINDOW_EV

The window of the rate refinement and of the rate gate (#291, #339): E - E_F from 2 to 20 eV. The upper end covers the steep rise of the rate above the Fermi level, which for the Al single pole peaks in error at 15 to 17 eV. The lower end lies below the transport cutoff of every default δ run of validation/experiments/se_yield.py (1 eV above the vacuum level, at least 5 eV above E_F); the cutoff-band-bottom runs read the rate down to 1 eV above E_F. Below 2 eV the rate starts from zero at E_F, and a relative tolerance there would bisect towards that node without end.

Functions

build_inelastic_table

Build the inelastic energy-loss table of penn on the grid of options.

build_inelastic_table_for_model

Build the inelastic energy-loss table of any of the three ELF-extension models (PennInelastic: single-pole Penn, full Penn or Mermin-ELF) on the grid of options, by the same procedure as build_inelastic_table (module docs), with that model’s DIIMFP and inverse IMFP in place of the single pole’s. For PennInelastic::SinglePole the rows are those of build_inelastic_table; only the model and provenance strings differ, and they carry PennInelastic::model_identity, so a table names the algorithm it was built with.

build_shell_channel_tables

Build the energy-loss tables of every channel of channels on the grid of options: the valence table, and for inner shell i the table of ShellResolvedChannels::shell_model (that shell’s own optical ELF with the valence model’s settings), wrapped with the shell’s identity, binding energy and binding_provenance (the origin of the binding energies, for example crate::electron::data::SubshellBindingTable::provenance).

default_probability_grid

The default starting probability grid: logistic_probability_grid with DEFAULT_PROBABILITY_TAIL and DEFAULT_LOGISTIC_STEP.

loss_density_node_count

The number of nodes of the loss-density model of the row at energy_ev (eV, on the model’s own axis) of a table of model built with density tolerance density_tolerance, reading the DIIMFP from grid where it covers the loss as a table build does (pass the grid of FullPenn::diimfp_grid for the table’s energies, or None). The count is bounded by MAX_NODES (a few more: the last refinement pass splits every failing panel); a row that reaches it is refused by the table build (InelasticTableError::DensityNodeLimit). Each node not covered by the grid costs one direct DIIMFP evaluation.

mean_loss_ev

∫_0^1 W(u) du by the trapezoid rule of the table’s linear interpolation.

stopping_power_ev_per_m

The stopping power λ⁻¹ ⟨W⟩ (eV/m) of row i of an energy-loss table, with ⟨W⟩ from mean_loss_ev; zero for a zero-rate row. None if i is out of range.