he-burn approximate rates need denominator guard

[BenWibking]

Approximate effective-rate helpers divide by (r_pg + r_pa) without denominator floor

Summary

Multiple approximate-rate helper functions in he-burn reaclib_rates.H files compute:

dd = 1.0_rt / (r_pg + r_pa)

with no protection against near-zero or zero denominator. Under extreme conditions (very small screened rates), this can produce inf/nan and destabilize RHS/Jacobian evaluation.

Severity

Medium

Affected Code

• networks/he-burn/*/reaclib_rates.H
  • families of functions named *_approx using dd = 1.0_rt / (r_pg + r_pa)
  • examples:
    • networks/he-burn/ase/reaclib_rates.H (e.g. around lines 1934, 1953, ...)
    • networks/he-burn/cno-he-burn-33a/reaclib_rates.H (e.g. around lines 3824, 3843, ...)

Why This Is a Bug

Even if uncommon, exact/underflow-zero denominator can occur in screened-rate arithmetic, and current code has no numerical guard. A single inf/nan can contaminate burn state and Jacobian entries.

Proposed Patch

Use a denominator floor (and derivative-consistent handling) in all affected helpers:

-    amrex::Real dd = 1.0_rt / (r_pg + r_pa);
+    constexpr amrex::Real denom_floor = 1.0e-300_rt;
+    const amrex::Real denom = amrex::max(r_pg + r_pa, denom_floor);
+    amrex::Real dd = 1.0_rt / denom;

For derivative expressions containing dd * dd, reuse denom consistently.

Validation

• Add a stress test at low-temperature/high-screening conditions where both channel rates become tiny.
• Verify no inf/nan in rate or drate_dT and no regressions in standard burn trajectories.

[zingale]

it looks like we have this sort of protection in aprox13. Maybe we should add it to pynucastro