Block-structured adaptive mesh refinement (AMR)

.github/scripts/check_caseopt_shared_slugs.py

@@ -0,0 +1,78 @@
+#!/usr/bin/env python3
+"""Assert every benchmark case shares benchmarks[0]'s slug for syscheck, pre_process and post_process.
+
+prebuild-case-optimization.sh builds those once from benchmarks[0] behind a done-marker; a case with a different
+slug (changed chemistry, eos_state_dependent or analytic ICs) would build its staging dir concurrently, unserialized.
+
+Usage: check_caseopt_shared_slugs.py <acc|mp|no> <case.py> [<case.py> ...]
+"""
+
+import os
+import sys
+
+SHARED_TARGETS = ("syscheck", "pre_process", "post_process")
+
+
+def main(argv):
+    if len(argv) < 3:
+        print(f"usage: {os.path.basename(argv[0])} <acc|mp|no> <case.py> [<case.py> ...]", file=sys.stderr)
+        return 2
+
+    interface, cases = argv[1], argv[2:]
+
+    sys.path.insert(0, os.path.join(os.getcwd(), "toolchain"))
+    from mfc import state
+    from mfc.state import MFCConfig
+
+    state.gCFG = MFCConfig(gpu=interface)
+    # Only keys load() forwards to case.py as DICT (benchmarks read gpu/nodes/tasks_per_node).
+    state.gARG = {
+        "gpu": interface != "no",
+        "nodes": 1,
+        "tasks_per_node": 1,
+        "case_optimization": True,
+        "input": None,
+        "--": [],
+        "rdma_mpi": False,
+    }
+
+    from mfc.build import get_target
+    from mfc.run import input as mfc_input
+
+    slugs = {}
+    for case in cases:
+        state.gARG["input"] = case
+        ifile = mfc_input.load(case, [], {}, do_print=False)
+        slugs[case] = {name: get_target(name).get_slug(ifile) for name in SHARED_TARGETS}
+
+    reference = cases[0]
+    mismatches = []
+    for case in cases[1:]:
+        for name in SHARED_TARGETS:
+            if slugs[case][name] != slugs[reference][name]:
+                mismatches.append((name, case, slugs[case][name], slugs[reference][name]))
+
+    if mismatches:
+        print("ERROR: case-optimization benchmarks no longer share one slug per shared target.", file=sys.stderr)
+        print(f"       Reference case (benchmarks[0]): {reference}", file=sys.stderr)
+        for name, case, got, want in mismatches:
+            print(f"         {name}: {case}", file=sys.stderr)
+            print(f"           got  {got}", file=sys.stderr)
+            print(f"           want {want}", file=sys.stderr)
+        print(
+            "       prebuild-case-optimization.sh builds these three targets once, from the reference\n"
+            "       case, and lets other shards skip them -- which serializes only the reference slug.\n"
+            "       Give every shard its own shared-target build, or drop the offending case, before\n"
+            "       relying on the marker.",
+            file=sys.stderr,
+        )
+        return 1
+
+    print(f"Shared-target slugs agree across {len(cases)} benchmark case(s):")
+    for name in SHARED_TARGETS:
+        print(f"  {name}: {slugs[reference][name]}")
+    return 0
+
+
+if __name__ == "__main__":
+    sys.exit(main(sys.argv))

.github/scripts/prebuild-case-optimization.sh

@@ -59,6 +59,14 @@ esac
 # Concurrent shards must not build those shared staging dirs simultaneously:
 # shard 1 builds them first and drops a done marker; other shards wait for it,
 # after which their builds no-op in the shared dirs.
+# The marker covers only the benchmarks[0] slug, so any other shared slug would build unserialized.
+# Slugs move with compile-time axes in Case.get_fpp (chemistry, eos_state_dependent, analytic ICs): check up front in every shard.
+(
+    . ./mfc.sh venv > /dev/null
+    python3 .github/scripts/check_caseopt_shared_slugs.py \
+        "$([ "$job_interface" = "omp" ] && echo mp || echo acc)" "${benchmarks[@]}"
+)
+
 if [ -n "$shard" ] && [ "$shard_count" -gt 1 ]; then
     shared_marker_done="build/.prebuild-shared-targets-done"
     shared_marker_failed="build/.prebuild-shared-targets-failed"

.github/scripts/submit-slurm-job.sh

@@ -153,11 +153,12 @@ elif [ "$device" = "gpu" ]; then
 
     case "$cluster" in
         phoenix)
+            # 8 tasks: the largest AMR test (3D, ppn=8) needs 8 MPI slots.
             # --exclude is rendered separately (see $node_exclude) so the
             # preflight can add a node to it and resubmit.
             sbatch_device_opts="\
 #SBATCH -p $gpu_partition
-#SBATCH --ntasks-per-node=4
+#SBATCH --ntasks-per-node=8
 #SBATCH -G${gpu_count}"
             node_exclude="atl1-1-03-007-29-0,atl1-1-03-007-31-0,atl1-1-01-002-28-0"
             ;;

.github/workflows/common/build.sh

@@ -54,16 +54,17 @@ if [ "$job_cluster" = "phoenix" ]; then
 fi
 
 # --- Variant selection ---
-# The suite needs two simulation binaries: the base one and a chemistry one (the
-# mechanism is compiled in). On amdflang each device link is ~1 hour, so building
-# both serially exceeds the 2 h walltime. $job_variant lets the caller split them
+# The suite needs three simulation binaries (base, chem, eos_state_dependent), each with its own build slug.
+# Serial amdflang device links exceed the 2 h walltime, so $job_variant lets the caller split them
 # into concurrent jobs that write to disjoint staging directories:
-#   base -> plain targets, what every non-chemistry test runs
+#   base -> empty case (chemistry=False, eos_state_dependent=False)
 #   chem -> the chemistry variant only
+#   eos  -> eos_state_dependent=True (Mie-Gruneisen, JWL, Vinet); one slug per target, so any such case builds it
 # Unset builds everything in this job (default for every other cluster).
 case "${job_variant:-}" in
     base) build_cmd=(./mfc.sh build -j 8 $build_opts) ;;
     chem) build_cmd=(./mfc.sh test -v --dry-run -a -j 8 -o Chemistry $build_opts) ;;
+    eos)  build_cmd=(./mfc.sh test -v --dry-run -a -j 8 -o eos=mie_gruneisen $build_opts) ;;
     "")   build_cmd=(./mfc.sh test -v --dry-run -a -j 8 $build_opts) ;;
     *)    echo "ERROR: unknown job_variant '$job_variant'"; exit 1 ;;
 esac

.github/workflows/test.yml

@@ -472,13 +472,11 @@ jobs:
 
       - name: Build (concurrent variants)
         if:   matrix.cluster == 'frontier_amd' && matrix.device == 'gpu'
-        # The suite needs a base and a chemistry simulation binary, and on amdflang
-        # each device LTO link is ~1 h -- building both serially overruns the 2 h
-        # walltime. They touch disjoint staging directories, so submit one SLURM job
-        # per variant and wait on both (mirrors the case-optimization pre-build).
+        # base, chem and eos_state_dependent simulation binaries (see common/build.sh) use disjoint staging dirs; submit one
+        # SLURM job per variant so amdflang LTO links do not serialize past the 2 h walltime, and wait on all.
         run:  |
           pids=""
-          for v in base chem; do
+          for v in base chem eos; do
             bash .github/scripts/submit-slurm-job.sh .github/workflows/common/build.sh ${{ matrix.device }} ${{ matrix.interface }} ${{ matrix.cluster }} "${{ matrix.shard }}" "$v" &
             pids="$pids $!"
           done

.gitignore

@@ -37,6 +37,10 @@ docs/documentation/parameters.md
 /tests/*/**
 !/tests/*/golden.txt
 !/tests/*/golden-metadata.txt
+# Written by ./mfc.sh test when cases fail and deleted again on a clean run, so a committed copy
+# flips in and out of every working tree depending on the last run's outcome. The pattern above
+# only covers tests/<uuid>/, not tests/ itself.
+/tests/failed_uuids.txt
 
 # NVIDIA Nsight Compute
 *.nsys-rep

.lychee.toml

@@ -33,6 +33,7 @@ exclude = [
     "https://code\\.visualstudio\\.com/?$",           # Root page returns 403 to automated requests
     "https://stackoverflow\\.com",                    # Returns 403 to automated requests
     "https://marketplace\\.visualstudio\\.com",        # Returns 503 to automated requests
+    "_8md\\.html$",  # Doxygen auto-links backticked *.md filenames in prose to per-file pages it never generates for markdown inputs; the real md_*.html page links are still checked
     "https://web\\.eng\\.ucsd\\.edu",                    # San Diego mechanism page has an untrusted SSL cert
     "https://www\\.openmp\\.org/spec-html/5\\.1/openmp\\.html",  # Times out from CI runners
     "https://openmp\\.llvm\\.org/design/Runtimes\\.html",           # Intermittent connection failures from CI runners

.typos.toml

@@ -7,6 +7,9 @@ extend-ignore-identifiers-re = [
 AttributeIDSupressMenu = "AttributeIDSupressMenu"
 
 [default.extend-words]
+# Cray CCE spells it this way in the lib-4425 runtime error; quoted verbatim in the AMR ledger so the
+# message stays greppable against what the machine actually prints.
+Unitialized = "Unitialized"
 INOUT = "INOUT"
 WRONLY = "WRONLY"
 nd = "nd"
@@ -22,6 +25,9 @@ TKE = "TKE"
 HSA = "HSA"
 infp = "infp"
 Sur = "Sur"
+thi = "thi"           # AMR clustering local: tagged-box hi index (tlo/thi)
+alo = "alo"           # AMR clustering local: accepted-box lo array (alo/ahi)
+thr = "thr"           # AMR clustering local: min-separation merge threshold
 equil = "equil"       # abbreviation for "equilibrium" (flamelet chemistry)
 chioces = "chioces"   # typo for "choices" - tests constraint key validation
 reqires = "reqires"   # typo for "requires" - tests dependency key validation

CMakeLists.txt

@@ -179,10 +179,11 @@ if (MFC_SIMULATION)
     MFC_SETUP_TARGET(TARGET  simulation
                      SOURCES "${simulation_SRCs}"
                      MPI FFTW OpenACC OpenMP)
-    # CCE 19.0.0 optcg ICE: IPA inlining triggers two distinct crashes:
+    # CCE 19.0.0 IPA inlining bugs:
     #   m_bubbles_EL:  castIsValid assertion (InstCombine/foldIntegerTypedPHI)
     #   m_phase_change, m_mpi_common, m_start_up: bring_routine_resident SIGSEGV
-    # GPU builds: apply -Oipa0 per-file to known crash sites. Whole-target
+    #   m_amr_regrid: kernel routines inlined into flattened host loops drop stores (garbage regrid indices)
+    # GPU builds: apply -Oipa0 per-file to the files above. Whole-target
     # -Oipa0 is not safe for GPU builds — it breaks cross-file !$acc routine /
     # !DIR$ INLINEALWAYS device-call registration, causing runtime aborts.
     # CPU builds: disable IPA for the whole target (safe — no GPU requirements,
@@ -195,12 +196,18 @@ if (MFC_SIMULATION)
                 "${CMAKE_BINARY_DIR}/fypp/simulation/m_phase_change.fpp.f90"
                 "${CMAKE_BINARY_DIR}/fypp/simulation/m_mpi_common.fpp.f90"
                 "${CMAKE_BINARY_DIR}/fypp/simulation/m_start_up.fpp.f90"
+                "${CMAKE_BINARY_DIR}/fypp/simulation/m_amr_regrid.fpp.f90"
                 PROPERTIES COMPILE_OPTIONS "-Oipa0"
             )
         else()
             target_compile_options(simulation PRIVATE -Oipa0)
         endif()
     endif()
+    # nvfortran 25.5-26.3 -mp=gpu (not 25.7): fort2 SIGSEGV while emitting -Minfo for m_amr_exchange
+    if (CMAKE_Fortran_COMPILER_ID STREQUAL "NVHPC" AND MFC_OpenMP)
+        set_source_files_properties("${CMAKE_BINARY_DIR}/fypp/simulation/m_amr_exchange.fpp.f90"
+                                    PROPERTIES COMPILE_OPTIONS "-Mnoinfo")
+    endif()
 endif()
 
 if (MFC_POST_PROCESS)

cmake/Fypp.cmake

@@ -120,6 +120,7 @@ macro(HANDLE_SOURCES target useCommon)
                                  -D MFC_COMPILER="${CMAKE_Fortran_COMPILER_ID}"
 				                 -D MFC_CASE_OPTIMIZATION=False
                                  -D chemistry=False
+                                 -D eos_state_dependent=False
                                  --line-numbering
                                  --no-folding
 								 --line-length=999

cmake/GPU.cmake

@@ -90,9 +90,11 @@ elseif (CMAKE_Fortran_COMPILER_ID STREQUAL "Cray")
     add_link_options("SHELL:-hkeepfiles")
 
     if (CMAKE_BUILD_TYPE STREQUAL "Debug")
+        # -h bounds: Cray array-bounds/pointer checking, matching the other compilers' debug branches above.
         add_compile_options(
                 "SHELL:-h acc_model=auto_async_none"
                 "SHELL: -h acc_model=no_fast_addr"
+                "SHELL: -h bounds"
                 "SHELL: -K trap=fp" "SHELL: -g" "SHELL: -O0"
         )
         add_link_options("SHELL: -K trap=fp" "SHELL: -g" "SHELL: -O0")

docs/documentation/amr.md

@@ -0,0 +1,132 @@
+@page amr Adaptive Mesh Refinement
+
+# Adaptive Mesh Refinement
+
+> **Experimental.** AMR is off by default (`amr = F`). The case validator rejects unsupported combinations by name.
+
+## Overview {#amr-overview}
+
+AMR adds refined rectangular blocks on top of the base grid (level 0) where the flow needs resolution: shocks,
+interfaces, wakes, bodies. The base grid is advanced by the ordinary solver; refined blocks run the same solver on
+their own grids, nested up to `amr_max_level` levels at `amr_ref_ratio`:1 per level.
+
+With `amr_regrid_int > 0`, every `amr_regrid_int` steps the blocks are rebuilt around the cells where the density
+gradient exceeds `amr_tag_eps`. With `amr_regrid_int = 0` a single block stays at `amr_block_beg`..`amr_block_end`.
+
+AMR lives in `simulation` (`src/simulation/m_amr*.fpp`); post_process can overlay the blocks. Internals:
+@subpage amr_implementation.
+
+## Algorithm {#amr-algorithm}
+
+- **Tagging and clustering.** A coarse cell is tagged when its normalized density jump exceeds `amr_tag_eps`
+  (`amr_tag_abs = T` uses the absolute jump instead). Tags are clustered by Berger-Rigoutsos bisection into boxes
+  that reach `amr_cluster_eff` tag efficiency, are at least `amr_blocking_factor` cells wide, are padded by `amr_buf`
+  cells and are capped at `amr_max_grid_size` cells per dimension. Levels nest properly (a level-2 box lies inside one
+  level-1 box).
+- **Prolongation.** New blocks and fine ghost cells are filled from the parent level by conservative minmod-limited
+  linear interpolation (`amr_prolong = 1`) or conservative quartic interpolation (`amr_prolong = 2`). The quartic fill
+  falls back to limited linear wherever it would give non-positive density or internal energy. Volume fractions are
+  renormalized to sum to one and species mass fractions to `sum(Y) = 1, Y >= 0`. A block that survives a regrid
+  keeps its own fine solution.
+- **Advance.** Lock-step (default): every level takes each step at `dt`, so `dt` must satisfy the finest level's CFL
+  limit. Subcycled (`amr_subcycle = T`): each refined level takes two steps per parent step, with ghost cells
+  interpolated in time; `dt` is the level-0 step.
+- **Refluxing and restriction.** Fine face fluxes (advective, viscous, diffusive) are accumulated in flux registers and
+  replace the coarse fluxes at coarse/fine faces, so mass, momentum and energy are conserved across the interface. Each
+  covered coarse cell is then set to the average of its fine children.
+
+## Conservation and accuracy {#amr-conservation}
+
+Conservation across coarse/fine faces holds to round-off for single- and multi-fluid, viscous and chemistry cases, and
+a uniform state stays uniform with regrid on. Interiors keep the WENO order; the coarse/fine ghost fill is second
+order (linear) or higher (quartic). Results do not depend on message timing, and agree across rank counts to round-off
+when `amr_max_grid_size` is set (with `amr_max_grid_size = 0` the block size cap, and therefore the box set, depends
+on the decomposition).
+
+## Parallelism {#amr-parallel}
+
+Each block is owned whole by one rank. Owners come from a cost-weighted space-filling-curve cut over the blocks (cost:
+one per coarse cell plus `K_ib` per immersed-boundary cell), recomputed at every regrid. `amr_place_owner = T` instead
+places each block on the rank that owns most of the coarse cells under it, capped at 5% above the mean load, which
+cuts communication. A block may straddle rank boundaries; its coarse-side coupling moves point-to-point between the
+owner and the ranks it overlaps. All exchanges are posted in a fixed order, so runs are deterministic.
+
+On GPUs every block array is device-resident. A rank advances its blocks of equal level in batches stacked along the
+last dimension, one solver call per batch; this is why AMR needs a uniform Cartesian grid. The toolchain enables
+`amr_device_pack` (fused exchange kernels) when `amr_max_grid_size <= 64` and sets `amr_snap = min(2, amr_buf - 2)`
+(regrid hysteresis) under dynamic regrid, unless they are set explicitly.
+
+`rank_time_wrt = T` prints a per-phase time budget of the AMR step and the block-balance report.
+
+## Supported physics {#amr-physics}
+
+Required: `model_eqns = 2` or `3`, `time_stepper = 3`, WENO (`recon_type = 1`) or IGR, a uniform Cartesian grid, and
+`mpp_lim = T` when `num_fluids > 1`.
+
+| Feature | Status |
+| :--- | :--- |
+| Single- and multi-fluid Euler, 6-equation model | Supported |
+| Viscous, hypoelastic (with continuum damage), body forces | Supported |
+| Chemistry (reactions, diffusion) | Lock-step only; no reaction substepping |
+| Immersed boundaries (fixed, prescribed or force-driven motion; not STL) | Supported; blocks always contain each body; moving bodies are lock-step only; multi-level with IB is single-rank |
+| Lagrangian bubbles | Lock-step; the cloud stays on the base grid (blocks are kept clear of it) |
+| Acoustic sources | On the base grid; blocks are kept clear of the source support |
+| MHD / RMHD | 1D only |
+| IGR | Single-rank, lock-step; coarse/fine conservation is truncation-order |
+| Euler-Euler bubbles, QBMM, phase change (`relax`), surface tension, synthetic turbulence | Not supported |
+| Cylindrical coordinates, stretched grids, `bc = -4` | Not supported |
+
+## Parameters {#amr-parameters}
+
+Full descriptions: @ref case.
+
+| Parameter | Default | Description |
+| :--- | :---: | :--- |
+| `amr` | F | Enable AMR |
+| `amr_block_beg(i)`, `amr_block_end(i)` | 0 | Initial block, level-0 cell indices (inclusive) |
+| `amr_regrid_int` | 0 | Steps between regrids; 0 keeps the initial block |
+| `amr_tag_eps` | 0.1 | Tagging threshold on the normalized density jump (> 0 with regrid) |
+| `amr_tag_abs` | F | Tag on the absolute density jump; `amr_tag_eps` is then a density |
+| `amr_buf` | 3 | Padding around tagged cells, coarse cells (>= 1 with regrid) |
+| `amr_cluster_eff` | 0.7 | Minimum tag efficiency of a clustered box |
+| `amr_blocking_factor` | 4 | Minimum box width, coarse cells |
+| `amr_max_grid_size` | 0 | Maximum box width, coarse cells; 0 derives it from the decomposition |
+| `amr_max_blocks` | 1024 | Maximum number of blocks |
+| `amr_max_level` | 1 | Number of refined levels |
+| `amr_ref_ratio` | 2 | Refinement ratio, 2 or 4 (4: one level only) |
+| `amr_subcycle` | F | Subcycle refined levels in time (ratio 2, up to two levels) |
+| `amr_prolong` | 1 | 1: limited linear, 2: quartic (ratio 2) |
+| `amr_place_owner` | F | Owner-aware block placement |
+| `amr_snap` | 0 | Regrid hysteresis, coarse cells (<= `amr_buf - 2`) |
+| `amr_device_pack` | F | Fused device pack/unpack of the coarse-patch exchange |
+
+## Example {#amr-example}
+
+```python
+# 2D shock tracking: one refined level, regrid every 5 steps
+"amr": "T",
+"amr_block_beg(1)": 48, "amr_block_end(1)": 79,
+"amr_block_beg(2)": 48, "amr_block_end(2)": 79,
+"amr_regrid_int": 5,
+"amr_tag_eps": 0.05,
+"amr_buf": 3,
+"amr_max_grid_size": 32,
+```
+
+See `examples/2D_amr_droplet` for a complete case.
+
+## Limitations {#amr-limitations}
+
+- **Block size.** A block spans at most about half a rank's subdomain per dimension unless `amr_max_grid_size` is set;
+  larger features are covered by several blocks.
+- **Restart.** `parallel_io` restarts work on any rank count when `amr_max_grid_size` is set; otherwise restart on no
+  more ranks than wrote the checkpoint. The serial format needs the rank count that wrote it.
+  The file header records `amr_ref_ratio`, `amr_max_level` and the storage precision, which must match (the level count may
+  only grow); a missing fine-level file aborts the restart.
+- **Nesting.** A level-2 block must lie inside one level-1 block; the run aborts otherwise (reduce
+  `amr_max_grid_size`).
+- **Regrid.** When no cell is tagged, regrid keeps the current blocks, so refinement does not retract until something is
+  tagged again. Level-2 and deeper blocks are re-prolonged from their parent at every regrid (only level-1 data is
+  carried over), so frequent regrids cost accuracy at those levels.
+- **Output.** Field output is at level-0 resolution with the fine solution restricted into it; post_process with
+  `amr = T` overlays the blocks as separate SILO domains.