Add convergence-rate CI for WENO5/3/1 and MUSCL2

[sbryngelson]

Closes #305

Summary

Adds a convergence CI workflow (.github/workflows/convergence.yml) with three jobs:

1D spatial convergence (run_convergence_1d.py) — single-fluid Euler density sine wave, one period advection; all schemes with 4 MPI ranks:

• WENO5: N ∈ [128, 512], rate ≥ 4.8
• WENO3: N ∈ [256, 1024], rate ≥ 1.8 (JS smooth-extremum asymptote = 2)
• WENO1: N ∈ [128, 1024], rate ≥ 0.95
• MUSCL2: N ∈ [128, 1024], rate ≥ 1.9 (unlimited slope)
• TENO5: N ∈ [128, 512], rate ≥ 4.8
• WENO7: N ∈ [64, 128], CFL=0.005, rate ≥ 6.5 (capped before round-off floor at N=256)
• TENO7: N ∈ [64, 128], CFL=0.005, rate ≥ 6.5

2D spatial convergence (run_convergence.py) — isentropic vortex eps=0.01 with Gauss-Legendre IC (hcid=283); N ∈ [32, 128], 4 MPI ranks:

• WENO5, WENO3, WENO1, MUSCL2, TENO5 (WENO7/TENO7 excluded — covariance floor dominates at these resolutions)

Temporal order (run_temporal_order.py) — fixed N=512/WENO5, vary CFL with 4 MPI ranks:

• RK1 (Forward Euler): CFL=[0.10, 0.05], rate ≥ 0.9
• RK2 (TVD Heun): CFL=[0.50, 0.25], rate ≥ 1.8
• RK3 (TVD Shu-Osher): CFL=[0.50, 0.25], rate ≥ 2.7

Also adds examples/1D_euler_convergence/ and examples/2D_isentropicvortex_convergence/ benchmark cases, and --time-stepper argument to the 1D case.

Test plan

• All three CI jobs pass locally with 4 MPI ranks
• 1D: WENO5≈5.0, WENO3≈1.9, WENO1≈1.0, MUSCL2≈2.0, TENO5≈5.0, WENO7≈7.0, TENO7≈7.0
• 2D: WENO5≈4.8, WENO3≈2.0, WENO1≈0.9, MUSCL2≈2.2, TENO5≈4.6
• Temporal: RK1≈1.0, RK2≈2.0, RK3≈2.95

[github-actions]

Claude Code Review

Head SHA: 3a069d2

Files changed:

• 13
• .github/workflows/convergence.yml
• examples/1D_advection_convergence/case.py
• examples/1D_euler_convergence/case.py
• examples/2D_isentropicvortex_convergence/case.py
• src/common/include/2dHardcodedIC.fpp
• src/simulation/m_muscl.fpp
• toolchain/mfc/case_validator.py
• toolchain/mfc/params/definitions.py
• toolchain/mfc/params/registry.py
• toolchain/mfc/test/cases.py

Findings

Fortran-side checker for muscl_lim bound likely not updated

Files: toolchain/mfc/case_validator.py:790, src/simulation/m_muscl.fpp:172

The PR extends muscl_lim's valid range from {1,2,3,4,5} to {0,1,2,3,4,5} in three places: case_validator.py (Python runtime bound changed from < 1 to < 0), definitions.py (choices list updated), and m_muscl.fpp (new if (muscl_lim == 0) branch added). However, no Fortran checker file (m_checker.fpp, m_checker_common.fpp, or similar) appears in the diff.

MFC duplicates physics-constraint validation in both Python (case_validator.py) and Fortran (m_checker*.fpp) via @:PROHIBIT() macros. The changed line in case_validator.py:

# before
self.prohibit(muscl_lim is not None and (muscl_lim < 1 or muscl_lim > 5), ...)
# after
self.prohibit(muscl_lim is not None and (muscl_lim < 0 or muscl_lim > 5), ...)

strongly implies a parallel Fortran check of the form @:PROHIBIT(muscl_lim < 1 .or. muscl_lim > 5, ...) exists in a checker file. If that check was not updated, any simulation with muscl_lim = 0 will abort at Fortran startup even though the Python validator accepts it, making the new unlimited-limiter mode completely non-functional at runtime despite appearing correct at the source level.

---

Hardcoded IC 283: kinetic energy back-conversion inconsistency

File: src/common/include/2dHardcodedIC.fpp (case 283)

The GL loop accumulates total energy from point-wise values:

Eq = pq/0.4_wp + 0.5_wp*rhoq*(uq**2 + vq**2)

The back-conversion to pressure then computes:

q_prim_vf(eqn_idx%E)%sf(i, j, 0) = (E_avg - 0.5_wp*(rhou_avg**2 + rhov_avg**2)/rho_avg)*0.4_wp

The kinetic energy subtracted during back-conversion, 0.5*(rhou_avg²+rhov_avg²)/rho_avg, is evaluated from the conserved averages (i.e., ⟨ρu⟩²/⟨ρ⟩), whereas the kinetic energy accumulated into E_avg during quadrature is ⟨ρu²⟩ (point-wise). By Jensen's inequality, ⟨ρu²⟩ ≥ ⟨ρu⟩²/⟨ρ⟩, so the recovered pressure p_avg is higher than the true GL average of pq. The error is O(h²) covariance, which the design comments (lines 1072–1077 of run_convergence.py) argue stays below scheme error for the weak-vortex (ε=0.01) regime at N≤46. At larger N this O(h²) floor sets the convergence ceiling for all high-order schemes. Confirm the documented crossover bound h > ε^(1/3) ≈ 0.22 (i.e., N < 46) actually keeps the covariance error subdominant at the tested resolutions (N=32–128 in CI) — if not, WENO5 will stall below the 4.0 threshold before reaching the regime where it shows true 5th-order rate.

---

CI workflow triggers on every push without branch filtering

File: .github/workflows/convergence.yml:9-11

on:
  push:
  workflow_dispatch:

There is no branches: filter. The 1D run at --resolutions 128 256 512 1024 with the default --cfl 0.02 requires ~112,000 time steps at N=1024; running all four schemes across four resolutions will likely approach or exceed the 60-minute job timeout. The plan document (Task 3 in docs/superpowers/plans/2026-05-05-convergence-ci.md) specified:

on:
  push:
    branches: [master]
  pull_request:
  workflow_dispatch:

The implemented workflow omits branches: [master] and pull_request:, diverging from the intended design and causing the expensive 1D job to fire on every push to every branch.