Merge branch 'master' into feature/non-newtonian-viscosity

Author: sbryngelson

src/simulation/m_ibm.fpp

@@ -193,7 +193,7 @@ contains
         type(ghost_point) :: gp
         type(ghost_point) :: innerp
 
-        ! set the Moving IBM interior Pressure Values
+        ! set the Moving IBM interior conservative variables
         $:GPU_PARALLEL_LOOP(private='[i,j,k,patch_id,rho]', copyin='[E_idx,momxb]', collapse=3)
         do l = 0, p
             do k = 0, n

toolchain/mfc/case.py

@@ -62,7 +62,7 @@ def get_inp(self, _target) -> str:
         cons.print(f"Generating [magenta]{target.name}.inp[/magenta]:")
         cons.indent()
 
-        MASTER_KEYS: list = case_dicts.get_input_dict_keys(target.name)
+        MASTER_KEYS = case_dicts.get_input_dict_keys(target.name)
 
         ignored = []
 

toolchain/mfc/case_validator.py

@@ -19,6 +19,7 @@
 
 from .common import MFCException
 from .params.definitions import CONSTRAINTS
+from .params.namelist_parser import get_fortran_constants
 from .state import CFG
 
 # Physics documentation for check methods.
@@ -572,7 +573,12 @@ def check_ibm(self):
         num_ibs = self.get("num_ibs", 0)
 
         self.prohibit(ib and n <= 0, "Immersed Boundaries do not work in 1D (requires n > 0)")
-        self.prohibit(ib and (num_ibs <= 0 or num_ibs > 1000), "num_ibs must be between 1 and num_patches_max (1000)")
+        self.prohibit(ib and num_ibs <= 0, "num_ibs must be >= 1 when ib is enabled")
+        num_patches_max = get_fortran_constants().get("num_patches_max", 1000)
+        self.prohibit(
+            ib and num_ibs > num_patches_max,
+            f"num_ibs must be <= {num_patches_max} (num_patches_max in m_constants.fpp)",
+        )
         self.prohibit(not ib and num_ibs > 0, "num_ibs is set, but ib is not enabled")
 
     def check_stiffened_eos(self):
@@ -1211,6 +1217,11 @@ def check_restart(self):
         self.prohibit(old_grid and t_step_old is None, "old_grid requires t_step_old to be set")
         self.prohibit(num_patches < 0, "num_patches must be non-negative")
         self.prohibit(num_patches == 0 and t_step_old is None, "num_patches must be positive for the non-restart case")
+        num_patches_max = get_fortran_constants().get("num_patches_max", 1000)
+        self.prohibit(
+            num_patches > num_patches_max,
+            f"num_patches must be <= {num_patches_max} (num_patches_max in m_constants.fpp)",
+        )
 
     def check_qbmm_pre_process(self):
         """Checks QBMM constraints for pre-process"""
@@ -1422,7 +1433,7 @@ def check_patch_physics(self):
     def check_bc_patches(self):
         """Checks boundary condition patch geometry (pre-process)"""
         num_bc_patches = self.get("num_bc_patches", 0)
-        num_bc_patches_max = self.get("num_bc_patches_max", 10)
+        num_bc_patches_max = get_fortran_constants().get("num_bc_patches_max", 10)
 
         if num_bc_patches <= 0:
             return

toolchain/mfc/params/__init__.py

@@ -23,11 +23,12 @@
 # and freezes it. It must come after REGISTRY is imported and must not be removed.
 from . import definitions  # noqa: F401
 from .definitions import CONSTRAINTS, DEPENDENCIES, get_value_label
-from .registry import REGISTRY, RegistryFrozenError
+from .registry import REGISTRY, IndexedFamily, RegistryFrozenError
 from .schema import ParamDef, ParamType
 
 __all__ = [
     "REGISTRY",
+    "IndexedFamily",
     "RegistryFrozenError",
     "ParamDef",
     "ParamType",

toolchain/mfc/params/definitions.py

@@ -8,11 +8,34 @@
 import re
 from typing import Any, Dict
 
-from .registry import REGISTRY
+from .namelist_parser import get_fortran_constants
+from .registry import REGISTRY, IndexedFamily
 from .schema import ParamDef, ParamType
 
-# Index limits
-NP, NF, NI, NA, NPR, NB = 10, 10, 1000, 4, 10, 10  # patches, fluids, ibs, acoustic, probes, bc_patches
+# Index limits — sourced from Fortran compile-time constants (m_constants.fpp).
+# These must stay in sync with Fortran; we error if the source can't be parsed.
+_FC = get_fortran_constants()
+
+
+def _fc(name: str) -> int:
+    """Get a required Fortran constant, raising if unavailable."""
+    if name not in _FC:
+        raise RuntimeError(
+            f"Fortran constant '{name}' not found in m_constants.fpp. "
+            f"Toolchain is out of sync with Fortran source."
+        )
+    return _FC[name]
+
+
+NF = _fc("num_fluids_max")  # fluid_pp
+NPR = _fc("num_probes_max")  # probe, acoustic, integral
+NB = _fc("num_bc_patches_max")  # patch_bc
+NUM_PATCHES_MAX = _fc("num_patches_max")  # patch_icpp, patch_ib (Fortran array bound)
+# Enumeration limits for families not yet converted to IndexedFamily.
+# These are smaller than the Fortran array bounds to keep the registry compact.
+# The CONSTRAINTS dict below uses the Fortran constants for validation.
+NP = 10  # patch_icpp: has per-index variations, can't easily be IndexedFamily
+NA = 4  # acoustic sources: enumerated individually
 
 
 # Auto-generated Descriptions
@@ -637,9 +660,9 @@ def get_value_label(param_name: str, value: int) -> str:
     "R0ref": {"min": 0},
     "sigma": {"min": 0},
     # Counts (must be positive)
-    "num_fluids": {"min": 1, "max": 10},
-    "num_patches": {"min": 0, "max": 10},
-    "num_ibs": {"min": 0, "max": 1000},
+    "num_fluids": {"min": 1, "max": NF},
+    "num_patches": {"min": 0, "max": NUM_PATCHES_MAX},
+    "num_ibs": {"min": 0},
     "num_source": {"min": 1},
     "num_probes": {"min": 1},
     "num_integrals": {"min": 1},
@@ -1139,26 +1162,37 @@ def _load():
     ]:
         _r(f"bub_pp%{a}", REAL, {"bubbles"}, math=sym)
 
-    # patch_ib (10 immersed boundaries)
-    for i in range(1, NI + 1):
-        px = f"patch_ib({i})%"
-        for a in ["geometry", "moving_ibm"]:
-            _r(f"{px}{a}", INT, {"ib"})
-        for a, pt in [("radius", REAL), ("theta", REAL), ("slip", LOG), ("c", REAL), ("p", REAL), ("t", REAL), ("m", REAL), ("mass", REAL)]:
-            _r(f"{px}{a}", pt, {"ib"})
-        for j in range(1, 4):
-            _r(f"{px}angles({j})", REAL, {"ib"})
-        for d in ["x", "y", "z"]:
-            _r(f"{px}{d}_centroid", REAL, {"ib"})
-            _r(f"{px}length_{d}", REAL, {"ib"})
-        for a, pt in [("model_filepath", STR), ("model_spc", INT), ("model_threshold", REAL)]:
-            _r(f"{px}{a}", pt, {"ib"})
-        for t in ["translate", "scale", "rotate"]:
-            for j in range(1, 4):
-                _r(f"{px}model_{t}({j})", REAL, {"ib"})
+    # patch_ib (immersed boundaries) — registered as indexed family for O(1) lookup.
+    # max_index is None so the parameter registry stays compact (no enumeration).
+    # The Fortran-side upper bound (num_patches_max in m_constants.fpp) is parsed
+    # and enforced by the case_validator, not by max_index here.
+    _ib_tags = {"ib"}
+    _ib_attrs: Dict[str, tuple] = {}
+    for a in ["geometry", "moving_ibm"]:
+        _ib_attrs[a] = (INT, _ib_tags)
+    for a, pt in [("radius", REAL), ("theta", REAL), ("slip", LOG), ("c", REAL), ("p", REAL), ("t", REAL), ("m", REAL), ("mass", REAL)]:
+        _ib_attrs[a] = (pt, _ib_tags)
+    for j in range(1, 4):
+        _ib_attrs[f"angles({j})"] = (REAL, _ib_tags)
+    for d in ["x", "y", "z"]:
+        _ib_attrs[f"{d}_centroid"] = (REAL, _ib_tags)
+        _ib_attrs[f"length_{d}"] = (REAL, _ib_tags)
+    for a, pt in [("model_filepath", STR), ("model_spc", INT), ("model_threshold", REAL)]:
+        _ib_attrs[a] = (pt, _ib_tags)
+    for t in ["translate", "scale", "rotate"]:
         for j in range(1, 4):
-            _r(f"{px}vel({j})", A_REAL, {"ib"})
-            _r(f"{px}angular_vel({j})", A_REAL, {"ib"})
+            _ib_attrs[f"model_{t}({j})"] = (REAL, _ib_tags)
+    for j in range(1, 4):
+        _ib_attrs[f"vel({j})"] = (A_REAL, _ib_tags)
+        _ib_attrs[f"angular_vel({j})"] = (A_REAL, _ib_tags)
+    REGISTRY.register_family(
+        IndexedFamily(
+            base_name="patch_ib",
+            attrs=_ib_attrs,
+            tags=_ib_tags,
+            max_index=NUM_PATCHES_MAX,
+        )
+    )
 
     # acoustic sources (4 sources)
     for i in range(1, NA + 1):

toolchain/mfc/params/namelist_parser.py

@@ -11,7 +11,7 @@
 
 import re
 from pathlib import Path
-from typing import Dict, Set
+from typing import Dict, Optional, Set
 
 # Fallback parameters for when Fortran source files are not available.
 # Generated from the namelist definitions in src/*/m_start_up.fpp.
@@ -464,6 +464,44 @@ def parse_all_namelists(mfc_root: Path) -> Dict[str, Set[str]]:
     return result
 
 
+def parse_fortran_constants(filepath: Path) -> Dict[str, int]:
+    """
+    Parse integer parameter constants from a Fortran source file.
+
+    Extracts lines like ``integer, parameter :: name = 123`` and returns
+    a dict mapping constant names to their integer values.
+    """
+    constants: Dict[str, int] = {}
+    pattern = re.compile(
+        r"integer\s*,\s*parameter\s*::\s*(\w+)\s*=\s*(\d+)", re.IGNORECASE
+    )
+    try:
+        text = filepath.read_text()
+    except FileNotFoundError:
+        return constants
+    for m in pattern.finditer(text):
+        constants[m.group(1)] = int(m.group(2))
+    return constants
+
+
+# Module-level cache for Fortran constants (None = not yet loaded)
+_FORTRAN_CONSTANTS_CACHE: Optional[Dict[str, int]] = None
+
+
+def get_fortran_constants() -> Dict[str, int]:
+    """
+    Get Fortran compile-time constants from m_constants.fpp.
+
+    Cached after first call. Returns empty dict if source unavailable.
+    """
+    global _FORTRAN_CONSTANTS_CACHE  # noqa: PLW0603
+    if _FORTRAN_CONSTANTS_CACHE is None:
+        root = get_mfc_root()
+        path = root / "src" / "common" / "m_constants.fpp"
+        _FORTRAN_CONSTANTS_CACHE = parse_fortran_constants(path)
+    return _FORTRAN_CONSTANTS_CACHE
+
+
 def get_mfc_root() -> Path:
     """Get the MFC root directory from this file's location."""
     # This file is at toolchain/mfc/params/namelist_parser.py