IGR integration

src/common/include/3dHardcodedIC.fpp

@@ -134,6 +134,7 @@
             q_prim_vf(eqn_idx%cont%beg)%sf(i, j, k) = (rho_th - rho_am)*rcut*xcut + rho_am
         else
             q_prim_vf(eqn_idx%adv%beg)%sf(i, j, k) = (1._wp - 2._wp*eps)*rcut*xcut + eps
+            q_prim_vf(eqn_idx%adv%end)%sf(i, j, k) = 1._wp - q_prim_vf(eqn_idx%adv%beg)%sf(i, j, k)
             q_prim_vf(eqn_idx%cont%beg)%sf(i, j, k) = rho_th*q_prim_vf(eqn_idx%adv%beg)%sf(i, j, k)
             q_prim_vf(eqn_idx%cont%end)%sf(i, j, k) = rho_am*(1._wp - q_prim_vf(eqn_idx%adv%beg)%sf(i, j, k))
         end if
@@ -160,6 +161,7 @@
             q_prim_vf(eqn_idx%cont%beg)%sf(i, j, k) = (rho_th - rho_am)*rcut*xcut + rho_am
         else
             q_prim_vf(eqn_idx%adv%beg)%sf(i, j, k) = (1._wp - 2._wp*eps)*rcut*xcut + eps
+            q_prim_vf(eqn_idx%adv%end)%sf(i, j, k) = 1._wp - q_prim_vf(eqn_idx%adv%beg)%sf(i, j, k)
             q_prim_vf(eqn_idx%cont%beg)%sf(i, j, k) = rho_th*q_prim_vf(eqn_idx%adv%beg)%sf(i, j, k)
             q_prim_vf(eqn_idx%cont%end)%sf(i, j, k) = rho_am*(1._wp - q_prim_vf(eqn_idx%adv%beg)%sf(i, j, k))
         end if

src/common/m_global_parameters_common.fpp

@@ -167,15 +167,10 @@ contains
             eqn_idx%mom%end = eqn_idx%cont%end + num_vels
             eqn_idx%E = eqn_idx%mom%end + 1
 
-            if (igr) then
-                ! IGR: volume fractions after energy (N-1 for N fluids; skipped when num_fluids=1)
-                eqn_idx%adv%beg = eqn_idx%E + 1
-                eqn_idx%adv%end = eqn_idx%E + num_fluids - 1
-            else
-                ! WENO/MUSCL + Riemann tracks a total of (N) volume fractions for N fluids
-                eqn_idx%adv%beg = eqn_idx%E + 1
-                eqn_idx%adv%end = eqn_idx%E + num_fluids
-            end if
+            ! Store all N volume fractions for N fluids. IGR originally used a compact N-1 layout,
+            ! but the shared Riemann path expects the final volume fraction to be available.
+            eqn_idx%adv%beg = eqn_idx%E + 1
+            eqn_idx%adv%end = eqn_idx%E + num_fluids
 
             sys_size = eqn_idx%adv%end
 

src/common/m_variables_conversion.fpp

@@ -820,12 +820,10 @@ contains
                     ! Obtaining the density, specific heat ratio function and the liquid stiffness function, respectively
                     call s_convert_to_mixture_variables(q_prim_vf, j, k, l, rho, gamma, pi_inf, qv, Re_K, G, fluid_pp(:)%G)
 
-                    if (.not. igr .or. num_fluids > 1) then
-                        ! Transferring the advection equation(s) variable(s)
-                        do i = eqn_idx%adv%beg, eqn_idx%adv%end
-                            q_cons_vf(i)%sf(j, k, l) = q_prim_vf(i)%sf(j, k, l)
-                        end do
-                    end if
+                    ! Transferring the advection equation(s) variable(s)
+                    do i = eqn_idx%adv%beg, eqn_idx%adv%end
+                        q_cons_vf(i)%sf(j, k, l) = q_prim_vf(i)%sf(j, k, l)
+                    end do
 
                     if (relativity) then
                         if (n == 0) then

src/pre_process/m_assign_variables.fpp

@@ -258,11 +258,9 @@ contains
 
         call s_convert_to_mixture_variables(q_prim_vf, j, k, l, orig_rho, orig_gamma, orig_pi_inf, orig_qv)
 
-        if (.not. igr .or. num_fluids > 1) then
-            do i = eqn_idx%adv%beg, eqn_idx%adv%end
-                q_prim_vf(i)%sf(j, k, l) = patch_icpp(patch_id)%alpha(i - eqn_idx%E)
-            end do
-        end if
+        do i = eqn_idx%adv%beg, eqn_idx%adv%end
+            q_prim_vf(i)%sf(j, k, l) = patch_icpp(patch_id)%alpha(i - eqn_idx%E)
+        end do
 
         if (mpp_lim .and. bubbles_euler) then
             ! adjust volume fractions, according to modeled gas void fraction
@@ -291,11 +289,9 @@ contains
             end do
         end if
 
-        if (.not. igr .or. num_fluids > 1) then
-            do i = eqn_idx%adv%beg, eqn_idx%adv%end
-                q_prim_vf(i)%sf(j, k, l) = patch_icpp(smooth_patch_id)%alpha(i - eqn_idx%E)
-            end do
-        end if
+        do i = eqn_idx%adv%beg, eqn_idx%adv%end
+            q_prim_vf(i)%sf(j, k, l) = patch_icpp(smooth_patch_id)%alpha(i - eqn_idx%E)
+        end do
 
         if (mpp_lim .and. bubbles_euler) then
             ! adjust volume fractions, according to modeled gas void fraction
@@ -356,11 +352,9 @@ contains
 
         q_prim_vf(eqn_idx%E)%sf(j, k, l) = (eta*patch_icpp(patch_id)%pres + (1._wp - eta)*orig_prim_vf(eqn_idx%E))
 
-        if (.not. igr .or. num_fluids > 1) then
-            do i = eqn_idx%adv%beg, eqn_idx%adv%end
-                q_prim_vf(i)%sf(j, k, l) = eta*patch_icpp(patch_id)%alpha(i - eqn_idx%E) + (1._wp - eta)*orig_prim_vf(i)
-            end do
-        end if
+        do i = eqn_idx%adv%beg, eqn_idx%adv%end
+            q_prim_vf(i)%sf(j, k, l) = eta*patch_icpp(patch_id)%alpha(i - eqn_idx%E) + (1._wp - eta)*orig_prim_vf(i)
+        end do
 
         if (mhd) then
             if (n == 0) then  ! 1D: By, Bz