feat(inelastic): refine band-bottom tables in lindhard run (#339)

CHANGELOG.md

@@ -12,6 +12,21 @@ one version).
 - `electron_psf_parameters.csv` gains `converged` (1/0) and `iterations` rows per
   fit, so the standalone CSV shows whether a fit hit its iteration budget (#347).
 
+### Changed
+
+- `lindhard run` builds the inelastic table of every material with a `band`
+  with the rate refinement (`RateRefinement::default()`, #339): rows are
+  added between 2 and 20 eV above the Fermi level until the interpolated
+  inelastic rate is within 1 % of the exact one. **Results of runs with a
+  `band` change**; the elastic table and the tables of materials without a
+  band are unchanged, bit for bit. The committed secondary-yield δ(E) tables
+  (`validation/experiments/se_yield_results.json`) predate it and are stale
+  until their re-run (#287). `electron_tables.csv` has a line per inelastic
+  row, with an empty elastic column at the added rows (the elastic value is
+  now matched by energy). The table-cache key gains `inelastic.rate_refinement`
+  (`KEY_VERSION` 4), so cached tables are rebuilt once, and a cached
+  inelastic table must contain the run's grid in order rather than equal it.
+
 ### Fixed
 
 - Python `RunResult.write()` now follows the command's output lifecycle: it
@@ -32,16 +47,17 @@ one version).
   of it): rows added to one band-bottom inelastic table from the Fermi
   level upward. Every grid cell overlapping `E - E_F` = 2 to 20 eV is
   bisected at its midpoint while the linear interpolation of the exact rate
-  misses it there by more than 1 %; more than `max_rows` added rows is an
+  misses it by more than 1 % at the cell's quartiles, midpoint or estimated
+  error peak (the vertex of the parabola through the largest sampled error
+  and its neighbours); more than `max_rows` added rows is an
   error, and so is a failing cell bisected down to rounding (a step in the
-  rate), each naming the cell. The grid is the same on any thread count. Off by default, so every
-  existing table is unchanged, and refused on the model's own axis.
+  rate), each naming the cell. The grid is the same on any thread count. Off by default, and
+  refused on the model's own axis.
   `RATE_TOLERANCE` and `RATE_WINDOW_EV` (the #291 gate, moved from the
   example into `electron::inelastic::table`) are its defaults. The `rate`
   mode of `lindhard-cli/examples/inelastic_low_energy.rs` also evaluates
   each coarse table refined (`docs/validation.md`, "Table grid above the
-  Fermi level"). `lindhard run` does not use it yet: the Al single pole at
-  the `se_yield.py` grid misses the 1 % gate by 0.0008 points.
+  Fermi level"); every refined table there passes the 1 % gate.
 - `[electron.materials.<name>.mermin_fit]` (#306): per-material options of
   the `mermin-melf` oscillator fit, `oscillators` (1 to 16, default 3),
   `weighting` (`"relative"`, the default, or `"uniform"`) and

docs/cli.md

@@ -512,7 +512,9 @@ 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.
+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).
 
 | Key | Default | Meaning |
 |---|---|---|
@@ -539,11 +541,15 @@ grid above the Fermi level"; reproduce with
 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`, but one of them (the Al single pole
-from 5 eV at 10 points per decade) misses 1 % by 0.0008 points, so
-`lindhard run` does not use it yet (#339). Until then a denser grid helps
-but does not close the gap: 40 points per decade cuts the error by a
-factor of 2 to 4 from 20.
+`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
@@ -816,8 +822,10 @@ skip them. Parameters from a fit with `converged` 0 need inspection of the fit
 diagnostics before use.
 
 `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 grid energy (the stopping power is
-`λ⁻¹ ⟨W⟩` of the stored loss distribution).
+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`)
 
@@ -851,7 +859,9 @@ the table depends on:
   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).
+  (`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
@@ -862,7 +872,8 @@ 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. A file found under the run's key that fails any of these checks is an
+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

docs/validation.md

@@ -1994,38 +1994,51 @@ result worse, not better; no Cu, Si or Au set changes status.
   rate needs rows from `E_F` upward, not a denser uniform grid.
 
   **Rate refinement (#339).** `InelasticTableOptions::rate_refinement`
-  (`electron::inelastic::table`, "Rate refinement"; off by default) gives
-  one band-bottom table those rows: every grid cell overlapping `E - E_F`
-  = 2 to 20 eV is bisected at its midpoint while the linear interpolation
-  of the exact rate misses the exact rate there by more than 1 %
-  (`RATE_TOLERANCE` and `RATE_WINDOW_EV`, the gate's own constants, now in
-  the library and read by the example). The `rate` mode builds each coarse
-  table again with it (`+refined`: the coarse rows and every row the
-  construction adds) and gates it the same way. Largest error over 2 to 20
-  eV (where it falls) and rows of the refined table (the coarse rows
+  (`electron::inelastic::table`, "Rate refinement") gives one band-bottom
+  table those rows: every grid cell overlapping `E - E_F` = 2 to 20 eV is
+  bisected at its midpoint while the linear interpolation of the exact
+  rate misses the exact rate by more than 1 % (`RATE_TOLERANCE` and
+  `RATE_WINDOW_EV`, the gate's own constants, in the library and read by
+  the example) at any of the cell's test points: its quartiles and
+  midpoint, and, where those pass, the vertex of the parabola through the
+  largest sampled error and its two neighbours (the estimated peak). The
+  `rate` mode builds each coarse table again with it (`+refined`: the
+  coarse rows and every row the construction adds) and gates it the same
+  way. The first version (#351) tested the midpoint only and missed in
+  four tables, by 1.0008 % (Al single pole, 5 eV, 10 per decade, at 6.1
+  eV) to 1.43 %: there the error peaks off the midpoint. Testing the
+  quartiles as well still left the first of them at 1.0008 % (its peak lies
+  between the first quartile and the midpoint of the cell 5.3 to 7.2 eV,
+  which miss by 0.84 and 0.99 %) and Al Mermin from 10 eV at 10 per decade
+  at 1.03 %; the estimated peak closes both. No constant was tuned: the
+  tolerance and window are #291's, and the test points are the next level
+  of the bisection and the parabolic vertex. Largest error over 2 to 20 eV
+  (where it falls) and rows of the refined table (the coarse rows
   bracketing 0.5 to 20 eV above `E_F` plus the rows added), at the two
   grids that matter, the `se_yield.py` grid (5 eV, 10 per decade) and the
-  library default (10 eV, 20 per decade):
+  library default (10 eV, 20 per decade), 2026-10-10:
 
   | Table | 5 eV, 10 / decade | rows (coarse + added) | 10 eV, 20 / decade | rows (coarse + added) |
   |---|---|---|---|---|
-  | Al single-pole | **1.0008 % (6.1), fails** | 6 + 22 | 0.84 % (3.1) | 11 + 24 |
-  | Al Mermin | 0.66 % (18.6) | 6 + 10 | 1.00 % (2.7) | 11 + 6 |
-  | Si single-pole | 0.84 % (15.2) | 6 + 25 | 0.96 % (3.9) | 10 + 21 |
-  | Si Mermin | 0.93 % (3.0) | 6 + 7 | 0.72 % (15.8) | 10 + 5 |
-
-  At the other grids of the mode (5 or 10 eV, 10 to 80 per decade) the
-  refined tables pass except the Al single pole from 5 eV at 20 per
-  decade (1.43 % at 17.6 eV), Al Mermin from 10 eV at 10 per decade
-  (1.03 %, 4.6 eV) and Si Mermin from 10 eV at 10 per decade (1.40 %,
-  18.6 eV); the uniform grids alone fail everywhere below 80 per decade.
-  The midpoint test is a criterion, not a bound: in those cells the error
-  peaks away from the midpoint the test samples. By the stop rule of #339
-  the criterion is not tightened by an empirical factor and neither the
-  tolerance nor the window is widened, so **`lindhard run` does not use
-  the refinement yet** and the grid of these runs is unchanged. #339 stays
-  open for that decision; #291 and #287 are not settled by it.
-  The committed results are unchanged.
+  | Al single-pole | 0.967 % (16.4) | 6 + 23 | 0.843 % (3.1) | 11 + 24 |
+  | Al Mermin | 0.656 % (18.6) | 6 + 10 | 0.880 % (4.8) | 11 + 7 |
+  | Si single-pole | 0.843 % (15.2) | 6 + 25 | 0.953 % (4.9) | 10 + 22 |
+  | Si Mermin | 0.934 % (3.0) | 6 + 7 | 0.721 % (15.8) | 10 + 5 |
+
+  Every refined table of the mode passes (5 and 10 eV, 10 to 80 per
+  decade, 32 tables), including the three that failed with the midpoint
+  test outside these grids: the Al single pole from 5 eV at 20 per decade
+  0.772 % (15.9; 11 + 23 rows), Al Mermin from 10 eV at 10 per decade
+  0.893 % (2.7; 7 + 11) and Si Mermin from 10 eV at 10 per decade 0.842 %
+  (18.2; 6 + 10). The uniform grids alone fail everywhere below 80 per
+  decade. The test is a criterion, not a bound; the dense gate is the
+  check. **`lindhard run` now builds the inelastic table of every material
+  with a band with the refinement** (`RateRefinement::default()`; the
+  elastic table and tables without a band keep the grid), which settles the
+  grid question of #291 for the δ runs: the `se_yield.py` grid (5 eV, 10
+  per decade) needs no change. The committed δ tables below were computed
+  before it, so they are stale until the re-run of #287; this change does
+  not re-run them.
 - Secondary generation: Kieft and Bosch (2008) as implemented
   (`secondaries = "kieft-bosch"`); transport cutoff 1 eV above the vacuum
   level (`cutoff_reference = "vacuum-level"`) unless a row says otherwise.

lindhard-cli/examples/inelastic_low_energy.rs

@@ -53,7 +53,8 @@
 //!   builder refuses is printed as refused. Each coarse table is followed by
 //!   the same rows built with the library's rate refinement
 //!   (`InelasticTableOptions::rate_refinement`, `RateRefinement::default()`,
-//!   #339; `+refined`; `lindhard run` does not use it yet): its grid is the
+//!   #339; `+refined`; `lindhard run` uses it for every table on the
+//!   band-bottom axis): its grid is the
 //!   coarse rows and every row the refinement adds, and `rows` counts them.
 //!   The largest error is printed to three decimals, so that a verdict at
 //!   the tolerance can be read. The verdict that 20 and 40 points per decade fail

lindhard-cli/src/electron.rs

@@ -22,6 +22,7 @@ use lindhard::electron::elastic::table::{
 use lindhard::electron::elastic::{SalvatDhfs, SolverOptions};
 use lindhard::electron::inelastic::table::{
     build_inelastic_table_for_model, mean_loss_ev, EnergyAxis, InelasticTableOptions,
+    RateRefinement,
 };
 use lindhard::electron::inelastic::{PennAlgorithm, PennInelastic};
 use lindhard::electron::transport::{LayerTables, RunMetadata, Transport};
@@ -236,13 +237,30 @@ fn inelastic_table(
         other => PennInelastic::try_new(other, m.optical_elf.clone())?,
     }
     .with_fermi_energy_ev(fermi_ev)?;
+    let mut options = InelasticTableOptions::new(r.table_energy_ev.clone()).with_axis(axis);
+    options.rate_refinement = rate_refinement(r, m);
     Ok(build_inelastic_table_for_model(
         &model,
         &m.material,
-        &InelasticTableOptions::new(r.table_energy_ev.clone()).with_axis(axis),
+        &options,
     )?)
 }
 
+/// The rate refinement of the inelastic table of `m` (#339): the library's
+/// default ([`RateRefinement::default`]) on the band-bottom axis of a
+/// material with a band, where the rate climbs from zero at the Fermi
+/// level; `None` on the model's own axis, whose table keeps exactly the
+/// rows of `electron.tables`. The refined table has the grid's rows and the
+/// rows the refinement adds between 2 and 20 eV above the Fermi level
+/// (`lindhard::electron::inelastic::table`, "Rate refinement"); the elastic
+/// table keeps the grid. Part of the table cache key.
+pub fn rate_refinement(
+    r: &ResolvedElectron,
+    m: &ResolvedElectronMaterial,
+) -> Option<RateRefinement> {
+    (inelastic_axis(r, m).1 == EnergyAxis::BandBottom).then(RateRefinement::default)
+}
+
 /// The model's Fermi energy and the table axis of the inelastic table of `m`.
 ///
 /// For a material with a band, the table is built on the band-bottom axis
@@ -740,14 +758,15 @@ pub fn tables_csv(r: &ResolvedElectron, sim: &ElectronSimulation) -> String {
          inelastic_mean_loss_ev,inelastic_stopping_ev_per_nm\n",
     );
     for (m, t) in r.materials.iter().zip(&sim.tables) {
-        // Both tables share the grid of electron.tables.
+        // One line per inelastic row. The elastic table has the rows of
+        // electron.tables, which the inelastic table contains; the rows the
+        // rate refinement adds (#339) have no elastic value.
+        let elastic = t.elastic.energy_ev();
         for (i, &e) in t.inelastic.energy_ev().iter().enumerate() {
-            let el = t
-                .elastic
-                .energy_ev()
-                .get(i)
-                .filter(|&&x| x == e)
-                .map(|_| t.elastic.inverse_mfp_per_m()[i] * NM);
+            let el = elastic
+                .binary_search_by(|x| x.total_cmp(&e))
+                .ok()
+                .map(|j| t.elastic.inverse_mfp_per_m()[j] * NM);
             let inv = t.inelastic.inverse_mfp_per_m()[i];
             let mean = t
                 .inelastic