diff --git a/docs/documentation/case.md b/docs/documentation/case.md index 98e75b326b..43de6842a1 100644 --- a/docs/documentation/case.md +++ b/docs/documentation/case.md @@ -993,6 +993,22 @@ When ``cyl_coord = 'T'`` is set in 2D the following constraints must be met: - `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 | +| `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 + +- `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 + + ## Enumerations ### Boundary conditions diff --git a/examples/1D_multispecies_diffusion/case.py b/examples/1D_multispecies_diffusion/case.py new file mode 100644 index 0000000000..ecacbd965c --- /dev/null +++ b/examples/1D_multispecies_diffusion/case.py @@ -0,0 +1,78 @@ +#!/usr/bin/env python3 +# References: +# + https://doi.org/10.1016/j.compfluid.2013.10.014: 4.4. Multicomponent diffusion test case + +import json +import argparse +import math +import cantera as ct + +ctfile = "gri30.yaml" +sol_L = ct.Solution(ctfile) +sol_L.TPX = 300, 8000, "O2:2,N2:2,H2O:5" + +L = 0.05 +Nx = 100 +dx = L / Nx +dt = 0.3e-6 +Tend = 0.05 + +NT = int(Tend / dt) +SAVE_COUNT = 2000 +NS = 2000 +case = { + "run_time_info": "T", + "x_domain%beg": 0, + "x_domain%end": +L, + "m": Nx, + "n": 0, + "p": 0, + "dt": float(dt), + "t_step_start": 0, + "t_step_stop": NT, + "t_step_save": NS, + "t_step_print": NS, + "parallel_io": "F", + "model_eqns": 2, + "num_fluids": 1, + "num_patches": 1, + "mpp_lim": "F", + "mixture_err": "F", + "time_stepper": 3, + "weno_order": 5, + "weno_eps": 1e-16, + "weno_avg": "F", + "mapped_weno": "T", + "mp_weno": "T", + "riemann_solver": 2, + "wave_speeds": 2, + "avg_state": 1, + "bc_x%beg": -1, + "bc_x%end": -1, + "viscous": "F", + "chemistry": "T", + "chem_params%diffusion": "T", + "chem_params%reactions": "F", + "chem_params%transport_model": 2, # Unity-Lewis + "format": 1, + "precision": 2, + "prim_vars_wrt": "T", + "chem_wrt_T": "T", + "patch_icpp(1)%geometry": 1, + "patch_icpp(1)%hcid": 182, + "patch_icpp(1)%x_centroid": L / 2, + "patch_icpp(1)%length_x": L, + "patch_icpp(1)%vel(1)": "0", + "patch_icpp(1)%pres": 1.01325e5, + "patch_icpp(1)%alpha(1)": 1, + "patch_icpp(1)%alpha_rho(1)": 1, + "fluid_pp(1)%gamma": 1.0e00 / (1.9326e00 - 1.0e00), + "fluid_pp(1)%pi_inf": 0, + "cantera_file": ctfile, +} + +for i in range(len(sol_L.Y)): + case[f"chem_wrt_Y({i + 1})"] = "T" + case[f"patch_icpp(1)%Y({i+1})"] = 0.0 +if __name__ == "__main__": + print(json.dumps(case)) diff --git a/src/common/m_chemistry.fpp b/src/common/m_chemistry.fpp index b43905dc7f..8339766813 100644 --- a/src/common/m_chemistry.fpp +++ b/src/common/m_chemistry.fpp @@ -13,7 +13,7 @@ module m_chemistry get_mole_fractions, get_species_binary_mass_diffusivities, & get_species_mass_diffusivities_mixavg, gas_constant, get_mixture_molecular_weight, & get_mixture_energy_mass, get_mixture_thermal_conductivity_mixavg, get_species_enthalpies_rt, & - get_mixture_viscosity_mixavg + get_mixture_viscosity_mixavg, get_mixture_specific_heat_cp_mass, get_mixture_enthalpy_mass use m_global_parameters @@ -174,10 +174,13 @@ contains real(wp), dimension(num_species) :: Xs_L, Xs_R, Xs_cell, Ys_L, Ys_R, Ys_cell real(wp), dimension(num_species) :: mass_diffusivities_mixavg1, mass_diffusivities_mixavg2 real(wp), dimension(num_species) :: mass_diffusivities_mixavg_Cell, dXk_dxi, h_l, h_r, h_k - real(wp), dimension(num_species) :: Mass_Diffu_Flux + real(wp), dimension(num_species) :: Mass_Diffu_Flux, dYk_dxi real(wp) :: Mass_Diffu_Energy real(wp) :: MW_L, MW_R, MW_cell, Rgas_L, Rgas_R, T_L, T_R, P_L, P_R, rho_L, rho_R, rho_cell, rho_Vic real(wp) :: lambda_L, lambda_R, lambda_Cell, dT_dxi, grid_spacing + real(wp) :: Cp_L, Cp_R + real(wp) :: diffusivity_L, diffusivity_R, diffusivity_cell + real(wp) :: hmix_L, hmix_R, dh_dxi integer :: x, y, z, i, n, eqn integer, dimension(3) :: offsets @@ -187,119 +190,218 @@ contains $:GPU_UPDATE(device='[isc1,isc2,isc3]') if (chemistry) then + ! Set offsets based on direction using array indexing offsets = 0 offsets(idir) = 1 - #:block UNDEF_AMD - $:GPU_PARALLEL_LOOP(collapse=3, private='[x,y,z,Ys_L, Ys_R, Ys_cell, Xs_L, Xs_R, mass_diffusivities_mixavg1, mass_diffusivities_mixavg2, mass_diffusivities_mixavg_Cell, h_l, h_r, Xs_cell, h_k, dXk_dxi,Mass_Diffu_Flux, Mass_Diffu_Energy, MW_L, MW_R, MW_cell, Rgas_L, Rgas_R, T_L, T_R, P_L, P_R, rho_L, rho_R, rho_cell, rho_Vic, lambda_L, lambda_R, lambda_Cell, dT_dxi, grid_spacing]', copyin='[offsets]') - do z = isc3%beg, isc3%end - do y = isc2%beg, isc2%end - do x = isc1%beg, isc1%end - ! Calculate grid spacing using direction-based indexing - select case (idir) - case (1) - grid_spacing = x_cc(x + 1) - x_cc(x) - case (2) - grid_spacing = y_cc(y + 1) - y_cc(y) - case (3) - grid_spacing = z_cc(z + 1) - z_cc(z) - end select - - ! Extract species mass fractions - $:GPU_LOOP(parallelism='[seq]') - do i = chemxb, chemxe - Ys_L(i - chemxb + 1) = q_prim_qp(i)%sf(x, y, z) - Ys_R(i - chemxb + 1) = q_prim_qp(i)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) - Ys_cell(i - chemxb + 1) = 0.5_wp*(Ys_L(i - chemxb + 1) + Ys_R(i - chemxb + 1)) - end do - - ! Calculate molecular weights and mole fractions - call get_mixture_molecular_weight(Ys_L, MW_L) - call get_mixture_molecular_weight(Ys_R, MW_R) - MW_cell = 0.5_wp*(MW_L + MW_R) - - call get_mole_fractions(MW_L, Ys_L, Xs_L) - call get_mole_fractions(MW_R, Ys_R, Xs_R) - - ! Calculate gas constants and thermodynamic properties - Rgas_L = gas_constant/MW_L - Rgas_R = gas_constant/MW_R - - P_L = q_prim_qp(E_idx)%sf(x, y, z) - P_R = q_prim_qp(E_idx)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) - - rho_L = q_prim_qp(1)%sf(x, y, z) - rho_R = q_prim_qp(1)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) - - T_L = P_L/rho_L/Rgas_L - T_R = P_R/rho_R/Rgas_R - - rho_cell = 0.5_wp*(rho_L + rho_R) - dT_dxi = (T_R - T_L)/grid_spacing - - ! Get transport properties - call get_species_mass_diffusivities_mixavg(P_L, T_L, Ys_L, mass_diffusivities_mixavg1) - call get_species_mass_diffusivities_mixavg(P_R, T_R, Ys_R, mass_diffusivities_mixavg2) - - call get_mixture_thermal_conductivity_mixavg(T_L, Ys_L, lambda_L) - call get_mixture_thermal_conductivity_mixavg(T_R, Ys_R, lambda_R) - - call get_species_enthalpies_rt(T_L, h_l) - call get_species_enthalpies_rt(T_R, h_r) - - ! Calculate species properties and gradients - $:GPU_LOOP(parallelism='[seq]') - do i = chemxb, chemxe - h_l(i - chemxb + 1) = h_l(i - chemxb + 1)*gas_constant*T_L/molecular_weights(i - chemxb + 1) - h_r(i - chemxb + 1) = h_r(i - chemxb + 1)*gas_constant*T_R/molecular_weights(i - chemxb + 1) - Xs_cell(i - chemxb + 1) = 0.5_wp*(Xs_L(i - chemxb + 1) + Xs_R(i - chemxb + 1)) - h_k(i - chemxb + 1) = 0.5_wp*(h_l(i - chemxb + 1) + h_r(i - chemxb + 1)) - dXk_dxi(i - chemxb + 1) = (Xs_R(i - chemxb + 1) - Xs_L(i - chemxb + 1))/grid_spacing + ! Model 1: Mixture-Average Transport + if (chem_params%transport_model == 1) then + #:block UNDEF_AMD + ! Note: Added 'i' and 'eqn' to private list. + $:GPU_PARALLEL_LOOP(collapse=3, private='[x,y,z,i,eqn,Ys_L, Ys_R, Ys_cell, Xs_L, Xs_R, mass_diffusivities_mixavg1, mass_diffusivities_mixavg2, mass_diffusivities_mixavg_Cell, h_l, h_r, Xs_cell, h_k, dXk_dxi,Mass_Diffu_Flux, Mass_Diffu_Energy, MW_L, MW_R, MW_cell, Rgas_L, Rgas_R, T_L, T_R, P_L, P_R, rho_L, rho_R, rho_cell, rho_Vic, lambda_L, lambda_R, lambda_Cell, dT_dxi, grid_spacing]', copyin='[offsets]') + do z = isc3%beg, isc3%end + do y = isc2%beg, isc2%end + do x = isc1%beg, isc1%end + ! Calculate grid spacing using direction-based indexing + select case (idir) + case (1) + grid_spacing = x_cc(x + 1) - x_cc(x) + case (2) + grid_spacing = y_cc(y + 1) - y_cc(y) + case (3) + grid_spacing = z_cc(z + 1) - z_cc(z) + end select + + ! Extract species mass fractions + $:GPU_LOOP(parallelism='[seq]') + do i = chemxb, chemxe + Ys_L(i - chemxb + 1) = q_prim_qp(i)%sf(x, y, z) + Ys_R(i - chemxb + 1) = q_prim_qp(i)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) + Ys_cell(i - chemxb + 1) = 0.5_wp*(Ys_L(i - chemxb + 1) + Ys_R(i - chemxb + 1)) + end do + + ! Calculate molecular weights and mole fractions + call get_mixture_molecular_weight(Ys_L, MW_L) + call get_mixture_molecular_weight(Ys_R, MW_R) + MW_cell = 0.5_wp*(MW_L + MW_R) + + call get_mole_fractions(MW_L, Ys_L, Xs_L) + call get_mole_fractions(MW_R, Ys_R, Xs_R) + + ! Calculate gas constants and thermodynamic properties + Rgas_L = gas_constant/MW_L + Rgas_R = gas_constant/MW_R + + P_L = q_prim_qp(E_idx)%sf(x, y, z) + P_R = q_prim_qp(E_idx)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) + + rho_L = q_prim_qp(1)%sf(x, y, z) + rho_R = q_prim_qp(1)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) + + T_L = P_L/rho_L/Rgas_L + T_R = P_R/rho_R/Rgas_R + + rho_cell = 0.5_wp*(rho_L + rho_R) + dT_dxi = (T_R - T_L)/grid_spacing + + ! Get transport properties + call get_species_mass_diffusivities_mixavg(P_L, T_L, Ys_L, mass_diffusivities_mixavg1) + call get_species_mass_diffusivities_mixavg(P_R, T_R, Ys_R, mass_diffusivities_mixavg2) + + call get_mixture_thermal_conductivity_mixavg(T_L, Ys_L, lambda_L) + call get_mixture_thermal_conductivity_mixavg(T_R, Ys_R, lambda_R) + + call get_species_enthalpies_rt(T_L, h_l) + call get_species_enthalpies_rt(T_R, h_r) + + ! Calculate species properties and gradients + $:GPU_LOOP(parallelism='[seq]') + do i = chemxb, chemxe + h_l(i - chemxb + 1) = h_l(i - chemxb + 1)*gas_constant*T_L/molecular_weights(i - chemxb + 1) + h_r(i - chemxb + 1) = h_r(i - chemxb + 1)*gas_constant*T_R/molecular_weights(i - chemxb + 1) + Xs_cell(i - chemxb + 1) = 0.5_wp*(Xs_L(i - chemxb + 1) + Xs_R(i - chemxb + 1)) + h_k(i - chemxb + 1) = 0.5_wp*(h_l(i - chemxb + 1) + h_r(i - chemxb + 1)) + dXk_dxi(i - chemxb + 1) = (Xs_R(i - chemxb + 1) - Xs_L(i - chemxb + 1))/grid_spacing + end do + + ! Calculate mixture-averaged diffusivities + $:GPU_LOOP(parallelism='[seq]') + do i = chemxb, chemxe + mass_diffusivities_mixavg_Cell(i - chemxb + 1) = & + (mass_diffusivities_mixavg2(i - chemxb + 1) + mass_diffusivities_mixavg1(i - chemxb + 1))/2.0_wp + end do + + lambda_Cell = 0.5_wp*(lambda_R + lambda_L) + + ! Calculate mass diffusion fluxes + rho_Vic = 0.0_wp + Mass_Diffu_Energy = 0.0_wp + + $:GPU_LOOP(parallelism='[seq]') + do eqn = chemxb, chemxe + Mass_Diffu_Flux(eqn - chemxb + 1) = rho_cell*mass_diffusivities_mixavg_Cell(eqn - chemxb + 1)* & + molecular_weights(eqn - chemxb + 1)/MW_cell*dXk_dxi(eqn - chemxb + 1) + rho_Vic = rho_Vic + Mass_Diffu_Flux(eqn - chemxb + 1) + Mass_Diffu_Energy = Mass_Diffu_Energy + h_k(eqn - chemxb + 1)*Mass_Diffu_Flux(eqn - chemxb + 1) + end do + + ! Apply corrections for mass conservation + $:GPU_LOOP(parallelism='[seq]') + do eqn = chemxb, chemxe + Mass_Diffu_Energy = Mass_Diffu_Energy - h_k(eqn - chemxb + 1)*Ys_cell(eqn - chemxb + 1)*rho_Vic + Mass_Diffu_Flux(eqn - chemxb + 1) = Mass_Diffu_Flux(eqn - chemxb + 1) - rho_Vic*Ys_cell(eqn - chemxb + 1) + end do + + ! Add thermal conduction contribution + Mass_Diffu_Energy = lambda_Cell*dT_dxi + Mass_Diffu_Energy + + ! Update flux arrays + flux_src_vf(E_idx)%sf(x, y, z) = flux_src_vf(E_idx)%sf(x, y, z) - Mass_Diffu_Energy + + $:GPU_LOOP(parallelism='[seq]') + do eqn = chemxb, chemxe + flux_src_vf(eqn)%sf(x, y, z) = flux_src_vf(eqn)%sf(x, y, z) - Mass_Diffu_Flux(eqn - chemxb + 1) + end do end do - - ! Calculate mixture-averaged diffusivities - $:GPU_LOOP(parallelism='[seq]') - do i = chemxb, chemxe - mass_diffusivities_mixavg_Cell(i - chemxb + 1) = & - (mass_diffusivities_mixavg2(i - chemxb + 1) + mass_diffusivities_mixavg1(i - chemxb + 1))/2.0_wp - end do - - lambda_Cell = 0.5_wp*(lambda_R + lambda_L) - - ! Calculate mass diffusion fluxes - rho_Vic = 0.0_wp - Mass_Diffu_Energy = 0.0_wp - - $:GPU_LOOP(parallelism='[seq]') - do eqn = chemxb, chemxe - Mass_Diffu_Flux(eqn - chemxb + 1) = rho_cell*mass_diffusivities_mixavg_Cell(eqn - chemxb + 1)* & - molecular_weights(eqn - chemxb + 1)/MW_cell*dXk_dxi(eqn - chemxb + 1) - rho_Vic = rho_Vic + Mass_Diffu_Flux(eqn - chemxb + 1) - Mass_Diffu_Energy = Mass_Diffu_Energy + h_k(eqn - chemxb + 1)*Mass_Diffu_Flux(eqn - chemxb + 1) - end do - - ! Apply corrections for mass conservation - $:GPU_LOOP(parallelism='[seq]') - do eqn = chemxb, chemxe - Mass_Diffu_Energy = Mass_Diffu_Energy - h_k(eqn - chemxb + 1)*Ys_cell(eqn - chemxb + 1)*rho_Vic - Mass_Diffu_Flux(eqn - chemxb + 1) = Mass_Diffu_Flux(eqn - chemxb + 1) - rho_Vic*Ys_cell(eqn - chemxb + 1) - end do - - ! Add thermal conduction contribution - Mass_Diffu_Energy = lambda_Cell*dT_dxi + Mass_Diffu_Energy - - ! Update flux arrays - flux_src_vf(E_idx)%sf(x, y, z) = flux_src_vf(E_idx)%sf(x, y, z) - Mass_Diffu_Energy - - $:GPU_LOOP(parallelism='[seq]') - do eqn = chemxb, chemxe - flux_src_vf(eqn)%sf(x, y, z) = flux_src_vf(eqn)%sf(x, y, z) - Mass_diffu_Flux(eqn - chemxb + 1) + end do + end do + $:END_GPU_PARALLEL_LOOP() + #:endblock UNDEF_AMD + + ! Model 2: Unity Lewis Number + else if (chem_params%transport_model == 2) then + #:block UNDEF_AMD + ! Note: Added ALL scalars and 'i'/'eqn' to private list to prevent race conditions. + $:GPU_PARALLEL_LOOP(collapse=3, private='[x,y,z,i,eqn,Ys_L, Ys_R, Ys_cell, dYk_dxi, Mass_Diffu_Flux, grid_spacing, MW_L, MW_R, MW_cell, Rgas_L, Rgas_R, P_L, P_R, rho_L, rho_R, rho_cell, T_L, T_R, Cp_L, Cp_R, hmix_L, hmix_R, dh_dxi, lambda_L, lambda_R, lambda_Cell, diffusivity_L, diffusivity_R, diffusivity_cell, Mass_Diffu_Energy]', copyin='[offsets]') + do z = isc3%beg, isc3%end + do y = isc2%beg, isc2%end + do x = isc1%beg, isc1%end + ! Calculate grid spacing using direction-based indexing + select case (idir) + case (1) + grid_spacing = x_cc(x + 1) - x_cc(x) + case (2) + grid_spacing = y_cc(y + 1) - y_cc(y) + case (3) + grid_spacing = z_cc(z + 1) - z_cc(z) + end select + + ! Extract species mass fractions + $:GPU_LOOP(parallelism='[seq]') + do i = chemxb, chemxe + Ys_L(i - chemxb + 1) = q_prim_qp(i)%sf(x, y, z) + Ys_R(i - chemxb + 1) = q_prim_qp(i)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) + Ys_cell(i - chemxb + 1) = 0.5_wp*(Ys_L(i - chemxb + 1) + Ys_R(i - chemxb + 1)) + end do + + ! Calculate molecular weights and mole fractions + call get_mixture_molecular_weight(Ys_L, MW_L) + call get_mixture_molecular_weight(Ys_R, MW_R) + MW_cell = 0.5_wp*(MW_L + MW_R) + + ! Calculate gas constants and thermodynamic properties + Rgas_L = gas_constant/MW_L + Rgas_R = gas_constant/MW_R + + P_L = q_prim_qp(E_idx)%sf(x, y, z) + P_R = q_prim_qp(E_idx)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) + + rho_L = q_prim_qp(1)%sf(x, y, z) + rho_R = q_prim_qp(1)%sf(x + offsets(1), y + offsets(2), z + offsets(3)) + + T_L = P_L/rho_L/Rgas_L + T_R = P_R/rho_R/Rgas_R + + rho_cell = 0.5_wp*(rho_L + rho_R) + + call get_mixture_specific_heat_cp_mass(T_L, Ys_L, Cp_L) + call get_mixture_specific_heat_cp_mass(T_R, Ys_R, Cp_R) + call get_mixture_enthalpy_mass(T_L, Ys_L, hmix_L) + call get_mixture_enthalpy_mass(T_R, Ys_R, hmix_R) + dh_dxi = (hmix_R - hmix_L)/grid_spacing + + ! Get transport properties + call get_mixture_thermal_conductivity_mixavg(T_L, Ys_L, lambda_L) + call get_mixture_thermal_conductivity_mixavg(T_R, Ys_R, lambda_R) + + ! Calculate species properties and gradients + $:GPU_LOOP(parallelism='[seq]') + do i = chemxb, chemxe + dYk_dxi(i - chemxb + 1) = (Ys_R(i - chemxb + 1) - & + Ys_L(i - chemxb + 1))/grid_spacing + end do + + ! Calculate mixture-averaged diffusivities + diffusivity_L = lambda_L/rho_L/Cp_L + diffusivity_R = lambda_R/rho_R/Cp_R + + lambda_Cell = 0.5_wp*(lambda_R + lambda_L) + diffusivity_cell = 0.5_wp*(diffusivity_R + diffusivity_L) + + ! Calculate mass diffusion fluxes + Mass_Diffu_Energy = 0.0_wp + + $:GPU_LOOP(parallelism='[seq]') + do eqn = chemxb, chemxe + Mass_Diffu_Flux(eqn - chemxb + 1) = rho_cell* & + diffusivity_cell* & + dYk_dxi(eqn - chemxb + 1) + end do + Mass_Diffu_Energy = rho_cell*diffusivity_cell*dh_dxi + + ! Update flux arrays + flux_src_vf(E_idx)%sf(x, y, z) = flux_src_vf(E_idx)%sf(x, y, z) - Mass_Diffu_Energy + + $:GPU_LOOP(parallelism='[seq]') + do eqn = chemxb, chemxe + flux_src_vf(eqn)%sf(x, y, z) = flux_src_vf(eqn)%sf(x, y, z) - Mass_Diffu_Flux(eqn - chemxb + 1) + end do end do end do end do - end do - $:END_GPU_PARALLEL_LOOP() - #:endblock UNDEF_AMD + $:END_GPU_PARALLEL_LOOP() + #:endblock UNDEF_AMD + end if end if end subroutine s_compute_chemistry_diffusion_flux diff --git a/src/common/m_derived_types.fpp b/src/common/m_derived_types.fpp index b26ecbcfa7..f3ef91e831 100644 --- a/src/common/m_derived_types.fpp +++ b/src/common/m_derived_types.fpp @@ -456,6 +456,7 @@ module m_derived_types !> gamma_method = 1: Ref. Section 2.3.1 Formulation of doi:10.7907/ZKW8-ES97. !> gamma_method = 2: c_p / c_v where c_p, c_v are specific heats. integer :: gamma_method + integer :: transport_model end type chemistry_parameters !> Lagrangian bubble parameters diff --git a/src/post_process/m_global_parameters.fpp b/src/post_process/m_global_parameters.fpp index 44000d7524..39e5c7800c 100644 --- a/src/post_process/m_global_parameters.fpp +++ b/src/post_process/m_global_parameters.fpp @@ -303,6 +303,7 @@ module m_global_parameters real(wp) :: rhoref, pref !> @} + type(chemistry_parameters) :: chem_params !> @name Bubble modeling variables and parameters !> @{ integer :: nb @@ -420,6 +421,9 @@ contains #:endfor #:endfor + chem_params%gamma_method = 1 + chem_params%transport_model = 1 + ! Fluids physical parameters do i = 1, num_fluids_max fluid_pp(i)%gamma = dflt_real diff --git a/src/pre_process/m_global_parameters.fpp b/src/pre_process/m_global_parameters.fpp index 4375086ed1..ea679ba6f5 100644 --- a/src/pre_process/m_global_parameters.fpp +++ b/src/pre_process/m_global_parameters.fpp @@ -220,6 +220,7 @@ module m_global_parameters real(wp) :: rhoref, pref !< Reference parameters for Tait EOS + type(chemistry_parameters) :: chem_params !> @name Bubble modeling !> @{ integer :: nb @@ -580,6 +581,9 @@ contains patch_ib(i)%rotation_matrix_inverse = patch_ib(i)%rotation_matrix end do + chem_params%gamma_method = 1 + chem_params%transport_model = 1 + ! Fluids physical parameters do i = 1, num_fluids_max fluid_pp(i)%gamma = dflt_real diff --git a/src/simulation/m_global_parameters.fpp b/src/simulation/m_global_parameters.fpp index e3b7bb696f..50967956ae 100644 --- a/src/simulation/m_global_parameters.fpp +++ b/src/simulation/m_global_parameters.fpp @@ -649,6 +649,7 @@ contains chem_params%diffusion = .false. chem_params%reactions = .false. chem_params%gamma_method = 1 + chem_params%transport_model = 1 num_bc_patches = 0 bc_io = .false. diff --git a/src/simulation/m_mpi_proxy.fpp b/src/simulation/m_mpi_proxy.fpp index f135cf8736..527e4b1aa4 100644 --- a/src/simulation/m_mpi_proxy.fpp +++ b/src/simulation/m_mpi_proxy.fpp @@ -125,7 +125,7 @@ contains call MPI_BCAST(chem_params%${VAR}$, 1, MPI_LOGICAL, 0, MPI_COMM_WORLD, ierr) #:endfor - #:for VAR in [ 'gamma_method' ] + #:for VAR in [ 'gamma_method', 'transport_model' ] call MPI_BCAST(chem_params%${VAR}$, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr) #:endfor end if diff --git a/toolchain/mfc/run/case_dicts.py b/toolchain/mfc/run/case_dicts.py index 1075477335..9240201527 100644 --- a/toolchain/mfc/run/case_dicts.py +++ b/toolchain/mfc/run/case_dicts.py @@ -355,7 +355,7 @@ def analytic(self): for var in [ 'diffusion', 'reactions' ]: SIMULATION[f'chem_params%{var}'] = ParamType.LOG -for var in [ 'gamma_method' ]: +for var in [ 'gamma_method', 'transport_model']: SIMULATION[f'chem_params%{var}'] = ParamType.INT for var in ["R0ref", "p0ref", "rho0ref", "T0ref", "ss", "pv", "vd", diff --git a/toolchain/mfc/test/cases.py b/toolchain/mfc/test/cases.py index 61008249b5..228c24ac1c 100644 --- a/toolchain/mfc/test/cases.py +++ b/toolchain/mfc/test/cases.py @@ -1031,7 +1031,7 @@ def foreach_example(): "2D_backward_facing_step", "2D_forward_facing_step", "1D_convergence", - "3D_IGR_33jet"] + "3D_IGR_33jet","1D_multispecies_diffusion"] if path in casesToSkip: continue name = f"{path.split('_')[0]} -> Example -> {'_'.join(path.split('_')[1:])}"