MHD branches in m_riemann_solver_lf.fpp are unreachable (validator forbids MHD + Lax-Friedrichs)

[sbryngelson]

Summary

The MHD branches in src/simulation/m_riemann_solver_lf.fpp are unreachable through any supported entry point, because the case validator rejects MHD with the Lax-Friedrichs solver:

https://github.com/MFlowCode/MFC/blob/master/toolchain/mfc/case_validator.py#L1048

self.prohibit(mhd and riemann_solver is not None and riemann_solver not in [1, 4],
              "MHD simulations require riemann_solver = 1 (HLL) or riemann_solver = 4 (HLLD)")

Attempting to run any MHD case with riemann_solver = "lax_friedrichs" fails before the solver is ever entered:

Error: Case validation failed for simulation:
Case parameter constraint violations:

  1. MHD simulations require riemann_solver = 1 (HLL) or riemann_solver = 4 (HLLD)

Yet m_riemann_solver_lf.fpp carries full ideal-MHD and relativistic-MHD state and flux code — the else if (mhd .and. relativity) and else if (mhd .and. .not. relativity) branches around lines 253-287, and the ideal-MHD energy flux around line 395.

The restriction dates back to #1066 ("Case checking moved to Python"), so this is long-standing rather than a recent regression.

Why it matters

Unreachable physics code silently accumulates bugs, because nothing can exercise it. Concretely: while reviewing #1658 I found that m_riemann_solver_lf.fpp:395 guarded the ideal-MHD energy flux with num_dims > 1 while every sibling guard used num_dims > 2. That inconsistency would have broken 1D ideal-MHD LF under --case-optimization — but it went unnoticed indefinitely, and could not have been caught by any test, because the path cannot be reached. #1658 corrects the guard, but the underlying exposure remains.

Secondary concern: the check is Python-side only

I could not find an equivalent constraint in src/simulation/m_checker.fpp or src/common/m_checker_common.fpp. A hand-written simulation.inp that bypasses the Python toolchain would therefore reach these branches with no diagnostic. If the restriction is intended to be a hard physics constraint, it likely belongs in the Fortran checker as well (or in addition).

Suggested resolution

Pick one:

1. Remove the MHD branches from m_riemann_solver_lf.fpp. Cheapest option; the code is not reachable and not tested. Reduces the surface area of guard-consistency bugs like the one above.
2. Lift the validator restriction and support MHD + LF properly — add at least one MHD LF case to the golden test suite, and one to the case-optimization CI, so the path is actually covered.

Either way, adding the Fortran-side check (or documenting why it is deliberately Python-only) would close the .inp bypass.

Related

• Fixes #1656: preserve MHD calculations with case optimization #1658 — fixes the case-optimization guards in both HLL and LF; the LF half of that diff touches this unreachable code.
• Case optimization compiles out required MHD calculations in 1D/2D HLL cases #1656 — the ideal-MHD case-optimization NaN bug that prompted the review.

---

Found in part via Claude Code, during a review of #1658. The validator block surfaced when an MHD + Lax-Friedrichs case was actually built and run rather than only read; the reachability conclusion and the Fortran-checker gap were then confirmed against the source.