Add TENO5, WENO7, TENO7 to convergence runners

Both 1D and 2D runners now test all 7 schemes: WENO5/3/1, MUSCL2, TENO5,
WENO7, TENO7.  Both case.py files gain --teno and --teno-ct flags.

Resolution bounds:
  1D: TENO5 [128,512], WENO7/TENO7 [128,256] (precision floor near N=512)
  2D: TENO5 [32,128],  WENO7/TENO7 [64,128]  (MFC stencil constraint: N>=35)

TENO CT values match the Shu-Osher examples: 1e-6 for order-5, 1e-9 for order-7.
weno_eps tightened to 1e-40 (was 1e-16) for all WENO/TENO schemes in case.py.

Author: sbryngelson

examples/1D_euler_convergence/case.py

@@ -18,8 +18,10 @@
 parser = argparse.ArgumentParser(description="1D Euler convergence case")
 parser.add_argument("--mfc", type=json.loads, default="{}", metavar="DICT")
 parser.add_argument("-N", type=int, default=64, help="Grid points (default: 64)")
-parser.add_argument("--order", type=int, default=5, help="WENO order: 1, 3, or 5")
+parser.add_argument("--order", type=int, default=5, help="WENO order: 1, 3, 5, or 7")
 parser.add_argument("--muscl", action="store_true", help="Use MUSCL-2 instead of WENO")
+parser.add_argument("--teno", action="store_true", help="Use TENO instead of WENO")
+parser.add_argument("--teno-ct", type=float, default=1e-6, help="TENO CT threshold (default: 1e-6)")
 parser.add_argument("--cfl", type=float, default=0.4, help="CFL number (default: 0.4)")
 parser.add_argument("--no-mapped", action="store_true", help="Disable mapped WENO")
 parser.add_argument("--muscl-lim", type=int, default=0, help="MUSCL limiter: 0=unlimited 1=minmod ... (default: 0)")
@@ -48,12 +50,14 @@
     scheme_params = {
         "recon_type": 1,
         "weno_order": args.order,
-        "weno_eps": 1.0e-16,
+        "weno_eps": 1.0e-40,
         "weno_Re_flux": "F",
         "weno_avg": "F",
-        "mapped_weno": "F" if (args.order == 1 or args.no_mapped) else "T",
+        "mapped_weno": "F" if (args.order == 1 or args.no_mapped or args.teno) else "T",
         "null_weights": "F",
         "mp_weno": "F",
+        "teno": "T" if args.teno else "F",
+        **({"teno_CT": args.teno_ct} if args.teno else {}),
     }
 
 print(

examples/2D_isentropicvortex_convergence/case.py

@@ -13,8 +13,10 @@
 parser = argparse.ArgumentParser(description="2D isentropic vortex convergence case")
 parser.add_argument("--mfc", type=json.loads, default="{}", metavar="DICT")
 parser.add_argument("-N", type=int, default=32, help="Grid points per dim (default: 32)")
-parser.add_argument("--order", type=int, default=5, help="WENO order: 1, 3, or 5 (default: 5)")
+parser.add_argument("--order", type=int, default=5, help="WENO order: 1, 3, 5, or 7 (default: 5)")
 parser.add_argument("--muscl", action="store_true", help="Use MUSCL instead of WENO")
+parser.add_argument("--teno", action="store_true", help="Use TENO instead of WENO")
+parser.add_argument("--teno-ct", type=float, default=1e-6, help="TENO CT threshold (default: 1e-6)")
 parser.add_argument("--muscl-lim", type=int, default=0, help="MUSCL limiter: 0=unlimited 1=minmod ... (default: 0)")
 args = parser.parse_args()
 
@@ -41,10 +43,12 @@
     scheme_params = {
         "recon_type": 1,
         "weno_order": args.order,
-        "weno_eps": 1.0e-16,
-        "mapped_weno": "F" if args.order == 1 else "T",
+        "weno_eps": 1.0e-40,
+        "mapped_weno": "F" if (args.order == 1 or args.teno) else "T",
         "null_weights": "F",
         "mp_weno": "F",
+        "teno": "T" if args.teno else "F",
+        **({"teno_CT": args.teno_ct} if args.teno else {}),
     }
 
 print(

toolchain/mfc/test/run_convergence.py

@@ -11,8 +11,12 @@
 L2(rho(T) - rho(0)) measures accumulated scheme error; the comparison to rho(0)
 (the numerical IC) eliminates IC discretisation error, isolating the scheme error.
 
+WENO7/TENO7 require min N=64 per dimension (MFC constraint: N >= 5 * weno_order = 35).
+With only N=64 and N=128 available from the default resolution set, the fitted rate
+is a single pairwise value; thresholds are set conservatively.
+
 Usage:
-    python toolchain/mfc/test/run_convergence.py [--no-build] [--resolutions 16 32 64]
+    python toolchain/mfc/test/run_convergence.py [--no-build] [--resolutions 32 64 128]
 """
 
 import argparse
@@ -30,19 +34,24 @@
 CASE = "examples/2D_isentropicvortex_convergence/case.py"
 MFC = "./mfc.sh"
 
-# (label, extra_args, expected_order, tolerance)
+# (label, extra_args, expected_order, tolerance, min_N, max_N)
 # With eps=0.01 and N=32..128 the prim->cons covariance error O(eps^3 h^2) is
 # well below the scheme's spatial error O(eps^2 h^p), so each scheme shows its
 # nominal rate.  The tolerance is the allowable shortfall from the nominal order.
 #
-# WENO3 note: at N=32-128 the rate is ~2.0-2.2 (pre-asymptotic; approaches 3
-# at finer grids) — still clearly above the 1D result (1.77, smooth-extremum
-# degradation), which is what this test is designed to show.  Threshold 1.8.
+# WENO3: at N=32-128 the rate is ~2.0-2.2 (pre-asymptotic; approaches 3 at
+#   finer grids).  Threshold 1.8.
+# WENO7/TENO7: require N >= 35 (MFC stencil constraint), so min_N=64.  With
+#   only N=64,128 in the default set the rate is a single pairwise value;
+#   thresholds set conservatively pending actual run data.
 SCHEMES = [
-    ("WENO5", ["--order", "5"], 5, 1.0),
-    ("WENO3", ["--order", "3"], 3, 1.2),
-    ("WENO1", ["--order", "1"], 1, 0.4),
-    ("MUSCL2", ["--muscl"], 2, 0.5),
+    ("WENO5", ["--order", "5"], 5, 1.0, 32, None),
+    ("WENO3", ["--order", "3"], 3, 1.2, 32, None),
+    ("WENO1", ["--order", "1"], 1, 0.4, 32, None),
+    ("MUSCL2", ["--muscl"], 2, 0.5, 32, None),
+    ("TENO5", ["--order", "5", "--teno", "--teno-ct", "1e-6"], 5, 1.0, 32, None),
+    ("WENO7", ["--order", "7"], 7, 3.0, 64, None),
+    ("TENO7", ["--order", "7", "--teno", "--teno-ct", "1e-9"], 7, 3.0, 64, None),
 ]
 
 
@@ -119,7 +128,11 @@ def run_case(tmpdir: str, N: int, extra_args: list):
     return Nt, os.path.join(tmpdir, f"N{N}")
 
 
-def test_scheme(label, extra_args, expected_order, tol, resolutions):
+def test_scheme(label, extra_args, expected_order, tol, resolutions, min_N=None, max_N=None):
+    if min_N is not None:
+        resolutions = [N for N in resolutions if N >= min_N]
+    if max_N is not None:
+        resolutions = [N for N in resolutions if N <= max_N]
     print(f"\n{'=' * 60}")
     print(f"  {label}  (need rate >= {expected_order - tol:.1f})")
     print(f"{'=' * 60}")
@@ -166,7 +179,7 @@ def main():
     parser = argparse.ArgumentParser(description="MFC convergence-rate verification")
     parser.add_argument("--no-build", action="store_true", help="Skip build step")
     parser.add_argument("--resolutions", type=int, nargs="+", default=[32, 64, 128], help="Grid resolutions (default: 32 64 128; N<32 unsupported for WENO5)")
-    parser.add_argument("--schemes", nargs="+", default=["WENO5", "WENO3", "WENO1", "MUSCL2"], help="Schemes to test (default: all)")
+    parser.add_argument("--schemes", nargs="+", default=["WENO5", "WENO3", "WENO1", "MUSCL2", "TENO5", "WENO7", "TENO7"], help="Schemes to test (default: all)")
     args = parser.parse_args()
 
     if not args.no_build:
@@ -180,11 +193,11 @@ def main():
             sys.exit(1)
 
     results = {}
-    for label, extra_args, expected_order, tol in SCHEMES:
+    for label, extra_args, expected_order, tol, min_N, max_N in SCHEMES:
         if label not in args.schemes:
             continue
         try:
-            passed = test_scheme(label, extra_args, expected_order, tol, args.resolutions)
+            passed = test_scheme(label, extra_args, expected_order, tol, args.resolutions, min_N, max_N)
         except Exception as e:
             print(f"  ERROR: {e}")
             passed = False

toolchain/mfc/test/run_convergence_1d.py

@@ -7,20 +7,25 @@
 L2(rho(T) - rho(0)) measures the accumulated scheme spatial truncation error.
 No non-conservative alpha equation — clean benchmark for all schemes.
 
-CFL=0.02 by default so that RK3 temporal error O(dt^3)~O(h^3) is negligible
-relative to the spatial error at all tested resolutions, allowing WENO5 to
-show its true 5th-order rate.
+CFL=0.02 by default so that RK3 temporal error O(dt^3) is negligible relative
+to the spatial error at all tested resolutions, allowing WENO5/7 to show their
+true spatial rates.
 
-WENO3-JS degrades to 2nd order at smooth extrema (Henrick et al. 2005, the
-Jiang-Shu smoothness indicators do not converge to optimal weights at critical
-points where f'=0).  The expected rate for WENO3 here is therefore 2, not 3;
-the 2D isentropic vortex test (run_convergence.py) verifies WENO3 rate 3 on
-a problem without smooth-extremum contamination.
+WENO3-JS degrades to 2nd order at smooth extrema (Henrick et al. 2005).
+The expected rate for WENO3 here is therefore 2, not 3; the 2D isentropic
+vortex test (run_convergence.py) verifies WENO3 rate 3.
 
 MUSCL2 uses muscl_lim=0 (unlimited central-difference) by default.  TVD
 limiters clip slopes to zero at smooth extrema and stall at 1st order on the
 sine wave; the unlimited limiter preserves 2nd-order convergence everywhere.
 
+TENO5 uses the same 5th-order stencil as WENO5 with threshold-based stencil
+selection (CT=1e-6).  On smooth problems all stencils are selected and the
+rate equals WENO5; TENO's advantage is sharper shock capturing.
+
+WENO7/TENO7: capped at N=256 — the machine-precision floor (~2e-15) is
+reached near N=512 for 7th-order schemes on this smooth problem.
+
 Usage:
     python toolchain/mfc/test/run_convergence_1d.py [--no-build] [--resolutions 32 64 128]
 """
@@ -43,18 +48,25 @@
 # (label, extra_args, expected_order, tolerance, min_N, max_N)
 # Per-scheme resolution bounds let each scheme run over the range where its
 # asymptotic order is cleanly visible:
-#   WENO5 : cap at N=512 — double-precision floor kills the rate at N=1024
-#           (error ~2.6e-12, rate collapses to 0.69); [128,512] gives 4.99.
-#   WENO3 : start at N=256 — skips the coarsest pre-asymptotic points;
-#           WENO3-JS degrades to 2nd order at smooth extrema (Henrick 2005),
-#           asymptote confirmed 1.99 at N=4096; [256,1024] gives ~1.87.
-#   WENO1 : full range [128,1024]; rate 0.97.
-#   MUSCL2: full range [128,1024]; unlimited slope, rate exactly 2.00.
+#   WENO5  : cap at N=512 — double-precision floor kills the rate at N=1024
+#            (error ~2.6e-12, rate collapses to 0.69); [128,512] gives 4.99.
+#   WENO3  : start at N=256 — skips the coarsest pre-asymptotic points;
+#            WENO3-JS degrades to 2nd order at smooth extrema (Henrick 2005),
+#            asymptote confirmed 1.99 at N=4096; [256,1024] gives ~1.87.
+#   WENO1  : full range [128,1024]; rate 0.97.
+#   MUSCL2 : full range [128,1024]; unlimited slope, rate exactly 2.00.
+#   TENO5  : same range as WENO5; CT=1e-6; rate matches WENO5 on smooth problems.
+#   WENO7  : cap at N=256 — machine-precision floor (~2e-15) reached near N=512
+#            for 7th-order schemes on this smooth problem.
+#   TENO7  : same range as WENO7; CT=1e-9.
 SCHEMES = [
     ("WENO5", ["--order", "5"], 5, 0.2, 128, 512),
     ("WENO3", ["--order", "3"], 2, 0.2, 256, None),
     ("WENO1", ["--order", "1"], 1, 0.05, 128, None),
     ("MUSCL2", ["--muscl"], 2, 0.1, 128, None),
+    ("TENO5", ["--order", "5", "--teno", "--teno-ct", "1e-6"], 5, 0.2, 128, 512),
+    ("WENO7", ["--order", "7"], 7, 1.0, 128, 256),
+    ("TENO7", ["--order", "7", "--teno", "--teno-ct", "1e-9"], 7, 1.0, 128, 256),
 ]
 
 
@@ -182,7 +194,7 @@ def main():
     parser.add_argument(
         "--schemes",
         nargs="+",
-        default=["WENO5", "WENO3", "WENO1", "MUSCL2"],
+        default=["WENO5", "WENO3", "WENO1", "MUSCL2", "TENO5", "WENO7", "TENO7"],
         help="Schemes to test",
     )
     parser.add_argument("--cfl", type=float, default=0.02, help="CFL number (default: 0.02; small so RK3 temporal error is negligible)")