JWL, and state-dependent EOS through Wood's law, characteristic BCs, probes and post-process output - #1810
Conversation
Any Mie-Gruneisen EOS p = p_ref + rho Gamma_G (e - e_ref) is MFC's rho e = Gamma p + Pi with Gamma = 1/Gamma_G and Pi = rho e_ref - p_ref/Gamma_G, so the existing operators need only the coefficients and dPi/drho. s_eos_coefficients is the single dispatch: the Mie-Gruneisen case supplies a linear-Hugoniot reference curve (u_s = c0 + s u_p, linear release) and one shared conversion produces Gamma, Pi and dPi/drho; a second family is one more case. Stiffened and ideal gas return the constants resolved at init, and no caller invokes the routine yet, so every existing answer is bit-identical - 27-case gate. Parameters mg_rho0, mg_c0, mg_s, mg_G0 follow the flat per-fluid style; the validator requires all four under mie_gruneisen, forbids them otherwise, and forbids qv there because e_ref carries the formation energy. Thirty pytest checks pin the maths against finite differences and a numerically integrated isentrope.
A Mie-Gruneisen fluid has neither gamma, pi_inf nor qv, and an ideal gas has no pi_inf: the validator now refuses them, and the suite's base fluid and every example drop the dead values. The ideal-gas half duplicates MFlowCode#1808 and drops out once it merges.
The mixture loop calls s_eos_coefficients per phase when any fluid's EOS is state dependent; the stiffened-gas arithmetic is untouched (27/27 bit-identical on the gate). The sound speed takes the partial densities as an optional argument and mixes frozen per-phase moduli, with the derivative term applied at that one site. Callers that cannot supply them are refused by the validator: 5-equation, HLL/HLLC/LF only, no bubbles, hypoelasticity, IBM, IGR, relativity, MHD, chemistry, acoustic source, probes, post-process c, characteristic BCs or Wood's law. Validation: symmetric-impact shock speed and plateau density match the Hugoniot to 1.4e-3, and a small pulse travels at the analytic sound speed to 1.9e-4 (the missing derivative term would give 24%).
gamma = gammas(1) and friends were assigned and never read; lit_gamma was assigned in both output routines and read in neither.
…eck passes from worktrees (cherry picked from commit 5ad5c85)
…ristic BCs, probes and post-process c
Isentropic release of a JWL fluid lands on the closed-form isentrope to 7.5e-5 at N = 200 and 9e-7 at N = 800; goldens for JWL beside an ideal gas, Mie-Gruneisen under Wood's law, and Mie-Gruneisen with characteristic walls. The acoustic pulse now comes from a rectangle patch and is tracked by its centroid: an analytic IC is compiled in and every distinct one costs the suite a full rebuild.
There was a problem hiding this comment.
Warning
Copilot couldn't run its full agentic review because it didn't start before the timeout. Make sure your repository has a runner available, or add a copilot-code-review.yml file specifying one with the runs-on attribute. See the docs for more details.
Pull request overview
This PR adds state-dependent EOS support (Mie–Grüneisen and JWL) and wires it through sound-speed / characteristic-BC / Wood’s-law paths, plus introduces validation and regression tests to lock in the new math and selector semantics.
Changes:
- Add Mie–Grüneisen and JWL EOS parameters, selectors, and coefficient evaluation (including derivative terms) across solver, post-process, and toolchain.
- Extend sound-speed computation to accept per-phase partial densities for state-dependent mixing (including Wood’s law).
- Add new unit + convergence tests and update example cases/validator behavior (notably: ideal gas no longer accepts
pi_inf=0).
Reviewed changes
Copilot reviewed 157 out of 157 changed files in this pull request and generated 1 comment.
Show a summary per file
| File | Description |
|---|---|
| toolchain/mfc/test_eos_mie_gruneisen.py | New manufactured-math tests for MG/JWL curve mapping, derivatives, and sound speed. |
| toolchain/mfc/test_case_validator.py | Adds validator tests for MG/JWL parameter ownership/constraints and ideal-gas stiffness rules. |
| toolchain/mfc/test/test_coverage_unit.py | Refactors git-env scrubbing helper to use shared implementation. |
| toolchain/mfc/test/coverage.py | Ensures git subprocess calls ignore inherited GIT_* vars (hook-safe). |
| toolchain/mfc/test/convergence.py | Adds convergence harness runners for MG wave speed, MG Hugoniot, and JWL release. |
| toolchain/mfc/test/cases.py | Registers new convergence cases; extends golden coverage and EOS selection cases. |
| toolchain/mfc/test/case.py | Updates post-process base config to stop setting pi_inf for ideal gas. |
| toolchain/mfc/params_tests/test_eos_selector.py | Updates selector tests to enforce “ideal gas owns no pi_inf”. |
| toolchain/mfc/params/definitions.py | Registers new EOS IDs and MG/JWL namelist parameters. |
| toolchain/mfc/eos.py | New Python-side EOS utilities used by tests and convergence harness. |
| toolchain/mfc/case_validator.py | Implements family-generic EOS selector validation + MG/JWL parameter rules and feature fences. |
| .github/scripts/check_coverage_map_health.py | Switches to shared _git() wrapper (hook-safe execution). |
| tests/ED75D01D/golden-metadata.txt | New/updated golden metadata for new EOS/paths coverage. |
| tests/A6846AD4/golden-metadata.txt | New/updated golden metadata for new EOS/paths coverage. |
| tests/5AC2F65D/golden-metadata.txt | New/updated golden metadata for new EOS/paths coverage. |
| tests/471270BB/golden-metadata.txt | New/updated golden metadata for new EOS/paths coverage. |
| src/simulation/m_time_steppers.fpp | Threads alpha_rho into state evaluation and sound-speed computation for dt. |
| src/simulation/m_sim_helpers.fpp | Extends s_compute_cell_state to output alpha_rho. |
| src/simulation/m_riemann_solver_lf.fpp | Passes alpha_rho_{L,R} into sound-speed computations. |
| src/simulation/m_riemann_solver_hllc.fpp | Passes alpha_rho_{L,R} into sound-speed computations (incl. avg-state). |
| src/simulation/m_riemann_solver_hll.fpp | Passes alpha_rho_{L,R} into sound-speed computations (incl. avg-state). |
| src/simulation/m_global_parameters.fpp | Initializes MG/JWL fluid parameter fields. |
| src/simulation/m_data_output.fpp | Threads alpha_rho into stability criteria and sound-speed computation. |
| src/simulation/m_cbc.fpp | Threads alpha_rho into sound speed + mixture-coefficient time derivatives for CBC. |
| src/pre_process/m_global_parameters.fpp | Initializes MG/JWL fluid parameter fields (pre-process). |
| src/pre_process/m_data_output.fpp | Removes unused isentrope exponent local and related computations. |
| src/post_process/m_start_up.fpp | Adds alpha_rho for post-process c output correctness under state-dependent EOS. |
| src/post_process/m_global_parameters.fpp | Initializes MG/JWL fluid parameter fields (post-process). |
| src/post_process/m_data_output.fpp | Threads alpha_rho into post-process speed-of-sound. |
| src/common/m_variables_conversion.fpp | Adds EOS parameter storage, s_eos_coefficients, state-dependent mixture coeffs, and sound-speed mixing changes. |
| src/common/m_global_parameters_common.fpp | Adds global arrays for EOS selectors + MG/JWL parameters + any_state_dependent_eos. |
| src/common/m_derived_types.fpp | Extends fluid_pp with MG/JWL parameters. |
| src/common/m_constants.fpp | Adds EOS IDs for MG and JWL. |
| examples/1D_mg_impact/case.py | New MG impact example for Hugoniot validation. |
| examples/1D_mg_acoustic/case.py | New MG acoustic pulse example for sound-speed validation. |
| examples/1D_jwl_release/case.py | New JWL release example for isentrope validation. |
| examples/*/case.py (many) | Removes pi_inf=0 from ideal-gas examples to satisfy new selector rules. |
Suppressed comments (8)
src/common/m_variables_conversion.fpp:1
alpha_rhois an OPTIONAL dummy argument ins_compute_speed_of_sound_avg, but it’s passed tos_compute_speed_of_soundunconditionally. In Fortran, an absent optional dummy may not be used as an actual argument, so this is undefined behavior and can crash/miscompile. Fix by guarding these calls withif (present(alpha_rho)) then ... call ...(..., alpha_rho) else ... call ...(...) end if, or by splitting into two call sites that do/don’t pass the argument.
src/common/m_variables_conversion.fpp:1alpha_rhois an OPTIONAL dummy argument ins_compute_speed_of_sound_avg, but it’s passed tos_compute_speed_of_soundunconditionally. In Fortran, an absent optional dummy may not be used as an actual argument, so this is undefined behavior and can crash/miscompile. Fix by guarding these calls withif (present(alpha_rho)) then ... call ...(..., alpha_rho) else ... call ...(...) end if, or by splitting into two call sites that do/don’t pass the argument.
src/common/m_variables_conversion.fpp:1alpha_rhois an OPTIONAL dummy argument ins_compute_speed_of_sound_avg, but it’s passed tos_compute_speed_of_soundunconditionally. In Fortran, an absent optional dummy may not be used as an actual argument, so this is undefined behavior and can crash/miscompile. Fix by guarding these calls withif (present(alpha_rho)) then ... call ...(..., alpha_rho) else ... call ...(...) end if, or by splitting into two call sites that do/don’t pass the argument.
toolchain/mfc/test/convergence.py:1- The report table includes a
"cells"column, but the row output always prints""for that column, which makes the output misleading and removes useful diagnostics (e.g., how many released cells were actually checked). Either populate this column (e.g., withsel.sum()per resolution) or remove the column entirely from the header/widths/row formatting.
toolchain/mfc/case_validator.py:1 - This turns
pi_inf=0.0for ideal gas from “accepted” into a hard error, which is a breaking change for existing input decks (many users setpi_inf=0explicitly). If the intent is “ideal gas ignores stiffness,” consider acceptingpi_inf=0with a warning (and still rejecting nonzero), or allowing it for one release cycle with a deprecation path, to preserve backwards compatibility while guiding users to stop setting it.
toolchain/mfc/case_validator.py:1 - The reverse mapping
{v: n for n, v in eos_names.items()}is rebuilt multiple times inside a per-fluid loop. Precompute it once (e.g.,eos_value_to_name = {v: n for n, v in eos_names.items()}) before the loop and reuse it for both the prohibit message andnamelookup—this reduces repeated work and makes the f-string much easier to read.
toolchain/mfc/case_validator.py:1 - The reverse mapping
{v: n for n, v in eos_names.items()}is rebuilt multiple times inside a per-fluid loop. Precompute it once (e.g.,eos_value_to_name = {v: n for n, v in eos_names.items()}) before the loop and reuse it for both the prohibit message andnamelookup—this reduces repeated work and makes the f-string much easier to read.
toolchain/mfc/test_eos_mie_gruneisen.py:1 - The final assertion relies on the subtle fact that
pis last computed at the pre-decrement density, so it compares againstjwl_isentrope(rho + dr, ...). This is easy to misread as an off-by-one bug. Consider storingrho_prevbeforerho -= dr, or recomputingponce after the loop at the finalrho, so the comparison is directlyjwl_isentrope(rho, ...)without the+ dradjustment.
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
|
|
||
| sys.path.insert(0, str(Path(__file__).resolve().parents[2] / "toolchain")) | ||
| from mfc.test.coverage import COVERAGE_MAP_PATH, load_map, map_health # noqa: E402 | ||
| from mfc.test.coverage import COVERAGE_MAP_PATH, _git, load_map, map_health # noqa: E402 |
Lines of Code
|
Stacked on #1805 and #1809 (and overlapping #1808); review the last two commits. Fourth step of #1638.
JWL. One more
caseins_eos_coefficients:p_ref = A e^{-R1 V} + B e^{-R2 V},V = rho0/rho, withGamma_G = omega. The curve is itself an isentrope, sode_ref/drho = p_ref/rho^2and the derivative pair is two lines; the dispatch tail now works inrho-derivatives so the conversion is family-agnostic. Six parameters (jwl_a, jwl_b, jwl_r1, jwl_r2, jwl_omega, jwl_rho0), synthetic values only. The validator's selector is family-generic: each family owns its parameters, a state-dependent fluid has nogamma,pi_inforqv, andR1 > R2 > 0.Fences lifted. Wood's law is an N-fluid loop over per-phase moduli (bit-identical for two stiffened gases); the characteristic-BC coefficient derivative applies the per-phase chain rule
d(alpha Pi)/dt = Pi dalpha/dt + Pi'(dalpha_rho/dt - rho dalpha/dt)(the alpha cancels, no division); probes and the post-processcoutput receive the partial densities. Still refused: bubbles, hypoelasticity, IBM, IGR, relativity, MHD, chemistry, acoustic sources, and the 6-equation model, whose relaxation step is written on the closed-form stiffened-gas isentrope and needs the phasic isentrope generalized first (its own PR).Verification.
1D_jwl_release): every released cell sits on the closed-form JWL isentrope through the initial state to 7.5e-5 at N = 200, 9e-7 at N = 800 (bound 1e-3).