diff --git a/docs/documentation/case.md b/docs/documentation/case.md index 512de0eeb..9d3480747 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 000000000..ecacbd965 --- /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 b43905dc7..833976681 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 28a83ca93..4bcab9166 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 44000d752..39e5c7800 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 4375086ed..ea679ba6f 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 951488d5d..6f5cde8bc 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 011c3b2e7..a30ea7246 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 e550d3317..05f349e50 100644 --- a/toolchain/mfc/run/case_dicts.py +++ b/toolchain/mfc/run/case_dicts.py @@ -356,7 +356,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 a297ed181..445d7bbcb 100644 --- a/toolchain/mfc/test/cases.py +++ b/toolchain/mfc/test/cases.py @@ -1030,7 +1030,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:])}"