feat: spherical (radial) geometry for VFP-1D + SNB/Spitzer/kinetic heat-flux comparison

[joglekara]

Summary

Adds a radial-coordinate (spherically symmetric) mode to the vfp-1d solver, file-based profile initialization, and Spitzer-Härm / SNB heat-flux diagnostics — so a kinetic vs SNB vs Spitzer heat-flow calculation can be run directly on tabulated ne(r), Te(r), and mass-density profiles.

Physics

For nl=1 (f0, f10), expanding the Vlasov equation in spherical symmetry gives exactly the Cartesian system except the f0 spatial term becomes the radial divergence:

$$\partial_t f_0 + \frac{v}{3}\frac{1}{r^2}\partial_r(r^2 f_{10}) = C[f_0], \qquad \partial_t f_{10} + v,\partial_r f_0 - E,\partial_v f_0 = C[f_{10}]$$

The angular coupling term $\propto (1-\mu^2)/r$ only feeds $f_2$ into the $f_{10}$ equation, which is dropped at nl=1 — velocity space and all collision operators are untouched. The divergence is discretized in conservative finite-volume form with exact cell volumes $(r_+^3 - r_-^3)/3$, so volume-weighted conservation telescopes; the flux through $r=0$ vanishes identically. Spherical geometry requires xmin: 0 and boundary: reflective ($f_{10} = E = 0$ at $r = 0$ by symmetry and at the outer wall), which the existing reflective BC machinery already provides.

Changes

• grid.geometry: cartesian | spherical — new Grid field with validation; OSHUN1D.div_e2c applies the radial divergence in the f0 equation.
• basis: file profiles — two-column (coordinate, value) CSV/text files for n and T, unit-converted with pint and interpolated onto the grid (clamped outside the tabulated range).
• density.ion block — a tabulated mass-density profile plus the ion mass number A sets $n_i(r) = \rho/(A m_u)$ and a spatially varying $Z(r) = n_e/n_i$ in the state (stored relative to units.Z, consistent with the FLMCollisions convention).
• adept/vfp1d/heat_flux.py — post-processing (pure numpy) diagnostics:
  • local Spitzer-Härm(-Epperlein-Haines) flux, consistent by construction with the existing kappa metric;
  • multigroup SNB flux (Schurtz et al. 2000; the "separated" variant of Brodrick et al. 2017 with r=2, $\xi(Z)$ e-e correction), geometry-aware tridiagonal solves, using the same normalized collision rates as the kinetic solver;
  • comparison netCDF + plot (plots/fields/heat-flux-comparison.png) and q_ratio_kinetic_over_sh, q_ratio_snb_over_sh metrics on every vfp-1d run.
• Example config configs/vfp-1d/spherical-heatflow.yaml set up for the profile-driven comparison (fill in the three CSV paths, reference units, and A).
• Docs — vfp1d config reference updated (geometry, file basis, ion block, diagnostics).

CD plasma specifics (second commit)

• units.logLambda: lee-more — Lee & More 1984 Coulomb logarithm: $\ln\Lambda = \max[2, \tfrac{1}{2}\ln(1 + b_{max}^2/b_{min}^2)]$ with Debye-Hückel screening (electron + ion) floored at the ion-sphere radius, and $b_{min}$ = max(classical closest approach, thermal de Broglie). At conduction-zone conditions (1e21/cc, 1 keV) it gives 7.09 vs NRL 6.73.
• density.ion.mixture — atomic mixtures represented by one effective ion preserving quasineutrality ($Z n_i = n_e$) and e-i collisionality ($Z^2 n_i = \sum n_i Z_i^2 = Z_{eff} n_e$). For fully ionized CD (40% C / 60% D): $\bar Z = 3$, $Z_{eff} = \langle Z^2\rangle/\langle Z\rangle = 5$, A = 6.304. units.Z is the transport $Z_{eff}$ (= 5), which is what enters the EH and SNB coefficients. The mass-density profile cross-checks full-ionization quasineutrality against $n_e(r)$ and warns on >5% mismatch.
• Example config updated for the CD conduction zone.

Tests

New tests/test_vfp1d/test_spherical.py (all pass on CPU):

• uniform plasma in spherical geometry stays uniform (geometric terms vanish without a radial flux);
• central hotspot: outward heat flow, volume-weighted energy conserved to 1e-3, SNB/SH within an order of magnitude of kinetic q at the peak;
• SNB → Spitzer-Härm in the collisional limit (2%), flux inhibition in the nonlocal limit;
• file-basis + ion-block normalizations (density moment of f0 matches ne(r); Z, ni correct);
• config validation (xmin≠0, periodic BC rejected).

Existing Cartesian benchmarks test_kappa_eh[False-1] and test_kappa_eh_reflective[False-1] still pass (and now exercise the heat-flux diagnostic path).

🤖 Generated with Claude Code