imas: exclude fast-ion charge from the impurity split - #46
Open
d-burg wants to merge 3 commits into
Open
Conversation
effective_impurity_charge uses (ne - ni) as the impurity charge. With a beam population that difference also contains the fast ions, so nz is inflated and Z_imp drifts off the real impurity charge. Use ne - sum_s Z_s n_s^fast instead. On DIII-D 174823 (6.3 MW NBI) the fast fraction reaches 28% through the L->H transition and the pressure-completeness gap falls from 12.6% to 0.6%, below the 2% guard that was rejecting those slices. Zeff is left on the full ne: it feeds the bootstrap and forward solve, so correcting it there is a separate change. That leaves Z_imp still low (4.2 vs 6 at the worst slice) -- the pressure is right, the impurity density is not.
… subtraction alone was half the fix
Zeff in the IDS convention is the thermal-species numerator over the FULL
ne, so passing the thermal ne_th with the unrenormalized Zeff recovers
only half the bias: measured Z_imp raw 1.95, half-applied 3.19, true 6.00
for C6 at 25 % fast fraction (adversarial review, verified numerically).
The correct inversion uses zeff * ne / ne_th, and it is LOCAL to the
Z_imp/p_imp derivation -- the Zeff consumed by the bootstrap and the
forward solve deliberately stays the full-ne one, so the previously
documented residual bias ('4.2 vs 6 at the worst slice') was never the
forced consequence the PR text attributed to that choice.
The inversion now lives in physics.impurity_charge_with_fast_ions (with
the z_fast >= ne surfaces excluded by the validity mask) and imas.py
wires it in; the tests assert EXACT recovery, pin the half-applied form
as a regression marker, and check the imas.py wiring by source -- the
first round only asserted 'closer than raw', which blessed the half-fix.
Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Baseline gains z_fast (the reader's fast-ion charge profile, fixed across draws like p_fast) and it is plumbed through generate_bouquet into perturb_kinetic_equilibrium. All impurity math now runs on ne - z_fast with the SAME thermal-derived Z_imp the reader stores: - run.py baseline forward solve and state-anchor assemblies; - the per-draw pres_tmp impurity term and the baseline reference assembly (sigma=0 draw pressure again equals the reader's p_recon -- the skew was e*z_fast*ti/Z_imp); - the Zeff-primary channel: Z_imp via impurity_charge_with_fast_ions, ni derived from thermal quasineutrality with the full-ne Zeff draw (ni = (Z ne_th - Zeff ne)/(Z-1), exact inverse of the baseline set), and the draw clip generalized to ne_th/ne <= Zeff <= Z_imp ne_th/ne (reduces to [1, Z_imp] at z_fast=0); - main_ion_density_from_zeff grows an optional z_fast (None = exactly the old behaviour; recon path unaffected). 4 new tests: exact ni/nz recovery through the thermal inversion, the z_fast=None reduction, the Baseline field, and source-level wiring pins for all four consumer sites. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
There was a problem hiding this comment.
🟢 Approval recommended
The functional changes are coherent end-to-end (reader → baseline → draw/solve paths) and are backed by targeted regression tests; remaining feedback is limited to documentation/test robustness.
Pull request overview
This PR corrects the IMAS impurity inversion so fast-ion charge is not mistakenly attributed to impurities, and propagates the resulting “thermal-electron” density (ne_th = ne - z_fast) through downstream pressure construction paths.
Changes:
- Add
physics.impurity_charge_with_fast_ions()and extendmain_ion_density_from_zeff(..., z_fast=...)to support fast-ion dilution while keeping IMAS Zeff normalized to fullne. - Wire IMAS baseline reading to compute
z_fast = Σ Z_s n_s^fast, compute(Z_imp, ne_th)via the helper, and usene_thfor impurity pressure. - Thread
z_fastthrough bouquet generation / equilibrium perturbation / forward-solve assembly and add regression tests for the corrected behavior and wiring.
File summaries
| File | Description |
|---|---|
| tests/test_fast_ion_dilution.py | Adds regression tests covering corrected impurity inversion and wiring across IMAS + solver paths. |
| bouquet/physics.py | Introduces impurity_charge_with_fast_ions and extends Zeff→ni inversion to account for z_fast. |
| bouquet/io/imas.py | Computes z_fast from IMAS ions and uses the helper to derive Z_imp and thermal ne_th for impurity pressure. |
| bouquet/baseline.py | Adds Baseline.z_fast to carry fast-ion charge density on the kinetic grid. |
| bouquet/TokaMaker_interface.py | Threads z_fast into perturb/draw plumbing, Zeff-primary ni derivation bounds, and impurity pressure calculation. |
| bouquet/run.py | Ensures forward-solve and sigma=0 consistency pressure assembly uses thermal ne when z_fast is present. |
Review details
- Files reviewed: 6/6 changed files
- Comments generated: 3
- Review effort level: Lite
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
Comment on lines
1613
to
1617
| max_li_iter=_MAX_LI_ITER, | ||
| psi_N_kinetic=None, | ||
| p_fast=None, | ||
| z_fast=None, | ||
| j_NBI=None, |
Comment on lines
+500
to
+516
| def impurity_charge_with_fast_ions(ne, ni, zeff, z_fast): | ||
| """``(Z_imp, ne_th)`` when fast ions carry part of the neutralization. | ||
|
|
||
| With a beam population, quasineutrality reads | ||
| ``ne = ni + Z_imp*nz + z_fast`` -- only ``ne_th = ne - z_fast`` is | ||
| neutralized by THERMAL ions. A ``zeff`` normalized to the FULL ``ne`` | ||
| (the IMAS convention: thermal-species numerator over total electron | ||
| density) must therefore be renormalized to ``zeff * ne / ne_th`` before | ||
| the single-impurity inversion; passing the thermal ``ne`` with the | ||
| full-``ne`` ``zeff`` recovers only half the bias (measured: raw 1.95, | ||
| half-applied 3.19, true 6.00 for C6 at 25 % fast fraction). Surfaces | ||
| where ``z_fast >= ne`` get a non-finite renormalized zeff and are | ||
| excluded by :func:`effective_impurity_charge`'s own validity mask. | ||
| """ | ||
| ne = np.asarray(ne, dtype=float) | ||
| z_fast = np.asarray(z_fast, dtype=float) | ||
| ne_th = np.maximum(ne - z_fast, 0.0) |
Comment on lines
+79
to
+83
| import inspect | ||
| import bouquet.io.imas as imas | ||
| src = inspect.getsource(imas) | ||
| assert "impurity_charge_with_fast_ions(ne, ni, Zeff, z_fast)" in src | ||
| assert "effective_impurity_charge(ne_th" not in src |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
effective_impurity_chargetreatedne - nias impurity charge. With a beam population that difference also contains the fast ions, sonzwas inflated andZ_impdrifted off the real impurity charge.The fix, in two parts
1. Subtract the fast-ion charge and renormalize Zeff. Only
ne_th = ne − Σ_s Z_s n_s^fastis neutralized by thermal ions. The IDS Zeff is the thermal-species numerator over the fullne, so the inversion must also useZeff · ne / ne_th; the subtraction alone recovers only half the bias (adversarial review, verified numerically: raw 1.95 → half-applied 3.19 → true 6.00 for C6 at 25 % fast fraction). The inversion lives inphysics.impurity_charge_with_fast_ions, wired into the IMAS reader. The Zeff consumed by the bootstrap and the forward solve deliberately stays the full-neone — the renormalization is local to theZ_imp/p_impderivation, so the wider blast radius that concern implied never materializes.2. Run every impurity consumer on the thermal
ne.Baselinegainsz_fast(fixed across draws, likep_fast), plumbed throughgenerate_bouquetintoperturb_kinetic_equilibrium, so run.py's forward-solve assemblies, the per-draw impurity pressure, the baseline reference assembly, and the Zeff-primary channel all usene − z_fastwith the same thermal-derivedZ_impthe reader stores. Previously the reader derivedZ_imp/p_impon thermalnewhile the solves and draws used the fullne— a σ=0 pressure skew ofe·z_fast·ti/Z_impbetween reader and solver. The Zeff-primaryniderivation uses the exact thermal inverseni = (Z·ne_th − Zeff·ne)/(Z−1)and its draw clip generalizes to[ne_th/ne, Z_imp·ne_th/ne]; both reduce exactly to the old forms whenz_fastis absent, so the reconstruction path is bit-untouched.Motivation
On a DIII-D discharge with 6.3 MW NBI the fast fraction reaches 28 % through the L→H transition, and the pressure-completeness gap tracked it closely with the first (subtraction-only) form of the fix:
The 2% guard was rejecting three of those slices outright. The residual gaps in that table are the signature of the half-applied
Z_imp(p_impstill ~1.9× inflated); the renormalization is expected to close most of them.No behaviour change when no ion carries
density_fast.Tests
tests/test_fast_ion_dilution.py(10 tests): exact recovery of the true impurity charge through the shipped helper, the half-applied form pinned as a regression marker, the fast-ion-free reduction, the pathologicalz_fast ≥ nedegradation, the imas.py wiring, exactni/nzrecovery through the thermal inversion, theBaseline.z_fastfield, and source-level pins for all four consumer sites. Fast suite: 399 passed.🤖 Generated with Claude Code