Skip to content

Make the NEMO total-energy under-relaxation check reachable - #2884

Open
bellonarts wants to merge 2 commits into
su2code:developfrom
babybluechips:fix_nemo_energy_underrelaxation
Open

Make the NEMO total-energy under-relaxation check reachable#2884
bellonarts wants to merge 2 commits into
su2code:developfrom
babybluechips:fix_nemo_energy_underrelaxation

Conversation

@bellonarts

Copy link
Copy Markdown

Proposed Changes

CNEMOEulerSolver::ComputeUnderRelaxationFactor limits the implicit update so that no monitored variable changes by more than MAX_UPDATE_FLOW in one nonlinear iteration. The total-energy check

/*--- Energy ---*/
if (iVar == (nVar-2)) { ... }

is nested inside the if (iVar < nSpecies) branch. In the NEMO conservative layout total energy sits at nVar - 2 = nSpecies + nDim, which can never satisfy iVar < nSpecies, so the energy clamp is unreachable and the solver never under-relaxes the total-energy update. This PR closes the species block before the energy check and evaluates the check independently. No new logic is added; the existing check simply runs. This is a correctness repair of an unreachable branch.

Behavior change: aggressive implicit updates to total energy are now bounded by the same fraction as the species updates, so the implicit update changes wherever that bound is active. On the public viscous cone (TestCases/nonequilibrium/visc_wedge/axi_visccone.cfg) the iteration-10 snapshot moves visibly. Against the vector stored in serial_regression.py, rms[Rho_0] goes from -5.215230 to -5.113745, rms[RhoE] from 1.262782 to 1.473443 and the total heat flux from 32641 to 42641; the run completes normally. The regression vectors of every affected NEMO case are rebaselined in this PR:

  • serial_regression.py: invwedge, visc_cone
  • parallel_regression.py: invwedge_a, invwedge_msw, invwedge_roe, invwedge_lax, invwedge_ss_inlet, visc_cone, super_cat, ion_gy

Every other NEMO case (thermalbath, thermalbath_frozen, invwedge_ap2, invwedge_ausm_m, partial_cat, ionized) prints the same digits with and without this hunk in both modes. The new x86-64 vectors come from an independent Linux build of develop at 07aa46b1 with this hunk applied, and the numbers quoted above are the values that build prints.

On the tolerance of these vectors: TestCases/TestCase.py sets self.tol = 0.0 by default and no NEMO case in either regression script overrides it, so the harness compares the printed digits exactly. A stock build of 07aa46b1 on the same Linux machine reproduces the stored x86 vectors exactly on invwedge_a (MPI), invwedge_lax, invwedge_ss_inlet, super_cat and ion_gy, but not on every case:

  • visc_cone: 4e-6 in the last digit of several entries, serial and MPI
  • invwedge (serial): 1e-6 in one entry (-18.573334 stored, -18.573333 on the machine)
  • invwedge_msw (MPI): 6.6e-4
  • invwedge_roe (MPI): 4.7e-3 (stock serial and MPI vectors for this case also differ from each other by more than the other wedges do)

Those differences are machine-to-CI noise on a stock build, so the digits this PR writes for those four cases will likely need a touch-up from the CI output, and invwedge_roe is the most sensitive of them. The final x86 digits are for CI to confirm. The test_vals_aarch64 vectors are not touched; the arm64 workflow only runs on push to develop, so those need a maintainer.

A unit test, UnitTests/SU2_CFD/solvers/CNEMOEulerSolver_tests.cpp, builds a NEMO Euler solver (AIR-5, SU2_NONEQ) on the unit box, writes a single linear-system update that exceeds MAX_UPDATE_FLOW by a factor of four for one variable at one point, and calls ComputeUnderRelaxationFactor. It checks that a species update is limited (the control, which also passes on stock), that a total-energy update is limited to the same factor (fails on stock, where the factor stays at one), that neighbouring points keep a factor of one, and that an in-bounds energy update is not limited. The test constructs CNEMOEulerSolver directly rather than through UnitQuadTestCase::InitSolver (which goes through CSolverFactory) so that it depends on the one class under test.

Scope, stated plainly: on our own Mach 7.95 viscous cone a bisection of a larger stabilization patch set showed that this hunk alone does not rescue the run (the energy-only arm diverged later than stock but still diverged). The PR therefore claims correctness of an unreachable branch and the matching residual signature (stock's rms[RhoE] runaway), not robustness.

Related Work

Related to the NEMO solver (SU2_NEMO). No open issue. This is one of a series of NEMO PRs from the same work (freestream mass-fraction tolerance, history MAX_/BGS_ indexing, consistent NEMO viscous Jacobian, coarse-grid viscous Jacobian edge length). The source hunks are independent and each applies to develop on its own. The regression-script edits are not independent: the two viscous Jacobian PRs also rebaseline the visc_cone, super_cat and ion_gy test_vals lines in parallel_regression.py (and visc_cone in serial_regression.py), so whichever of those PRs merges after this one will conflict on those lines, and the vectors for the combined code differ from any single-PR set. Suggested handling: merge these PRs one at a time, rebase the next one on the new develop, and take the digits for the overlapping cases from that PR's CI run rather than from any of the individual PR bodies.

PR Checklist

  • I am submitting my contribution to the develop branch.
  • My contribution generates no new compiler warnings (try with --warnlevel=3 when using meson).
  • My contribution is commented and consistent with SU2 style (https://su2code.github.io/docs_v7/Style-Guide/).
  • I used the pre-commit hook to prevent dirty commits and used pre-commit run --all to format old commits.
  • I have added a test case that demonstrates my contribution, if necessary.
  • I have updated appropriate documentation (Tutorials, Docs Page, config_template.cpp), if necessary.

bellonarts and others added 2 commits September 3, 2026 13:25
In CNEMOEulerSolver::ComputeUnderRelaxationFactor the total-energy check
(iVar == nVar-2) was nested inside the species block (iVar < nSpecies).
Total energy sits at nSpecies+nDim, so the check could never run and the
energy update was never under-relaxed. Close the species block first and
evaluate the energy check on its own; no new logic is added.

This changes the implicit update wherever the energy bound is active, so
the affected NEMO regression vectors are rebaselined (serial: invwedge,
visc_cone; parallel: invwedge_a, invwedge_msw, invwedge_roe, invwedge_lax,
invwedge_ss_inlet, visc_cone, super_cat, ion_gy). The aarch64 vectors are
left for a maintainer. A unit test builds a NEMO Euler solver on the unit
box and checks that an over-large energy update is limited (it stays
unlimited on the previous code) alongside a species control.

Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
Update the branch-specific serial and two-rank aarch64 expectations from GitHub Actions run 33998938091 at exact pre-update head a495681. Artifact digest: sha256:4c9d814dd9a373592fce4c9f33f74e565e26016ee855b2df8660bbcc99e4d560. No x86 expectation or solver source changes.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant