Add physics constraint validation and auto-generated docs

.github/workflows/lint-toolchain.yml

@@ -42,5 +42,20 @@ jobs:
 
     - name: Validate example cases
       run: |
-        ./mfc.sh validate examples/1D_sodshocktube/case.py
-        ./mfc.sh validate examples/2D_shockbubble/case.py
+        failed=0
+        passed=0
+        for case in examples/*/case.py; do
+          output=$(./mfc.sh validate "$case" 2>&1)
+          if [ $? -eq 0 ]; then
+            passed=$((passed + 1))
+          else
+            echo "FAIL: $case"
+            echo "$output"
+            failed=$((failed + 1))
+          fi
+        done
+        echo ""
+        echo "Results: $passed passed, $failed failed"
+        if [ "$failed" -ne 0 ]; then
+          exit 1
+        fi

.gitignore

@@ -44,6 +44,7 @@ docs/*/result*
 docs/documentation/*-example.png
 docs/documentation/examples.md
 docs/documentation/case_constraints.md
+docs/documentation/physics_constraints.md
 
 examples/*batch/*/
 examples/**/D/*

CMakeLists.txt

@@ -660,16 +660,19 @@ if (MFC_DOCUMENTATION)
         VERBATIM
     )
 
-    # Generate case_constraints.md from case_validator.py and examples/
+    # Generate case_constraints.md and physics_constraints.md together.
+    # Both are produced by gen_constraints.sh, so a single command avoids races.
     add_custom_command(
         OUTPUT  "${CMAKE_CURRENT_SOURCE_DIR}/docs/documentation/case_constraints.md"
+                "${CMAKE_CURRENT_SOURCE_DIR}/docs/documentation/physics_constraints.md"
         DEPENDS "${CMAKE_CURRENT_SOURCE_DIR}/toolchain/mfc/gen_case_constraints_docs.py"
+                "${CMAKE_CURRENT_SOURCE_DIR}/toolchain/mfc/gen_physics_docs.py"
                 "${CMAKE_CURRENT_SOURCE_DIR}/toolchain/mfc/case_validator.py"
                 "${CMAKE_CURRENT_SOURCE_DIR}/toolchain/mfc/params/definitions.py"
                 "${examples_DOCs}"
         COMMAND "bash" "${CMAKE_CURRENT_SOURCE_DIR}/docs/gen_constraints.sh"
                        "${CMAKE_CURRENT_SOURCE_DIR}"
-        COMMENT "Generating case_constraints.md"
+        COMMENT "Generating case_constraints.md and physics_constraints.md"
         VERBATIM
     )
 
@@ -732,11 +735,13 @@ if (MFC_DOCUMENTATION)
 
         set(opt_example_dependency "")
         set(opt_constraints_dependency "")
+        set(opt_physics_dependency "")
         set(opt_cli_reference_dependency "")
         set(opt_parameters_dependency "")
         if (${target} STREQUAL documentation)
             set(opt_example_dependency "${CMAKE_CURRENT_SOURCE_DIR}/docs/documentation/examples.md")
             set(opt_constraints_dependency "${CMAKE_CURRENT_SOURCE_DIR}/docs/documentation/case_constraints.md")
+            set(opt_physics_dependency "${CMAKE_CURRENT_SOURCE_DIR}/docs/documentation/physics_constraints.md")
             set(opt_cli_reference_dependency "${CMAKE_CURRENT_SOURCE_DIR}/docs/documentation/cli-reference.md")
             set(opt_parameters_dependency "${CMAKE_CURRENT_SOURCE_DIR}/docs/documentation/parameters.md")
         endif()
@@ -749,6 +754,7 @@ if (MFC_DOCUMENTATION)
             DEPENDS "${CMAKE_CURRENT_BINARY_DIR}/${target}-Doxyfile"
                     "${opt_example_dependency}"
                     "${opt_constraints_dependency}"
+                    "${opt_physics_dependency}"
                     "${opt_cli_reference_dependency}"
                     "${opt_parameters_dependency}"
                     "${${target}_SRCs}" "${${target}_DOCs}"

docs/custom.css

@@ -3,6 +3,22 @@
  * Overrides for doxygen-awesome theme
  */
 
+/* Seamless split <code> tags for Fortran % struct accessors.
+ * Doxygen consumes %<word> even inside code spans, so we split around %
+ * into adjacent <code> elements and remove internal borders/padding. */
+code.f90l {
+    border-right: 0;
+    border-top-right-radius: 0;
+    border-bottom-right-radius: 0;
+    padding-right: 0;
+}
+code.f90r {
+    border-left: 0;
+    border-top-left-radius: 0;
+    border-bottom-left-radius: 0;
+    padding-left: 0;
+}
+
 /* Fix inline code visibility in colored admonition blocks (warning, attention, important, note, etc.) */
 
 /* Warning/Attention/Important blocks (red/pink background) */

docs/documentation/case.md

@@ -65,7 +65,7 @@ To run such a case, use the following format:
 For example, to run the `scaling` case in "weak-scaling" mode:
 
 ```shell
-./mfc.sh run examples/scaling/case.py -t pre_process -j 8 -- --scaling weak
+./mfc.sh run examples/scaling/benchmark.py -t pre_process -j 8 -- --scaling weak
 ```
 
 ## Parameters
@@ -297,7 +297,7 @@ These physical parameters must be consistent with fluid material's parameters de
 
 #### Elliptic Smoothing
 
-Initial conditions in which not all patches support the `patch_icpp(j)%smoothen` parameter can still be smoothed by applying iterations of the heat equation to the initial condition.
+Initial conditions in which not all patches support the `patch_icpp(j)%%smoothen` parameter can still be smoothed by applying iterations of the heat equation to the initial condition.
 This is enabled by adding `'elliptic_smoothing': "T",` and `'elliptic_smoothing_iters': N,` to the case dictionary, where `N` is the number of smoothing iterations to apply.
 
 ### 4. Immersed Boundary Patches {#sec-immersed-boundary-patches}
@@ -620,8 +620,8 @@ To restart the simulation from $k$-th time step, see @ref running "Restarting Ca
 | `num_probes`            | Integer | Number of probes	|
 | `probe(i)%[x,y,z]`      | Real	  | Coordinates of probe $i$	|
 | `output_partial_domain` | Logical | Output part of the domain |
-| `[x,y,z]_output%beg`    | Real    | Beginning of the output domain in the [x,y,z]-direction |
-| `[x,y,z]_output%end`    | Real    | End of the output domain in the [x,y,z]-direction |
+| `[x,y,z]_output%%beg`    | Real    | Beginning of the output domain in the [x,y,z]-direction |
+| `[x,y,z]_output%%end`    | Real    | End of the output domain in the [x,y,z]-direction |
 | `lag_txt_wrt`           | Logical | Write Lagrangian bubble data to `.dat` files |
 | `lag_header`            | Logical | Write header to Lagrangian bubble `.dat` files |
 | `lag_db_wrt`            | Logical | Write Lagrangian bubble data to silo/hdf5 database files |
@@ -667,7 +667,7 @@ If `file_per_process` is true, then pre_process, simulation, and post_process mu
 
 - `probe_wrt` activates the output of state variables at coordinates specified by `probe(i)%[x;y,z]`.
 
-- `output_partial_domain` activates the output of part of the domain specified by `[x,y,z]_output%beg` and `[x,y,z]_output%end`.
+- `output_partial_domain` activates the output of part of the domain specified by `[x,y,z]_output%%beg` and `[x,y,z]_output%%end`.
 This is useful for large domains where only a portion of the domain is of interest.
 It is not supported when `precision = 1` and `format = 1`.
 It also cannot be enabled with `flux_wrt`, `heat_ratio_wrt`, `pres_inf_wrt`, `c_wrt`, `omega_wrt`, `ib`, `schlieren_wrt`, `qm_wrt`, or 'liutex_wrt'.
@@ -991,30 +991,30 @@ Note: For relativistic flow, the conservative and primitive densities are differ
 
 When ``cyl_coord = 'T'`` is set in 3D the following constraints must be met:
 
-- `bc_y%beg = -14`  enables the axis boundary condition
+- `bc_y%%beg = -14`  enables the axis boundary condition
 
-- `bc_z%beg = bc_z%end = -1`  enables periodic boundary conditions in the azimuthal direction
+- `bc_z%%beg = bc_z%%end = -1`  enables periodic boundary conditions in the azimuthal direction
 
-- `z_domain%beg = 0`  sets the azimuthal starting point to 0
+- `z_domain%%beg = 0`  sets the azimuthal starting point to 0
 
-- `z_domain%end = 2*math.pi` to set the azimuthal ending point to \f$2\pi\f$ (note, requires `import math` in the case file)
+- `z_domain%%end = 2*math.pi` to set the azimuthal ending point to \f$2\pi\f$ (note, requires `import math` in the case file)
 
 When ``cyl_coord = 'T'`` is set in 2D the following constraints must be met:
 
-- `bc_y%beg = -2` to enable reflective boundary conditions
+- `bc_y%%beg = -2` to enable reflective boundary conditions
 
 ### 17. Chemistry
 
 | Parameter                     | Type    | Description                                              |
 | ---:                          | :---:   | :---                                                     |
 | `chemistry`                   | Logical | Enable chemistry simulation                              |
-| `chem_params%diffusion`       | Logical | Enable multispecies diffusion                            |
-| `chem_params%reactions`       | Logical | Enable chemical reactions                                |
-| `chem_params%gamma_method`    | Integer | Methodology for calculating the heat capacity ratio      |
-| `chem_params%transport_model` | Integer | Methodology for calculating the diffusion coefficients   |
+| `chem_params%%diffusion`       | Logical | Enable multispecies diffusion                            |
+| `chem_params%%reactions`       | Logical | Enable chemical reactions                                |
+| `chem_params%%gamma_method`    | Integer | Methodology for calculating the heat capacity ratio      |
+| `chem_params%%transport_model` | Integer | Methodology for calculating the diffusion coefficients   |
 | `cantera_file`                | String  | Cantera-format mechanism file (e.g., .yaml)              |
 
-- `chem_params%transport_model` specifies the methodology for calculating diffusion coefficients and other transport properties, `1` for mixture-average, `2` for Unity-Lewis
+- `chem_params%%transport_model` specifies the methodology for calculating diffusion coefficients and other transport properties, `1` for mixture-average, `2` for Unity-Lewis
 
 - `cantera_file` specifies the chemical mechanism file. If the file is part of the standard Cantera library, only the filename is required. Otherwise, the file must be located in the same directory as your `case.py` file
 
@@ -1051,15 +1051,15 @@ The entries labeled "Characteristic." are characteristic boundary conditions bas
 
 | Parameter                     | Type    | Description |
 | ---:                          | :----:  | :--- |
-| `bc_[x,y,z]%grcbc_in`         | Logical | Enable grcbc for subsonic inflow |
-| `bc_[x,y,z]%grcbc_out`        | Logical | Enable grcbc for subsonic outflow (pressure)|
-| `bc_[x,y,z]%grcbc_vel_out`    | Logical | Enable grcbc for subsonic outflow (pressure + normal velocity) |
-| `bc_[x,y,z]%vel_in`           | Real Array | Inflow velocities in x, y and z directions |
-| `bc_[x,y,z]%vel_out`          | Real Array | Outflow velocities in x, y and z directions |
-| `bc_[x,y,z]%pres_in`          | Real    | Inflow pressure |
-| `bc_[x,y,z]%pres_out`         | Real    | Outflow pressure |
-| `bc_[x,y,z]%alpha_rho_in`     | Real Array | Inflow density |
-| `bc_[x,y,z]%alpha_in`         | Real Array | Inflow void fraction |
+| `bc_[x,y,z]%%grcbc_in`         | Logical | Enable grcbc for subsonic inflow |
+| `bc_[x,y,z]%%grcbc_out`        | Logical | Enable grcbc for subsonic outflow (pressure)|
+| `bc_[x,y,z]%%grcbc_vel_out`    | Logical | Enable grcbc for subsonic outflow (pressure + normal velocity) |
+| `bc_[x,y,z]%%vel_in`           | Real Array | Inflow velocities in x, y and z directions |
+| `bc_[x,y,z]%%vel_out`          | Real Array | Outflow velocities in x, y and z directions |
+| `bc_[x,y,z]%%pres_in`          | Real    | Inflow pressure |
+| `bc_[x,y,z]%%pres_out`         | Real    | Outflow pressure |
+| `bc_[x,y,z]%%alpha_rho_in`     | Real Array | Inflow density |
+| `bc_[x,y,z]%%alpha_in`         | Real Array | Inflow void fraction |
 
 This boundary condition can be used for subsonic inflow (`bc_[x,y,z]%[beg,end]` = -7) and subsonic outflow (`bc_[x,y,z]%[beg,end]` = -8) characteristic boundary conditions. These are based on \cite Pirozzoli13. This enables to provide inflow and outflow conditions outside the computational domain.
 

docs/documentation/contributing.md

@@ -64,8 +64,8 @@ q_cons_vf (conservative variables: density, momentum, energy, volume fractions)
 ```
 
 Key data structures (defined in `src/common/m_derived_types.fpp`):
-- `scalar_field` — wraps a 3D `real(stp)` array (`%sf(i,j,k)`)
-- `vector_field` — array of `scalar_field` (`%vf(1:sys_size)`)
+- `scalar_field` — wraps a 3D `real(stp)` array (`%%sf(i,j,k)`)
+- `vector_field` — array of `scalar_field` (`%%vf(1:sys_size)`)
 - `q_cons_vf` / `q_prim_vf` — conservative and primitive state vectors
 
 ### Build Toolchain
@@ -136,7 +136,7 @@ Both human reviewers and AI code reviewers reference this section.
 
 ### Array Bounds and Indexing
 
-- MFC uses **non-unity lower bounds** (e.g., `idwbuff(1)%beg:idwbuff(1)%end` with negative ghost-cell indices). Always verify loop bounds match array declarations.
+- MFC uses **non-unity lower bounds** (e.g., `idwbuff(1)%%beg:idwbuff(1)%%end` with negative ghost-cell indices). Always verify loop bounds match array declarations.
 - **Riemann solver indexing:** Left states at `j`, right states at `j+1`. Off-by-one here corrupts fluxes.
 
 ### Precision and Type Safety
@@ -164,7 +164,7 @@ Both human reviewers and AI code reviewers reference this section.
 - **Pressure formula** must match `model_eqns` value. Model 2/3 (multi-fluid), model 4 (bubbles), MHD, and hypoelastic each use different EOS formulations. Wrong formula = wrong physics.
 - **Conservative-primitive conversion:** Density recovery, kinetic energy, and pressure each have model-specific paths. Verify the correct branch is taken.
 - **Volume fractions** must sum to 1. `alpha_rho_K` must be non-negative. Species mass fractions should be clipped to [0,1].
-- **Boundary conditions:** Periodic BCs must match at both ends (`bc_x%beg` and `bc_x%end`). Cylindrical coordinates have special requirements (`bc_y%beg = -14` for axis in 3D).
+- **Boundary conditions:** Periodic BCs must match at both ends (`bc_x%%beg` and `bc_x%%end`). Cylindrical coordinates have special requirements (`bc_y%%beg = -14` for axis in 3D).
 - **Parameter constraints:** New parameters or physics features must be validated in `toolchain/mfc/case_validator.py`. New features should add corresponding validation.
 
 ### Python Toolchain
@@ -274,6 +274,8 @@ def check_my_feature(self):
         self.errors.append("my_param requires other_param to be set")
 ```
 
+If your check enforces a physics constraint, also add a `PHYSICS_DOCS` entry (see [How to Document Physics Constraints](#how-to-document-physics-constraints) below).
+
 **Step 5: Declare in Fortran** (`src/<target>/m_global_parameters.fpp`)
 
 Add the variable declaration in the appropriate target's global parameters module. Choose the target(s) where the parameter is used:
@@ -659,6 +661,41 @@ Checklist:
 
 See @ref troubleshooting for debugging workflows, profiling tools, GPU diagnostic environment variables, common build/runtime errors, and fixes.
 
+### How to Document Physics Constraints {#how-to-document-physics-constraints}
+
+When adding a new `check_` method to `case_validator.py`, document its physics by adding an entry to the `PHYSICS_DOCS` dict at the top of the file:
+
+```python
+PHYSICS_DOCS = {
+    ...
+    "check_my_feature": {
+        "title": "My Feature Constraint",          # Required: human-readable title
+        "category": "Thermodynamic Constraints",    # Required: groups the constraint in docs
+        "explanation": "Why this constraint exists.", # Required: plain English
+        "math": r"\alpha > 0",                      # Optional: LaTeX formula
+        "references": ["Wilfong26"],                # Optional: BibTeX keys from references.bib
+        "exceptions": ["IBM cases"],                # Optional: when constraint doesn't apply
+    },
+}
+```
+
+The @ref physics_constraints "Physics Constraints" page is **auto-generated** — run `./mfc.sh generate` to rebuild it.
+The generator merges your `PHYSICS_DOCS` entry with the AST-extracted `prohibit()`/`warn()` calls,
+so stage, severity, and parameter information appear automatically.
+
+**Fields:**
+
+| Field | Required | Description |
+|-------|----------|-------------|
+| `title` | Yes | Section heading in generated docs |
+| `category` | Yes | Grouping category (e.g., "Mixture Constraints") |
+| `explanation` | Yes | Plain English description of the physics |
+| `math` | No | LaTeX formula (rendered by Doxygen's MathJax) |
+| `references` | No | List of BibTeX cite keys from `docs/references.bib` |
+| `exceptions` | No | List of cases where the constraint doesn't apply |
+
+**Categories:** Thermodynamic Constraints, Mixture Constraints, Domain and Geometry, Velocity and Dimensional Consistency, Model Equations, Boundary Conditions, Bubble Physics, Feature Compatibility, Numerical Schemes, Acoustic Sources, Post-Processing.
+
 ## Testing
 
 MFC has 500+ regression tests. See @ref testing for the full guide.

docs/documentation/equations.md

@@ -777,7 +777,7 @@ Used for viscous fluxes and velocity gradients.
 | `-15` | Slip wall |
 | `-16` | No-slip wall |
 
-### 16.2 Characteristic BCs (\cite Thompson87, \cite Thompson90; `bc_x%beg = -5` to `-12`)
+### 16.2 Characteristic BCs (\cite Thompson87, \cite Thompson90; `bc_x%%beg = -5` to `-12`)
 
 **Characteristic decomposition:**
 

docs/documentation/readme.md

@@ -15,6 +15,7 @@ Welcome to the Multi-component Flow Code (MFC) documentation.
 - @ref parameters "Case Parameters" - All ~3,400 parameters
 - @ref cli-reference "CLI Reference" - Command line options
 - @ref case_constraints "Case Creator Guide" - Feature compatibility
+- @ref physics_constraints "Physics Constraints" - Mathematical basis for validation rules
 
 ## Examples & Visualization
 

docs/documentation/running.md

@@ -309,8 +309,8 @@ If you want to restart a simulation,
 1. For the Computational Domain Parameters
     - Have the following removed __except__ `m`, `n`, and `p`:
 		- All domain/mesh information
-			- `(xyz)_domain%beg`
-			- `(xyz)_domain%end`
+			- `(xyz)_domain%%beg`
+			- `(xyz)_domain%%end`
 			- `stretch_(xyz)`
 			- `a_(xyz)`
 			- `(xyz)_a`
@@ -334,11 +334,11 @@ If you want to restart a simulation,
 
 3. For Patches
 	- Have all information about old patches (used in the `case.py` file) removed.
-		- `patch_icpp(1)%all variables`
-		- `patch_icpp(2)%all variables`
-		- `patch_icpp(num_patches)%all variables`
+		- `patch_icpp(1)%%all variables`
+		- `patch_icpp(2)%%all variables`
+		- `patch_icpp(num_patches)%%all variables`
 	- Add information about new patches that will be introduced, if any. The parameter num_patches should reflect this addition.
-		- e.g. `patch_icpp(1)%some variables of interest`
+		- e.g. `patch_icpp(1)%%some variables of interest`
 
 4. For Fluid properties
 	- Keep information about the fluid properties

docs/gen_constraints.sh

@@ -14,4 +14,5 @@ echo "Generating case constraints documentation..."
 python3 "$REPO_ROOT/toolchain/mfc/gen_case_constraints_docs.py" > "$REPO_ROOT/docs/documentation/case_constraints.md"
 echo "✓ Generated docs/documentation/case_constraints.md"
 
-
+python3 "$REPO_ROOT/toolchain/mfc/gen_physics_docs.py" > "$REPO_ROOT/docs/documentation/physics_constraints.md"
+echo "✓ Generated docs/documentation/physics_constraints.md"

examples/2D_bubbly_steady_shock/case.py

@@ -155,6 +155,7 @@
             "fluid_pp(1)%pi_inf": gam_l * (pi_inf_l) / (gam_l - 1.0),
             # Bubbles
             "bubbles_euler": "T",
+            "nb": nb,
             "bubble_model": 2,
             "polytropic": "T",
             "polydisperse": "F",