Skip to content

Anchored Dual-pass HLLD for Hypoelasticity (+ HLLC and interface-consistent HLL)#1414

Open
ChrisZYJ wants to merge 113 commits into
MFlowCode:masterfrom
ChrisZYJ:hypo_hlld
Open

Anchored Dual-pass HLLD for Hypoelasticity (+ HLLC and interface-consistent HLL)#1414
ChrisZYJ wants to merge 113 commits into
MFlowCode:masterfrom
ChrisZYJ:hypo_hlld

Conversation

@ChrisZYJ

@ChrisZYJ ChrisZYJ commented May 9, 2026

Copy link
Copy Markdown
Contributor

Description

Adds:

  1. Hypoelasticity: Anchored Dual-pass HLLD
  2. Hypoelasticity: HLLC
  3. Hypoelasticity: HLL option (interface-consistent)
  4. HLL option (alpha div U) so non-conservative treatment aligns with HLLC

Key Design Choices

Separate HLLD Riemann Solvers

At a glance it might be tempting to combine HLLD MHD with dual-pass hypoelasticity HLLD, but keeping them separate makes the code cleaner and much easier to maintain because:

  1. Unlike HLL or HLLC, HLLD is a class of HLLD-type solvers, with formulas and states dependent on the eigenstructure of the governing equations, so the inner states' equations are completely different for MHD vs Hypoelasticity.
  2. HLLD hypoelasticity has a newly developed dual-pass anchored form, making it different from any convenional HLLD Riemann solver. The anchored forms are necessary for the non-conservative hypoelasticity terms, which MHD does not have.
  3. MHD and Hypoelasticity deal with completely different physical regimes with different governing equations, and any changes or new physical models added in the future will not apply to both modules at once.

Riemann Source Terms

For the non-conservative terms, unlike the usual governing equations that only need div U i.e. du/dx, dv/dy, dw/dz (alpha div U, K div U, etc.), Hypoelasticity has cross terms like du/dy, so we must also pass those Riemann-consistent traces from Riemann solver to the rhs. (The old Hypoelasticity code with the HLL Riemann solver uses finite difference for non-conservative rhs, which provides enough stability given that HLL smears the interface immediately, so there wasn't a need to pass the du/dy traces before this PR. But that does not work for HLLC/HLLD for Hypoelasticity.)

Also grouped/named the condition branches (with lots of comments within the code):

Branch Face quantity read RHS formula per $\alpha_k$ K*div(u) velocity source
adv_src_alpha_iface flux_src_n(dir)%vf(j_adv) = per-fluid $\Psi_{\alpha_k}$ $u_\text{cell} \cdot \Delta\Psi_\alpha / \Delta x$ nc_iface_vel_n(dir)%vf(1)
adv_src_vel_iface flux_src_n(dir)%vf(adv\%beg) = shared $\Psi_u$ $\alpha_k \cdot \Delta\Psi_u / \Delta x$ Same flux_src_n slot (already $\Psi_u$)
adv_src_none Skipped (HLLD handles internally)

The derivations, meanings, and usage of the Riemann source variables are not straightforward. I've added some hopefully very helpful notes in misc/dev_notes for future developers (or AI agents; directing them to my notes should help them make fewer mistakes with the source terms) in terms of the understanding and derivations for the HLL/HLLC non-conservative fluxes, and their variable mapping for Riemann solvers and RHS.

Backwards Compatibiilty

  • All default behaviors preserved exactly (newly added features as options)
    • The only exception is the removal of an incorrect ad-hoc fluids-limit guard that affects only Hypoelasticity HLL
  • All existing usage of Riemann and rhs source terms are preserved. No refactor is done to keep the scope of this PR limited (any refactoring would touch most of the existing HLLC functionalities)

Type of change

  • New feature

Testing

  • All tests passed locally on CPU and Nvidia GPU, and on Frontier.

  • Smooth Eigenmode Convergence

image
  • Weak Solution Comparison (Rodriguez & Johnsen (2019) §5.3(b))
image
  • Weak Scaling on Frontier
image

Checklist

  • I added or updated tests for new behavior
  • I updated documentation if user-facing behavior changed
  • GPU results match CPU results
  • Tested on NVIDIA GPU or AMD GPU

AI code reviews

Reviews are not triggered automatically. To request a review, comment on the PR:

  • @coderabbitai review — incremental review (new changes only)
  • @coderabbitai full review — full review from scratch
  • /review — Qodo review
  • /improve — Qodo code suggestions
  • @claude full review — Claude full review (also triggers on PR open/reopen/ready)
  • Add label claude-full-review — Claude full review via label

@ChrisZYJ

ChrisZYJ commented Jul 2, 2026

Copy link
Copy Markdown
Contributor Author

And thanks @sbryngelson for the suggestions! I addressed the following cleanup items:

  • Added checker + Python validator guard for mhd + hypoelasticity + HLLD.
  • Deduped repeated hypo elastic-energy and HLL/HLLC wave-speed formulas via Fypp inline macros.
  • Fixed an HLLC 2D-axisym elastic-energy bug: HLLC now doubles only shear_indices, matching m_variables_conversion.
  • Kept the intrusive s_hypo_hlld_riemann_solver helper extraction as follow-up scope. The solve kernel has been optimized for Frontier, and prior profiling showed it is register/occupancy-sensitive, so moving large blocks behind helpers could be done separately with performance profiling.

Also fixed the spelling CI failure: docs/documentation/case.md now uses patch(es) instead of patche(s). I think
this was exposed by a newer typos version on CI.

@ChrisZYJ

ChrisZYJ commented Jul 2, 2026

Copy link
Copy Markdown
Contributor Author

@sbryngelson If this looks good after the CI rerun, would it be possible to merge this PR soon? I’m happy to continue improving the code in follow-up PRs after this lands. The branch has already gone through several upstream re-merges and CI/review cleanup rounds, so merging the core implementation now would make future cleanup work smaller and easier to review. Thank you!

@sbryngelson

Copy link
Copy Markdown
Member

I'm no longer merging messy/bloated/WIP code on the promise that it will be fixed later. I've been burned on this far too many times (perhaps not by you), and it only ratchets up, not down.

Please clean your code of smells, refactor into nice helpers, and minimize SLOC while maintaining speed if you want the PR merged. Thanks.

@ChrisZYJ

ChrisZYJ commented Jul 3, 2026

Copy link
Copy Markdown
Contributor Author

Okay. I'll refactor the solver into helpers according to your reviews, and minimize SLOC as you requested. Will post the changes here soon.

@ChrisZYJ

ChrisZYJ commented Jul 4, 2026

Copy link
Copy Markdown
Contributor Author

I've done the refactor. Here is what changed:

  • Factored the 25 repeated per-component wave-fan flux folds into small Fypp component tables driving one shared template.
  • Factored the 12 outer-wave stress star-state updates into one formula family with a table over component, sign, and elastic coefficient.
  • Used a single definition for the per-component HLL flux formula in both the degenerate-fan fallback and the ADC blend.
  • Factored s_finalize_riemann_solver_hatR and the generic finalizer into one template in m_riemann_state.fpp, so they will not diverge in the future.
  • Used a single [seq] helper f_hlld_wave_zone for the five-wave fan-zone selection, now shared by the flux fold, the NC face-velocity export, and the axisymmetric face-state pick.
  • Changed the three mutually exclusive hypo_nc_* booleans into one hypo_nc_mode enum, following the same pattern as adv_src_mode.

This removes 95 lines across the touched files. m_riemann_solver_hypo_hlld.fpp drops from 1084 to 968 lines, the main subroutine from 938 to 859, and the GPU loop body from ~820 to ~690.

What is left in the main subroutine is dominated by the per-cell star-state algebra, which appears only once, so there is no duplication left to remove. Extracting it behind an interface would take ~50 scalar arguments or regrouping the privatized locals into derived types, and either would change register allocation on the GPU kernel. Since the goal was fewer lines at unchanged speed, that block stays inline, and the helper extractions went where the interfaces are narrow.

I've taken care to ensure that most changes do not affect the generated code, and that they have minimal impact on performance. Also merged in the current master.

@ChrisZYJ

ChrisZYJ commented Jul 4, 2026

Copy link
Copy Markdown
Contributor Author

One note on CI. The last two runs each failed a single Frontier job at the Fetch Dependencies step, before anything gets built or tested (gpu-omp 1/2, then gpu-acc 2/2). Both failures show the same error on the same runner, frontier-5:

error: Failed to install: pandas-3.0.3-cp313-cp313-manylinux_2_24_x86_64.manylinux_2_28_x86_64.whl
  Caused by: failed to open file `/tmp/uv-cache-sbryngelson/archive-v0/h07_8qgFw9kMCAsSlH4jZ/pandas-3.0.3.dist-info/METADATA`: No such file or directory

Maybe a corrupted uv cache entry on that runner's node? Could you help me re-run the one failed job? All the tests should pass.

sbryngelson added a commit that referenced this pull request Jul 4, 2026
Self-hosted Frontier/Frontier-AMD matrix legs (acc/omp/cpu x shards) run
their "Fetch Dependencies" step directly on the same login node as the
same OS user, all pointed at the same UV_CACHE_DIR (introduced in #1385
to dodge NFS file-lock errors on ~/.cache/uv). uv's own cache lock
guards individual entries, but concurrent installs from separate uv
processes can still race while one extracts/prunes the shared
archive-v0 store, leaving a corrupted entry behind (e.g. a missing
dist-info METADATA file) that fails every subsequent install until the
cache is cleared by hand -- as happened on PR #1414's Frontier gpu-acc
[2/2] job.

Serialize the actual `uv pip install` call with flock so only one
process touches a given cache dir at a time, while keeping the cache
itself shared and warm across runs.
@sbryngelson

Copy link
Copy Markdown
Member

One note on CI. The last two runs each failed a single Frontier job at the Fetch Dependencies step, before anything gets built or tested (gpu-omp 1/2, then gpu-acc 2/2). Both failures show the same error on the same runner, frontier-5:

error: Failed to install: pandas-3.0.3-cp313-cp313-manylinux_2_24_x86_64.manylinux_2_28_x86_64.whl
  Caused by: failed to open file ``/tmp/uv-cache-sbryngelson/archive-v0/h07_8qgFw9kMCAsSlH4jZ/pandas-3.0.3.dist-info/METADATA``: No such file or directory

Maybe a corrupted uv cache entry on that runner's node? Could you help me re-run the one failed job? All the tests should pass.

#1630

sbryngelson added a commit that referenced this pull request Jul 4, 2026
Addresses Copilot review on #1630, plus a live recurrence caught on
this PR's own CI run (job 85133634699, "Frontier (AMD) cpu [1/2]"):

- UV_LOCK_DIR's guard only checked -w, so a writable non-directory
  TMPDIR would pass and then get used as a directory prefix, breaking
  flock. Add a -d check alongside -w, per Copilot's suggestion.

- That same CI run hit a *new* corruption symptom ("The wheel is
  invalid: Missing .dist-info directory" for pandas) on the same
  physical login node (login05) as the original incident on #1414,
  even with the new lock in place. Root cause: a cache entry corrupted
  before the lock existed (or by any other cause) just fails forever
  until someone manually clears it -- which is exactly what had
  happened here; login05's ~1.2GiB cache had never actually been
  cleared (an earlier `uv cache clean` in this investigation was run
  from a different login node's session and never touched login05).

  Since self-hosted runners are spread across login nodes we can't all
  individually SSH into and inspect every time this happens, make the
  script self-heal instead: on install failure, clear the uv cache and
  retry once before giving up.
# Conflicts:
#	toolchain/mfc/test/cases.py
# Conflicts:
#	src/simulation/m_riemann_solver_hll.fpp
#	src/simulation/m_riemann_solver_hllc.fpp
#	src/simulation/m_riemann_solver_lf.fpp
#	src/simulation/m_riemann_state.fpp
@sbryngelson

Copy link
Copy Markdown
Member

Hi @ChrisZYJ#1572 (Phase-3 Riemann hot-path decomposition) restructured/removed the exact solver code this PR builds on, which is what put m_riemann_solver_hll/hllc/lf.fpp and m_riemann_state.fpp into conflict.

I resolved the merge of current master into this branch for you. It's pushed to:

sbryngelson:pr1414-hypo-hlld-master-merge

Its first parent is your current head (815f1e2f4), so you can take it as a clean fast-forward if you want — no obligation:

git fetch https://github.com/sbryngelson/MFC.git pr1414-hypo-hlld-master-merge
git checkout hypo_hlld
git merge --ff-only FETCH_HEAD
git push

(Or just git pull https://github.com/sbryngelson/MFC.git pr1414-hypo-hlld-master-merge.) It's your PR, so please review the resolution before adopting it — the choices below involve your physics, and I have not compiled or run it (fypp preprocesses cleanly and the generated Fortran is block-balanced, but goldens are untested).

The core tension

#1572 deleted the HLLC/LF hypoelastic-energy blocks as dead code, justified by case_validator.py restricting hypoelasticity to HLL (riemann_solver == 1). This PR lifts that restriction to [1, 2, 4] (HLL/HLLC/HLLD), so those blocks are live again. The resolution keeps your feature and layers master's refactors underneath.

Choices made (10 conflict hunks)

File / region Resolution Rationale
LF hypoelastic block took master (deleted) LF isn't in [1,2,4], so the block is dead on both sides; merged LF is identical to master
HLL inline energy macro took master's helper s_compute_hypoelastic_interface_energy (called just after the loop) already reproduces your inline physics including the shear_indices doubling — keeping the inline macro too would double-count
HLLC declarations took yours (superset) your ADC/star-state locals are needed downstream
HLLC 6-eq & 5-eq private lists took yours, and added master's alpha_lim_L/R to the 5-eq list GPU-safe superset; the merged 5-eq body uses both your ADC vars and master's alpha_lim arrays
HLLC 6-eq & 5-eq hypoelastic energy took yours (live feature) validator now permits HLLC hypoelasticity
HLLC 5-eq mixture loads took master (alpha_lim_L/R = ...) the auto-merged helper call accumulates mixture props from alpha_lim; your inline accumulation would double-count against it
HLLC flux else (not-hypoelasticity) took yours your else holds the full base momentum+energy+elastic flux; master's conflict-side tail alone would leave momentum/energy flux uncomputed
m_riemann_state union kept master's three shared helpers and your templated s_finalize_riemann_solver (+_hatR); dropped master's now-redundant plain s_finalize_riemann_solver decl since your template supersedes it

Two things worth your eyes

  1. hypo_energy_guard gating. Your HLLC path uses @:compute_hypo_elastic_energy (which honors hypo_energy_guard), while HLL now goes through master's helper (which hard-gates on G > verysmall). Both are live, and I preserved each side's own choice — but HLL and HLLC now treat the soft-G case differently. You may want them consistent.
  2. Cosmetic: an unused Re_max local survives in HLLC (Phase-3: Riemann hot-path decomposition into shared GPU device helpers #1572 had removed it). Harmless; drop at your leisure.

Happy to adjust any of these if you'd prefer a different call. cc @sbryngelson

ChrisZYJ added 7 commits July 11, 2026 18:50
…er (drop dead six-eqn block, unused Re_max); fix HLL right-state damage index
…n hypo feature constraints with validator parity, constraint tests, and damage/three-fluid regressions
…near-zero-G elastic energy guards; regenerate the 12 affected hypo goldens
Mirror the Python validator's model and solver restrictions with Fortran PROHIBIT checks, use named solver constants in the touched checker logic, and remove the one-off parameter pytest.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

enhancement New feature or request

Development

Successfully merging this pull request may close these issues.

3 participants