Reject chem_params%diffusion without chemistry

src/simulation/m_rhs.fpp

@@ -2242,7 +2242,7 @@ contains
                         end do
                     end if
 
-                    if (chem_params%diffusion .and. .not. (viscous .or. surface_tension .or. heat_conduction)) then
+                    if (chemistry .and. chem_params%diffusion .and. .not. (viscous .or. surface_tension .or. heat_conduction)) then
                         @:DEALLOCATE(flux_src_n(i)%vf(eqn_idx%E)%sf)
                     end if
 

toolchain/mfc/case_validator.py

@@ -2156,6 +2156,10 @@ def check_chemistry(self):
         # inconsistent and the simulation NaNs. See MFlowCode/MFC#1470.
         self.prohibit(chemistry and num_fluids is not None and num_fluids != 1, "chemistry is only supported for single-component flows (num_fluids = 1)")
 
+        # The diffusion flux arrays are allocated only under chemistry, but the Riemann-state and
+        # RHS loops test chem_params%diffusion alone, so without chemistry they write unallocated arrays.
+        self.prohibit(diffusion and not chemistry, "chem_params%diffusion requires chemistry = T")
+
         # Chemistry with Euler bubbles is not currently supported: the IBM image-point
         # interpolation branch selects the bubbles/QBMM path before the chemistry path, so
         # the species state is not carried when both are enabled. Disallow the combination

toolchain/mfc/test_case_validator.py

@@ -168,6 +168,14 @@ def test_accepts_2d(self):
         self.assertAccepts({**BASE_2D, "bf_spatial_support": "T"})
 
 
+class TestChemistryDiffusion(ConstraintTestCase):
+    def test_rejects_diffusion_without_chemistry(self):
+        self.assertRejects({**BASE, "chem_params%diffusion": "T"}, "chem_params%diffusion requires chemistry")
+
+    def test_accepts_diffusion_with_chemistry(self):
+        self.assertAccepts({**CHEMISTRY, "chem_params%diffusion": "T"})
+
+
 class TestChemistrySubstepping(ConstraintTestCase):
     def test_rejects_negative_substeps(self):
         self.assertRejects({**CHEMISTRY, "chem_params%reaction_substeps": -1}, "reaction_substeps must be >= 0")