Centralize the equation-of-state expressions in the Riemann solvers - #1762
Open
sbryngelson wants to merge 7 commits into
Open
Centralize the equation-of-state expressions in the Riemann solvers#1762sbryngelson wants to merge 7 commits into
sbryngelson wants to merge 7 commits into
Conversation
s_compute_speed_of_sound took H, |u|^2 and qv and then undid them: for a real state H = ((Gamma+1)p + Pi + qv)/rho + |u|^2/2, so c^2 = (H - |u|^2/2 - qv/rho)/Gamma reduces to ((Gamma+1)p + Pi)/(Gamma rho). The three arguments cancel. Callers therefore no longer supply them, and the invisible 'H must include qv' contract cannot be stated. That contract is what produced MFlowCode#1707: five sites open-coded H and three dropped the qv this routine went on to subtract. Those three are fixed here by construction. An average of two states is not a state - its enthalpy is a free input - so the four interface-averaged sites use s_compute_speed_of_sound_avg, whose arithmetic is unchanged. Falls out as dead: qv_sf (a full-domain array whose only reader was one of the open-coded enthalpies), the H dummy threaded through s_save_data and p_main, and a per-cell qv accumulation in post-process. Claude-Session: https://claude.ai/code/session_011BUQDw64EtzvDzTvtWTNj4
|
Claude Code Review Head SHA: 9974b68 Files changed:
Findings:
|
sbryngelson
marked this pull request as ready for review
August 25, 2026 15:56
sbryngelson
requested
a lite review from Copilot
and removed request for
Copilot
August 25, 2026 15:56
Contributor
There was a problem hiding this comment.
Pull request overview
Derives real-state sound speed directly from thermodynamic state while preserving enthalpy-based calculations for interface averages.
Changes:
- Simplifies real-state sound-speed calls across simulation and post-processing.
- Adds a dedicated interface-average sound-speed routine.
- Removes obsolete enthalpy threading and
qv_sfstorage.
Reviewed changes
Copilot reviewed 12 out of 12 changed files in this pull request and generated 1 comment.
Show a summary per file
| File | Description |
|---|---|
src/common/m_variables_conversion.fpp |
Implements split sound-speed APIs and removes qv_sf. |
src/simulation/m_time_steppers.fpp |
Uses state-derived sound speed. |
src/simulation/m_riemann_solver_lf.fpp |
Simplifies LF sound-speed calls. |
src/simulation/m_riemann_solver_hypo_hlld.fpp |
Simplifies hypoelastic HLLD calls. |
src/simulation/m_riemann_solver_hlld.fpp |
Simplifies HLLD calls. |
src/simulation/m_riemann_solver_hllc.fpp |
Separates real and averaged-state calculations. |
src/simulation/m_riemann_solver_hll.fpp |
Separates real and averaged-state calculations. |
src/simulation/m_data_output.fpp |
Corrects diagnostic sound-speed computation. |
src/simulation/m_cbc.fpp |
Removes unnecessary enthalpy calculation. |
src/post_process/p_main.fpp |
Removes obsolete enthalpy argument. |
src/post_process/m_start_up.fpp |
Simplifies sound-speed output generation. |
src/post_process/m_data_output.fpp |
Removes redundant qv and enthalpy accumulation. |
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
| ! Compute mixture sound Speed | ||
| call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, ((gamma + 1._wp)*pres + pi_inf)/rho, alpha, 0._wp, & | ||
| & 0._wp, c, qv) | ||
| call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c) |
Codecov Report❌ Patch coverage is Additional details and impacted files@@ Coverage Diff @@
## master #1762 +/- ##
==========================================
+ Coverage 61.67% 61.71% +0.04%
==========================================
Files 84 84
Lines 21619 21596 -23
Branches 3196 3197 +1
==========================================
- Hits 13334 13329 -5
+ Misses 6093 6074 -19
- Partials 2192 2193 +1 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
s_accumulate_mixture_properties and the else branch of s_convert_species_to_mixture_variables_kernel contained the same four-line accumulation loop. Move the routine from m_riemann_state.fpp into m_variables_conversion.fpp and have the kernel call it, below the mpp_lim clipping so it still sees the clipped volume fractions. The two routines are not merged. s_accumulate_mixture_properties is a parameterised subset accumulator - callers pass num_fluids or num_fluids - 1, and raw or limited volume fractions - while the kernel always spans num_fluids, clips in place, special-cases num_fluids == 1 with bubbles_euler, and optionally emits Re_K and G_K. Merging would need a clipping flag and optional dummies on a [seq] device routine, which is not portable across the offload backends. This is the single place per-fluid EOS dispatch will enter the stiffened-gas mixture path. Claude-Session: https://claude.ai/code/session_011BUQDw64EtzvDzTvtWTNj4
Adds s_compute_energy to m_variables_conversion: the stiffened-gas total energy, thermodynamic terms only. Magnetic and elastic energies stay at their call sites because they are not equation-of-state terms, and the chemistry and relativistic forms are left open-coded because they are not this relation. Converts the ten plain sites in hll, hllc, hypo_hlld and lf. The four MHD sites and the two hllc bubbles sites that omit qv are deliberately left for follow-up commits; the latter changes behaviour. Also deletes a dead energy assignment in m_cbc: E is read only inside the chemistry branch of the flux update, so the non-chemistry assignment was never used. Measured on MI210 / amdflang: all 60 Riemann kernel resource profiles identical (scratch, VGPR, AGPR); 5eq_rk3_weno3_hllc 2.290 / 2.385 / 2.407 against a base of 2.374 / 2.422 / 2.361. Claude-Session: https://claude.ai/code/session_011BUQDw64EtzvDzTvtWTNj4
hll, hllc, hlld, hypo_hlld and lf each derived the stiffened-gas mixture coefficients their own way: two called a shared accumulator, two hand-rolled the loop inline, and hllc carried a three-way branch. They now all call s_compute_mixture_coefficients, which owns the rule including its one special case. That special case is required, not incidental. Under bubbles_euler with num_fluids == 1 the sole advection slot aliases the void fraction (eqn_idx%alf == eqn_idx%adv%end), so alpha is not a composition there and the mixture rule does not apply; the coefficients are the liquid's. Clipping stays with the callers because it is genuinely not uniform - hll clips local arrays while hllc clips the shared reconstruction buffers in place - and it can never coincide with the special case, since case_validator prohibits mpp_lim with num_fluids == 1. Behaviour fix, not a refactor: hll and lf previously lacked the special case and computed gamma = alf*gammas(1) under bubbles_euler, which is physically meaningless. No test or example exercises that combination, but it runs today, so results change for anyone using it. hllc bubbles at num_fluids == 2 likewise moves from pure-liquid to standard accumulation; also unexercised. s_accumulate_mixture_properties is merged in and deleted: after hllc's three-way branch collapsed, its num_fluids - 1 call site was gone and it had exactly one caller. The gpuParallelization.md example that used it is updated. s_compute_energy now takes composition rather than coefficients, deriving them through the same rule, so no call site holds equation-of-state knowledge. Measured on MI210 / amdflang: hypo_hlld drops 24 VGPR and 24 AGPR because the caller no longer keeps gamma, pi_inf and qv live across the call; hllc and lf gain 8 VGPR at 118-142, far from any occupancy limit. 5eq_rk3_weno3_hllc 2.343 / 2.319 / 2.330 against 2.304 / 2.283 / 2.366 before. Claude-Session: https://claude.ai/code/session_011BUQDw64EtzvDzTvtWTNj4
sbryngelson
changed the title
The operator returns the thermodynamic energy; pres_mag is added at the call site, because magnetic energy is not an equation-of-state term and a second EOS backend must not have to know about it. Also trims the operator doc comments to the facts a reader cannot derive from the code.
The bubbles_euler branch computed E without qv while every other energy site in the file includes it, and s_compute_pressure's bubbles branch subtracts qv when inverting. Forward and inverse transforms disagreed, so E was inconsistent with the pressure it came from whenever fluid_pp(1)%qv is non-zero. E feeds the energy flux directly, so the error was not confined to diagnostics. Routing the site through s_compute_energy restores qv. No test or example combines bubbles_euler with a non-zero qv, so the suite neither demonstrates the bug nor the fix; the argument is the forward/inverse inconsistency above. Measured on MI210 / amdflang: HLLC drops 226 to 214 VGPR; 5eq_rk3_weno3_hllc 2.339 / 2.343 / 2.283.
sbryngelson
changed the title
…e shared rule m_viscous carried the same three-way branch four times over and m_acoustic_src once more, under Tait naming: B_tait is sum(alpha*pi_infs) and small_gamma is sum(alpha*gammas). Both now call s_compute_mixture_coefficients, leaving one implementation of the mixture rule in the codebase. m_acoustic_src also had a latent double-accumulation: with bubbles_euler, mpp_lim and num_fluids > 2 both of its blocks ran, adding onto variables the second block does not reset. Unreachable today, since bubbles cases use num_fluids = 1, and gone now. As with the solvers, hll and lf aside, this gives these paths the bubbles special case they lacked. No test or example combines bubbles_euler with num_fluids > 1, so the num_fluids > 2 arms were unexercised.
Lines of Code
|
sbryngelson
changed the title
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.
Gives the mixture coefficients, the energy and the sound speed one implementation each, so that adding a second EOS (#1638) is a change inside one module rather than an edit to every solver.
Net -202 lines.
The core of it
s_compute_speed_of_soundused to takeH,|u|^2andqvand then undo them. For a real stateH = ((Gamma+1)p + Pi + qv)/rho + |u|^2/2, soc^2 = (H - |u|^2/2 - qv/rho)/Gammareduces to((Gamma+1)p + Pi)/(Gamma rho). All three cancel. Removing them makes #1707 unrepresentable: thereis no
Hargument to get wrong.The same idea applied to the other two families:
s_compute_energytakes composition, not coefficients, so no solver holds any EOS knowledge.Magnetic and elastic energy stay at the call sites - they are not EOS terms.
Behaviour changes (not just refactoring)
hllandlflacked thebubbles_eulerspecial case and computedgamma = alf*gammas(1),with
alfa void fraction. Underbubbles_eulerwithnum_fluids == 1the sole advection slotaliases
alf, so the mixture rule does not apply there.qvwhiles_compute_pressure's bubbles branch subtracts it.Forward and inverse disagreed;
Efeeds the energy flux directly.Neither combination is exercised by any test or example, but both run today.
Verification
MI210, AFAR amdflang, OpenMP offload. Every commit measured; none slower.
5eq_rk3_weno3_hllchypo_hlld-24 VGPR / -24 AGPRSuite: 686 passed / 0 failed on the first two commits, no goldens regenerated. The later commits
have not been run locally - CI is their first correctness gate.
Scope
Closes three of roughly ten EOS expression families. Thirteen files and 96 parameter references
remain, the largest being
m_phase_change(39, and the only consumer ofcvs/qvps),m_igr(20)and
m_pressure_relaxation(14). Plan for the rest:docs/superpowers/plans/2026-08-25-finish-eos-centralisation.md.Two constraints found here bound the JWL follow-up more than any interface question:
m_pressure_relaxation's Newton solve needs a closed-form invertible isentrope, which JWL lacks;and phase change needs entropy and temperature, which JWL does not define without a caloric
extension.
Related
Supersedes #1714 (same goal; its
type(eos_state)is the +20%, because a derived-type dummy on a[seq]device routine forces a memory ABI). Closes #1707. Partially addresses #1708. Unblocks #1700and #1638.
Caveat on the measurements: #1759 documents amdflang regenerating the whole device image at link
time, so kernel counts changing can degrade untouched kernels. Counts are constant on both sides of
every comparison above (471 kernels), so it is not in play - but these builds do not carry that
workaround.