Fix chemistry‑related bugs (#829)

Co-authored-by: Spencer Bryngelson <[email protected]>
Co-authored-by: Spencer Bryngelson <[email protected]>

Author: DimAdam-01

src/common/m_chemistry.fpp

@@ -9,7 +9,7 @@
 module m_chemistry
 
     use m_thermochem, only: &
-        num_species, mol_weights, get_temperature, get_net_production_rates
+        num_species, molecular_weights, get_temperature, get_net_production_rates
 
     use m_global_parameters
 
@@ -91,7 +91,7 @@ contains
                     !$acc loop seq
                     do eqn = chemxb, chemxe
 
-                        omega_m = mol_weights(eqn - chemxb + 1)*omega(eqn - chemxb + 1)
+                        omega_m = molecular_weights(eqn - chemxb + 1)*omega(eqn - chemxb + 1)
 
                         rhs_vf(eqn)%sf(x, y, z) = rhs_vf(eqn)%sf(x, y, z) + omega_m
 

src/post_process/m_mpi_proxy.fpp

@@ -165,7 +165,7 @@ contains
             & 'alt_soundspeed', 'hypoelasticity', 'mhd', 'parallel_io',        &
             & 'rho_wrt', 'E_wrt', 'pres_wrt', 'gamma_wrt', 'sim_data',         &
             & 'heat_ratio_wrt', 'pi_inf_wrt', 'pres_inf_wrt', 'cons_vars_wrt', &
-            & 'prim_vars_wrt', 'c_wrt', 'qm_wrt','schlieren_wrt',              &
+            & 'prim_vars_wrt', 'c_wrt', 'qm_wrt','schlieren_wrt','chem_wrt_T', &
             & 'bubbles_euler', 'qbmm', 'polytropic', 'polydisperse',           &
             & 'file_per_process', 'relax', 'cf_wrt',                           &
             & 'adv_n', 'ib', 'cfl_adap_dt', 'cfl_const_dt', 'cfl_dt',          &

src/pre_process/m_mpi_proxy.fpp

@@ -128,6 +128,9 @@ contains
 
             call MPI_BCAST(patch_icpp(i)%model_spc, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
 
+            if (chemistry) then
+                call MPI_BCAST(patch_icpp(i)%Y, size(patch_icpp(i)%Y), mpi_p, 0, MPI_COMM_WORLD, ierr)
+            end if
             ! Broadcast IB variables
             call MPI_BCAST(patch_ib(i)%geometry, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
             call MPI_BCAST(patch_ib(i)%model_filepath, len(patch_ib(i)%model_filepath), MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr)

src/simulation/include/inline_riemann.fpp

@@ -37,8 +37,8 @@
         call get_species_enthalpies_rt(T_L, h_iL)
         call get_species_enthalpies_rt(T_R, h_iR)
 
-        h_iL = h_iL*gas_constant/mol_weights*T_L
-        h_iR = h_iR*gas_constant/mol_weights*T_R
+        h_iL = h_iL*gas_constant/molecular_weights*T_L
+        h_iR = h_iR*gas_constant/molecular_weights*T_R
         call get_species_specific_heats_r(T_L, Cp_iL)
         call get_species_specific_heats_r(T_R, Cp_iR)
 
@@ -48,17 +48,17 @@
 
         if (abs(T_L - T_R) < eps) then
             ! Case when T_L and T_R are very close
-            Cp_avg = sum(Yi_avg(:)*(0.5_wp*Cp_iL(i) + 0.5_wp*Cp_iR(:))*gas_constant/mol_weights(:))
-            Cv_avg = sum(Yi_avg(:)*((0.5_wp*Cp_iL(i) + 0.5_wp*Cp_iR(:))*gas_constant/mol_weights(:) - gas_constant/mol_weights(:)))
+            Cp_avg = sum(Yi_avg(:)*(0.5_wp*Cp_iL(:) + 0.5_wp*Cp_iR(:))*gas_constant/molecular_weights(:))
+            Cv_avg = sum(Yi_avg(:)*((0.5_wp*Cp_iL(:) + 0.5_wp*Cp_iR(:))*gas_constant/molecular_weights(:) - gas_constant/molecular_weights(:)))
         else
             ! Normal calculation when T_L and T_R are sufficiently different
             Cp_avg = sum(Yi_avg(:)*(h_iR(:) - h_iL(:))/(T_R - T_L))
-            Cv_avg = sum(Yi_avg(:)*((h_iR(:) - h_iL(:))/(T_R - T_L) - gas_constant/mol_weights(:)))
+            Cv_avg = sum(Yi_avg(:)*((h_iR(:) - h_iL(:))/(T_R - T_L) - gas_constant/molecular_weights(:)))
         end if
 
         gamma_avg = Cp_avg/Cv_avg
 
-        Phi_avg(:) = (gamma_avg - 1._wp)*(vel_avg_rms/2.0_wp - h_avg_2(:)) + gamma_avg*gas_constant/mol_weights(:)*T_avg
+        Phi_avg(:) = (gamma_avg - 1._wp)*(vel_avg_rms/2.0_wp - h_avg_2(:)) + gamma_avg*gas_constant/molecular_weights(:)*T_avg
         c_sum_Yi_Phi = sum(Yi_avg(:)*Phi_avg(:))
     end if
 

src/simulation/m_riemann_solvers.fpp

@@ -513,8 +513,8 @@ contains
                                 call get_mixture_molecular_weight(Ys_L, MW_L)
                                 call get_mixture_molecular_weight(Ys_R, MW_R)
 
-                                Xs_L(:) = Ys_L(:)*MW_L/mol_weights(:)
-                                Xs_R(:) = Ys_R(:)*MW_R/mol_weights(:)
+                                Xs_L(:) = Ys_L(:)*MW_L/molecular_weights(:)
+                                Xs_R(:) = Ys_R(:)*MW_R/molecular_weights(:)
 
                                 R_gas_L = gas_constant/MW_L
                                 R_gas_R = gas_constant/MW_R
@@ -1029,8 +1029,8 @@ contains
                                     Y_L = qL_prim_rs${XYZ}$_vf(j, k, l, i)
                                     Y_R = qR_prim_rs${XYZ}$_vf(j + 1, k, l, i)
 
-                                    flux_rs${XYZ}$_vf(j, k, l, i) = (s_M*Y_R*rho_R*vel_R(dir_idx(norm_dir)) &
-                                                                     - s_P*Y_L*rho_L*vel_L(dir_idx(norm_dir)) &
+                                    flux_rs${XYZ}$_vf(j, k, l, i) = (s_M*Y_R*rho_R*vel_R(dir_idx(1)) &
+                                                                     - s_P*Y_L*rho_L*vel_L(dir_idx(1)) &
                                                                      + s_M*s_P*(Y_L*rho_L - Y_R*rho_R)) &
                                                                     /(s_M - s_P)
                                     flux_src_rs${XYZ}$_vf(j, k, l, i) = 0._wp
@@ -2573,8 +2573,8 @@ contains
                                     call get_mixture_molecular_weight(Ys_L, MW_L)
                                     call get_mixture_molecular_weight(Ys_R, MW_R)
 
-                                    Xs_L(:) = Ys_L(:)*MW_L/mol_weights(:)
-                                    Xs_R(:) = Ys_R(:)*MW_R/mol_weights(:)
+                                    Xs_L(:) = Ys_L(:)*MW_L/molecular_weights(:)
+                                    Xs_R(:) = Ys_R(:)*MW_R/molecular_weights(:)
 
                                     R_gas_L = gas_constant/MW_L
                                     R_gas_R = gas_constant/MW_R

toolchain/mfc/run/input.py

@@ -76,10 +76,10 @@ def generate_fpp(self, target) -> None:
         # Write the generated Fortran code to the m_thermochem.f90 file with the chosen precision
         common.file_write(
             os.path.join(modules_dir, "m_thermochem.f90"),
-            pyro.codegen.fortran90.gen_thermochem_code(
+            pyro.FortranCodeGenerator().generate(
+                "m_thermochem",
                 self.get_cantera_solution(),
-                module_name="m_thermochem",
-                real_type=real_type
+                pyro.CodeGenerationOptions(scalar_type = real_type)
             ),
             True
         )

toolchain/pyproject.toml

@@ -38,8 +38,8 @@ dependencies = [
     "matplotlib",
 
     # Chemistry
-    "cantera==3.0.1",
-    "pyrometheus==1.0.2"
+    "cantera==3.1.0",
+    "pyrometheus == 1.0.3"
 ]
 
 [tool.hatch.metadata]