Fix three bugs flagged in code review

1. check_non_newtonian() was never called from validate_simulation(),
   so all non-Newtonian parameter constraints were silently skipped.
   Also add required-field checks for K, nn, hb_m and an IGR prohibition.

2. Re_visc_nn was missing from private(...) in all four GPU_PARALLEL_LOOP
   calls in m_viscous.fpp, causing thread race conditions on GPU builds.

3. igr + non_newtonian is now prohibited at validation time: q_prim_vf
   sub-arrays are only allocated in the non-IGR path, so s_compute_dt
   would dereference unallocated storage with any_non_newtonian=.true.

Author: sbryngelson

src/simulation/m_viscous.fpp

@@ -74,7 +74,7 @@ contains
         #:if not MFC_CASE_OPTIMIZATION or num_dims > 1
             if (shear_stress) then  ! Shear stresses
                 $:GPU_PARALLEL_LOOP(collapse=3, private='[i, j, k, l, rho_visc, gamma_visc, pi_inf_visc, alpha_visc_sum , &
-                                    & alpha_visc, alpha_rho_visc, Re_visc, tau_Re]')
+                                    & alpha_visc, alpha_rho_visc, Re_visc, Re_visc_nn, tau_Re]')
                 do l = is3_viscous%beg, is3_viscous%end
                     do k = -1, 1
                         do j = is1_viscous%beg, is1_viscous%end
@@ -168,7 +168,7 @@ contains
         #:if not MFC_CASE_OPTIMIZATION or num_dims > 1
             if (bulk_stress) then  ! Bulk stresses
                 $:GPU_PARALLEL_LOOP(collapse=3, private='[i, j, k, l, rho_visc, gamma_visc, pi_inf_visc, alpha_visc_sum , &
-                                    & alpha_visc, alpha_rho_visc, Re_visc, tau_Re]')
+                                    & alpha_visc, alpha_rho_visc, Re_visc, Re_visc_nn, tau_Re]')
                 do l = is3_viscous%beg, is3_viscous%end
                     do k = -1, 1
                         do j = is1_viscous%beg, is1_viscous%end
@@ -255,7 +255,7 @@ contains
         if (p == 0) return
         #:if not MFC_CASE_OPTIMIZATION or num_dims > 2
             if (shear_stress) then  ! Shear stresses
-                $:GPU_PARALLEL_LOOP(collapse=3, private='[alpha_visc, alpha_rho_visc, Re_visc, tau_Re]')
+                $:GPU_PARALLEL_LOOP(collapse=3, private='[alpha_visc, alpha_rho_visc, Re_visc, Re_visc_nn, tau_Re]')
                 do l = is3_viscous%beg, is3_viscous%end
                     do k = -1, 1
                         do j = is1_viscous%beg, is1_viscous%end
@@ -345,7 +345,7 @@ contains
             end if
 
             if (bulk_stress) then  ! Bulk stresses
-                $:GPU_PARALLEL_LOOP(collapse=3, private='[alpha_visc, alpha_rho_visc, Re_visc, tau_Re]')
+                $:GPU_PARALLEL_LOOP(collapse=3, private='[alpha_visc, alpha_rho_visc, Re_visc, Re_visc_nn, tau_Re]')
                 do l = is3_viscous%beg, is3_viscous%end
                     do k = -1, 1
                         do j = is1_viscous%beg, is1_viscous%end

toolchain/mfc/case_validator.py

@@ -883,8 +883,12 @@ def check_viscosity(self):
     def check_non_newtonian(self):
         """Checks constraints on non-Newtonian fluid parameters"""
         viscous = self.get("viscous", "F") == "T"
+        igr = self.get("igr", "F") == "T"
         num_fluids = self.get("num_fluids", 1)
 
+        any_nn = any(self.get(f"fluid_pp({i})%non_newtonian", "F") == "T" for i in range(1, num_fluids + 1))
+        self.prohibit(any_nn and igr, "non_newtonian is not supported with igr (primitive variable storage is not allocated in IGR mode)")
+
         for i in range(1, num_fluids + 1):
             nn_flag = self.get(f"fluid_pp({i})%non_newtonian", "F") == "T"
             if not nn_flag:
@@ -899,6 +903,10 @@ def check_non_newtonian(self):
             mu_max = self.get(f"fluid_pp({i})%mu_max")
             hb_m = self.get(f"fluid_pp({i})%hb_m")
 
+            self.prohibit(K is None, f"fluid_pp({i})%K must be set when non_newtonian = T")
+            self.prohibit(nn is None, f"fluid_pp({i})%nn must be set when non_newtonian = T")
+            self.prohibit(hb_m is None, f"fluid_pp({i})%hb_m must be set when non_newtonian = T")
+
             self.prohibit(K is not None and K <= 0, f"fluid_pp({i})%K (consistency index) must be > 0")
             self.prohibit(nn is not None and nn <= 0, f"fluid_pp({i})%nn (flow behavior index) must be > 0")
             self.prohibit(tau0 is not None and tau0 < 0, f"fluid_pp({i})%tau0 (yield stress) must be >= 0")
@@ -2111,6 +2119,7 @@ def validate_simulation(self):
         self.check_bubbles_euler_simulation()
         self.check_body_forces()
         self.check_viscosity()
+        self.check_non_newtonian()
         self.check_mhd_simulation()
         self.check_igr_simulation()
         self.check_acoustic_source()