Fix missing qv term in the Lagrange bubble initial pressure - #1709
Fix missing qv term in the Lagrange bubble initial pressure#1709sbryngelson wants to merge 2 commits into
Conversation
There was a problem hiding this comment.
Pull request overview
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.
Fixes incorrect initial Lagrangian bubble pressure when a nonzero heat-of-formation term (qv) is configured, and adds a targeted regression test/golden to prevent recurrence.
Changes:
- Include the missing
qvterm in the stiffened-gas pressure inversion used to seed Lagrangian bubble pressure. - Add a new lag-bubble test configuration with
fluid_pp(1)%qv = 0.01to exercise the previously untested path. - Add a new golden (
43EA05B4) to validate the behavior.
Reviewed changes
Copilot reviewed 3 out of 4 changed files in this pull request and generated 2 comments.
| File | Description |
|---|---|
src/simulation/m_bubbles_EL.fpp |
Fixes pressure inversion by subtracting qv in the bubble initialization path. |
toolchain/mfc/test/cases.py |
Adds a constrained regression case that sets qv nonzero for lag bubbles. |
tests/43EA05B4/golden-metadata.txt |
Adds metadata for the new golden tied to the new regression case. |
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
|
Claude Code Review Head SHA: b3431b4 Files changed:
Findings:
|
s_add_bubbles open-coded the stiffened-gas pressure inversion and omitted the
qv (heat of formation) term:
pliq = (E - dynP - pi_inf)/gamma
The canonical inversion in s_compute_pressure (m_variables_conversion.fpp:71)
is (energy - dyn_p - pi_inf - qv)/gamma. The omission is clearly accidental:
qv is an output argument of the s_convert_to_mixture_variables call one line
above and is never read anywhere else in that scope.
pliq seeds gas_p(bub_id, 1), so every Lagrangian bubble started from a wrong
internal pressure. This is physics, not diagnostics. It also propagates to the
initial gas mass and the bubble natural frequency, which are both derived from
gas_p in the same routine.
Consistency with the running solver settles which quantity pliq should be:
s_get_pinf, the Maeda and Colonius (2018) subgrid closure that supplies the
bubble driving pressure at every later step, interpolates q_prim_vf(eqn_idx%E),
and that field is filled by s_compute_pressure with qv included. Initializing
from an inversion that omits qv left each bubble out of equilibrium with its
own driving pressure by qv/gamma.
Only s_add_bubbles needs the change; the restart path reads gas_p straight from
the restart file rather than recomputing it.
Dormant when qv = 0, the default, and every existing lag-bubble test leaves it
there, which is why this survived.
Adds a 2D one-way-coupled case with qv /= 0. The golden lands in a follow-up
commit so that it is generated from a clean tree.
Fixes #1706
Generated from a clean checkout of the preceding commit so the recorded provenance corresponds exactly to committed sources. The golden captures beta, the Lagrangian void fraction, alongside the conservative variables. beta follows the bubble radius, which is driven by the initial gas pressure, so the defect is observable even under one-way coupling where the bubbles do not feed back into the Eulerian field.
b3431b4 to
9191c9a
Compare
|
Both Copilot points were valid and are addressed in the force-push ( Label Dirty working tree. Also correct, and the recorded SHA which is the parent commit in this PR. One thing the restructure caught. After committing the fix, my original verification method silently stopped working — Worth stating explicitly because a regression test that cannot fail is worse than no test, and the false pass looked identical to a real one. Also updated the description with two things I had missed on the first pass: the fix can newly trigger the vapour-pressure and negative-gas-mass aborts for |
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## master #1709 +/- ##
=======================================
Coverage 60.77% 60.77%
=======================================
Files 83 83
Lines 20872 20872
Branches 3101 3101
=======================================
Hits 12685 12685
Misses 6121 6121
Partials 2066 2066 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
Fixes #1706.
The defect
src/simulation/m_bubbles_EL.fppopen-codes the stiffened-gas pressure inversion and drops theqv(heat of formation) term:The canonical inversion,
s_compute_pressureinm_variables_conversion.fpp:71, isThe omission looks accidental rather than a deliberately different quantity:
qvis an output argument of thes_convert_to_mixture_variablescall on the line above, and it is never read anywhere else in that scope. It is computed and discarded.Impact
pliqis too large byqv/gamma, and it seeds the initial bubble gas pressure:So every Lagrangian bubble starts from the wrong internal pressure. This is physics, not diagnostic output.
Dormant when
qv = 0, which is the default. It triggers for any case settingfluid_pp(i)%qv— the phase-change and reactive-burn configurations.Why it survived
No test exercised the combination. Exactly one golden case sets
bubbles_lagrange, and it leavesqvat zero; the cases that setqvare phase-change and reactive-burn, none of which use Lagrangian bubbles. The suite could not have caught this.Test
Adds a 2D one-way-coupled Lagrange bubble case with
fluid_pp(1)%qv = 0.01, restricted to a single configuration so it contributes one golden (F428FDC0).The golden captures
beta, the Lagrangian void fraction, alongside the conservative variables.betafollows the bubble radius, which is driven by the initial gas pressure, so the error is observable even under one-way coupling.Verification
Built with GNU 15.2.0, MPI, on macOS. The negative result is the one that matters — without it, a green suite would prove nothing here, since the pre-existing cases are all
qv = 0.No existing golden moved: the change is confined to a path only reachable with
bubbles_lagrange, and no existing lag case has nonzeroqv.Behaviour change worth flagging
The correction shifts
pliqby-qv/gamma, andgamma > 0always, so the direction follows the sign ofqv. Forqv > 0the initial bubble pressure drops, which brings two existing guards in the same routine within reach:A case with
qv > 0and Lagrangian bubbles that previously ran can now abort at startup. That is the correct outcome — it was running on an inflated bubble pressure that masked a setup the model rejects — but it converts a silently wrong answer into a hard stop, so it is a user-visible change rather than a pure numerical shift. Forqv < 0(as in the phase-change examples, which useqv = -1.167e6for fluid 1) the pressure rises and the guards become less reachable.gas_palso feeds the initial gas mass and the bubble natural frequency in the same routine, so both are corrected by the same change.Adjacent issue found while reviewing
s_add_bubblesalso assigns the module-level globalpreffrom the last-added bubble'sgas_p(line 383).prefis a documented user-settable case parameter that the four-equation EOS branch reads. This PR changes the value it gets overwritten with, but the overwrite itself is pre-existing and looks wrong independently. Filed as #1710 rather than addressed here.Context
First of three PRs against #1708. The root cause is that the EOS algebra has several independent definitions, and
s_compute_pressureis awkward enough to call — it demands anum_speciesarray and anintent(inout)temperature — that call sites open-code the expression instead. This PR fixes only the bug; deduplication follows separately so that any behavior change lands in a PR whose purpose is to change behavior.Handover notes
Context for anyone picking this up on another machine.
Branch and commits
Two commits deliberately: the golden in the second is generated from a clean checkout of the first, so
golden-metadata.txtrecords a real committed SHA rather than a dirty tree. Copilot flagged the dirty-tree provenance on the original single-commit version; that is why it is split. Preserve this if you amend — regenerate the golden only from a committed, clean tree.Environment
Verified with GNU 15.2.0 (gfortran), MPI on, no GPU, macOS arm64. A first
./mfc.sh buildbootstraps a Python venv and builds hdf5/silo; budget ~20 min.Committing runs a pre-commit hook that executes the full CI lint-gate (7 checks). Use
--no-verifyonly when you have already run it.The exact defect
src/simulation/m_bubbles_EL.fpp, ins_add_bubbles:Canonical inversion is
s_compute_pressureinsrc/common/m_variables_conversion.fpp:pres = (energy - dyn_p - pi_inf - qv)/gamma.qvis an output argument of the call one line above and is read nowhere else in that scope — computed and discarded. That is the evidence the omission is accidental rather than a deliberately different quantity.Why
pliqmust be the full thermodynamic pressureSettled by internal consistency, not by the Maeda paper (which I did not read — it is paywalled JCP;
docs/documentation/equations.mdsection 6.2 documents the model):s_get_pinfin the same file — "Compute the bubble driving pressure p_inf", the Maeda & Colonius (2018) subgrid closure — interpolatesq_prim_vf(eqn_idx%E)%sf. That slot is filled atm_variables_conversion.fppbys_compute_pressurewithqv_Kpassed. So at every step after t=0 the pressure driving the bubble includesqv. Initialising from an inversion that omits it left each bubble out of equilibrium with its own driving pressure byqv/gamma.Only
s_add_bubblesneeds the change: the restart path (s_restart_bubbles) readsgas_pstraight from the restart file rather than recomputing it.Test case: why it is shaped the way it is
Added in
toolchain/mfc/test/cases.py, insidealter_lag_bubbles, guarded tondims == 2 and couplingMethod == 1 and adap_dt == "F"so it contributes exactly one golden (F428FDC0).Two constraints that are easy to get wrong:
beta. Under one-way coupling the bubbles do not feed back into the Eulerian field, so the conservative variables alone would not move. It works becausebeta, the Lagrangian void fraction, follows bubble radius and therefore the initial gas pressure. If you restructure the case, confirmbetais still intests/F428FDC0/golden.txt.qvto the existing lag-bubble stack. That would churn every existing lag golden. A new label is what keeps it to one.The case label feeds the golden UUID —
crc32(sha1(str(trace)))attoolchain/mfc/test/case.py. Renaming the label changes the directory name. Avoid!in labels (history expansion in interactive bash);qv_nonzeroreplaced an earlierqv!=0for that reason.How the fix was verified, and a trap
Trap: once the fix is committed,
git stash push -- src/simulation/m_bubbles_EL.fpphas nothing to stash and silently leaves the fixed binary in place — the "reverted" run then reports a pass that looks identical to a real one. Revert withgit checkout master -- src/simulation/m_bubbles_EL.fppinstead, rebuild, run, thengit checkout HEAD -- <file>to restore.This matters because no pre-existing golden covers the path: only one case sets
bubbles_lagrangeand it hasqv = 0, while theqv /= 0cases are phase-change and reactive-burn without Lagrangian bubbles. A green suite proves nothing here without the negative check.Related work
prefoverwrite in the same routine. Already resolved: Remove the four-equation model (model_eqns = 4) #1713 (merged) deleted the global.eos_staterefactor; covers the enthalpy half of EOS algebra is duplicated across call sites; needs a single source of truth #1708.m_pressure_relaxation.fpp,m_data_output.fpp(x2), and this file all open-code the same inversion.Purpose in the series
A correctness fix, and the first concrete evidence for #1708.
The EOS algebra has nine hand-written copies of two expressions. Four are wrong. This is one of them, and the only one that affects physics rather than diagnostics — the others are output paths. It stands alone as a bug fix, but the reason it exists is structural:
s_compute_pressureis awkward enough to call (it demands anum_speciesarray and anintent(inout)temperature) that call sites open-code the one line instead, and then drift.Fixing the instance without fixing the pressure is deliberate — the deduplication is tracked separately so any behaviour change lands in a PR whose stated purpose is to change behaviour.
Working conventions and hazards (shared across this series)
Collected from the work that produced #1705, #1709, #1712, #1713, #1714, #1716. Every one of these cost real time or produced a wrong result before being caught.
Testing
A regression test that cannot fail is worse than no test. Always verify the negative: revert the fix, rebuild, confirm the case fails, restore. Two ways this silently broke here:
git stash push -- <file>has nothing to stash once the fix is committed, so the "reverted" run tests the fixed binary and reports a pass identical to a real one. Usegit checkout master -- <file>, rebuild, test, thengit checkout HEAD -- <file>.model_eqns = 3cannot detect a sound-speed defect, because the six-equation branch ofs_compute_speed_of_soundtouches neitherHnorqv.The golden packer discards data.
toolchain/mfc/packer/pack.pytreated every.datunderD/as<x> [<y> <z>] <value>and kept only the last column of each row. Probe and integral output are multi-column time series, so most columns were never compared (#1711, fixed in #1712). Before asserting that a golden covers something, check it is actually ingolden.txt.Case labels are load-bearing. The golden UUID is
crc32(sha1(str(trace)))— the label chain determines the directory name. Renaming a label renames the golden. Avoid!in labels (history expansion in interactive bash).Local suite runs are flaky at high
-j. Non-reproducible failures appeared on several unrelated branches at-j 12–16(chemistry cases, probe cases) that passed individually and in clean reruns. Re-run before investigating.Removing parameters or features
Deregistering a parameter breaks things that are not the source tree. Removing
pref/rhoreffrom the registry broke the entire suite becauseBASE_CFGintoolchain/mfc/test/case.pyset them for every case. Also checkfp_stability.py,params_tests/mutation_tests.py, and lint fixtures that use real parameter names as examples.Grep the generated artifacts, not just the sources. A stale
TYPED_DECLSentry naming a deleted type survived removal and did not break the build only because the parameter had also left every target's namelist vars, so it was never emitted. Checkgenerated_decls.fpp,generated_constants.fpp,SIM_GPU_DECL_VARS, and the MPI broadcast generators.Dead-local tell: after removing a block, a local with exactly one remaining occurrence in its file is almost certainly its own declaration. Two occurrences often means declaration plus a
private()entry.Fortran is case-insensitive. A local
pRefshadowed the module globalprefin the hardcoded-IC files; the read site was spelledprefand looked like a reference to the global. It is not. Confirm scope before concluding a global is live.GPU
A CPU test run cannot catch a missing
private(). It is a silent device race. Audit by hand or by script when adding per-cell state.Do not match
GPU_PARALLEL_LOOPnaively —END_GPU_PARALLEL_LOOPcontains the same substring and will register as a loop start, producing false positives. Exclude it explicitly.Derived-type components cannot have runtime extents.
dimension(num_fluids)in a type fails to compile outside case-optimized builds, wherenum_fluidsis aparameter.Benchmarking
Run
./mfc.sh benchin the foreground on an idle machine. Running it in the background while a pre-commit precheck ran at-j 12produced a bogus +50% regression on a case that executes none of the changed code. Baseline noise here is ±3.5%; run the baseline twice before trusting any delta, and sanity-check that the regressing cases actually execute the modified code.GitHub mechanics
--force-with-leaseneeds an explicit SHA (--force-with-lease=<branch>:<sha>) when the ref has not been fetched in the current clone; the bare form fails with "stale info".gh run view --log-failedcan miss the real output entirely. On the Frontier jobs the failing step carried only a non-zero exit while the actual test output lived in a separatePrint Logsstep that succeeded. Fetch the full log.gh run rerun <id> --failedrefuses while the workflow is still running; retry later.file INSTALL cannot set modification time ... No such file or directory, exit 143). Check for a real error before assuming a code fault.