Skip to content

Commit 24a324c

Browse files
committed
Lodi for Multicomponent
1 parent 1182c2c commit 24a324c

File tree

8 files changed

+51
-261
lines changed

8 files changed

+51
-261
lines changed

examples/1D_Boundary_Lido/case.py

Lines changed: 0 additions & 113 deletions
This file was deleted.
-48.5 KB
Binary file not shown.
-83.1 KB
Binary file not shown.

examples/1D_reactive_shocktube/viz.py

Lines changed: 0 additions & 80 deletions
This file was deleted.

src/common/m_variables_conversion.fpp

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1529,14 +1529,14 @@ contains
15291529
!Computing the energy from the internal energy of the mixture
15301530
call get_mixture_molecular_weight(Y_k, mix_mol_weight)
15311531
R_gas = gas_constant/mix_mol_weight
1532-
T_K= pres_K/rho_K/R_gas
1532+
T_K = pres_K/rho_K/R_gas
15331533
call get_mixture_energy_mass(T_K, Y_K, E_K)
15341534
E_K = rho_K*E_K + 5e-1_wp*rho_K*vel_K_sum
15351535
else
15361536
! Computing the energy from the pressure
15371537
E_K = gamma_K*pres_K + pi_inf_K &
1538-
+ 5e-1_wp*rho_K*vel_K_sum + qv_K
1539-
end if
1538+
+ 5e-1_wp*rho_K*vel_K_sum + qv_K
1539+
end if
15401540

15411541
! mass flux, this should be \alpha_i \rho_i u_i
15421542
!$acc loop seq
@@ -1556,9 +1556,9 @@ contains
15561556
FK_vf(j, k, l, E_idx) = vel_K(dir_idx(1))*(E_K + pres_K)
15571557

15581558
! Species advection Flux, \rho*u*Y
1559-
if (chemistry) then
1559+
if (chemistry) then
15601560
!$acc loop seq
1561-
do i=1, num_species
1561+
do i = 1, num_species
15621562
FK_vf(j, k, l, i - 1 + chemxb) = vel_K(dir_idx(1))*(rho_K*Y_K(i))
15631563
end do
15641564
end if

src/pre_process/include/1dHardcodedIC.fpp

Lines changed: 0 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -8,30 +8,6 @@
88
select case (patch_icpp(patch_id)%hcid)
99
case (100)
1010
! Put your variable assignments here
11-
!
12-
q_prim_vf(momxb)%sf(i,0,0)=8*exp(-200*((x_cc(i)-0.0054/2)/0.0054)**2)
13-
q_prim_vf(E_idx)%sf(i,0,0)=1.01325*10**5+8.0d0*0.3d0*800.74d0*exp(-200*((x_cc(i)-0.0054/2)/0.0054)**2)
14-
! q_prim_vf(1)%sf(i,0,0)=(1.01325*10**5+2.2d0*0.227d0*800.74d0*exp(-400*((x_cc(i)-0.002/2)/0.002)**2))/(8.314*1000*((0.4+x_cc(i)/0.002*0.2)/2.016+(0.6-x_cc(i)/0.002*0.2)/32) *(300+30*exp(-400*((x_cc(i)-0.002/2)/0.002)**2)))
15-
! q_prim_vf(15)%sf(i,0,0)=300+50*exp(-300*((x_cc(i)-0.002/2)/0.002)**2)
16-
! q_prim_vf(1)%sf(i,0,0)=(1.01325*10**5)/(8.314*1000*((0.4+x_cc(i)/0.002*0.2)/2.016+(0.6-x_cc(i)/0.002*0.2)/32) *(300+50*exp(-400*((x_cc(i)-0.002/2)/0.002)**2)))
17-
q_prim_vf(1)%sf(i,0,0)=0.227+0.227/800.74*(8*exp(-200*((x_cc(i)-0.0054/2)/0.0054)**2))
18-
! q_prim_vf(5)%sf(i,0,0)=0.4+x_cc(i)/0.004*0.2
19-
q_prim_vf(7)%sf(i,0,0)=1.0
20-
! q_prim_vf(8)%sf(i,0,0)=0.6-x_cc(i)/0.004*0.2
21-
! q_prim_vf(10)%sf(i,0,0)=0.1
22-
!q_prim_vf(18)%sf(i,0,0)=0.1
23-
!q_prim_vf(19)%sf(i,0,0)=0.1
24-
!q_prim_vf(20)%sf(i,0,0)=0.1
25-
!q_prim_vf(52)%sf(i,0,0)=x_cc(i)/0.0054*0.4
26-
!q_prim_vf(5)%sf(i,0,0)=0.4
27-
! q_prim_vf(9)%sf(i,0,0)=0.6-x_cc(i)/0.04*0.2
28-
! q_prim_vf(9)%sf(i,0,0)=1.0d0
29-
! q_prim_vf(15)%sf(i,0,0)=300.0d0
30-
!q_prim_vf(1)%sf(1,0,0)=(1.01325*10**5+8.0d0*0.227d0*800.74d0*exp(-400*((x_cc(i)-0.002/2)/0.002)**2))/(8314/17.008*300)
31-
32-
! ! (0.4+x_cc(i)/0.002*0.2)/2.016+(0.6-x_cc(i)/0.002*0.2)/32)!
33-
! Put your variable assignments here
34-
! (0.4+x_cc(i)/0.002*0.2)/2.016+(0.6-x_cc(i)/0.002*0.2)/32)
3511
case default
3612
call s_int_to_str(patch_id, iStr)
3713
call s_mpi_abort("Invalid hcid specified for patch "//trim(iStr))

src/simulation/m_cbc.fpp

Lines changed: 25 additions & 25 deletions
Original file line numberDiff line numberDiff line change
@@ -31,11 +31,11 @@ module m_cbc
3131
use m_compute_cbc
3232

3333
use m_thermochem, only: &
34-
get_mixture_energy_mass,get_mixture_specific_heat_cv_mass, &
35-
get_mixture_specific_heat_cp_mass,gas_constant, &
36-
get_mixture_molecular_weight,get_species_enthalpies_rt, &
37-
molecular_weights, get_species_specific_heats_r, &
38-
get_mole_fractions,get_species_specific_heats_r
34+
get_mixture_energy_mass, get_mixture_specific_heat_cv_mass, &
35+
get_mixture_specific_heat_cp_mass, gas_constant, &
36+
get_mixture_molecular_weight, get_species_enthalpies_rt, &
37+
molecular_weights, get_species_specific_heats_r, &
38+
get_mole_fractions, get_species_specific_heats_r
3939

4040
implicit none
4141

@@ -134,7 +134,7 @@ contains
134134
135135
if (chemistry) then
136136
flux_cbc_index = sys_size
137-
else
137+
else
138138
flux_cbc_index = adv_idx%end
139139
end if
140140
!$acc update device(flux_cbc_index)
@@ -670,7 +670,6 @@ contains
670670
real(wp) :: Cv, Cp, e_mix, Mw, R_gas
671671
real(wp), dimension(num_species) :: Ys, h_k, dYs_dt, dYs_ds, Xs, Gamma_i, Cp_i
672672
673-
674673
real(wp) :: vel_K_sum, vel_dv_dt_sum
675674
676675
integer :: i, j, k, r, q !< Generic loop iterators
@@ -818,7 +817,7 @@ contains
818817
if (chemistry) then
819818
!$acc loop seq
820819
do i = chemxb, chemxe
821-
Ys(i - chemxb + 1) = q_prim_rs${XYZ}$_vf(0, k, r, i)
820+
Ys(i - chemxb + 1) = q_prim_rs${XYZ}$_vf(0, k, r, i)
822821
end do
823822
824823
call get_mixture_molecular_weight(Ys, Mw)
@@ -836,7 +835,6 @@ contains
836835
else if (chem_params%gamma_method == 2) then
837836
!> gamma_method = 2: c_p / c_v where c_p, c_v are specific heats.
838837
call get_mixture_specific_heat_cv_mass(T, Ys, Cv)
839-
call get_mixture_specific_heat_cp_mass(T, Ys, Cp)
840838
gamma = 1.0_wp/(Cp/Cv - 1.0_wp)
841839
end if
842840
else
@@ -866,10 +864,12 @@ contains
866864
dadv_ds(i) = 0._wp
867865
end do
868866
869-
!$acc loop seq
870-
do i=1, num_species
871-
dYs_ds(i) = 0._wp
872-
end do
867+
if (chemistry) then
868+
!$acc loop seq
869+
do i = 1, num_species
870+
dYs_ds(i) = 0._wp
871+
end do
872+
end if
873873
874874
!$acc loop seq
875875
do j = 0, buff_size
@@ -899,12 +899,12 @@ contains
899899
900900
if (chemistry) then
901901
!$acc loop seq
902-
do i=1, num_species
902+
do i = 1, num_species
903903
dYs_ds(i) = q_prim_rs${XYZ}$_vf(j, k, r, chemxb - 1 + i)* &
904-
fd_coef_${XYZ}$ (j, cbc_loc) + &
905-
dYs_ds(i)
904+
fd_coef_${XYZ}$ (j, cbc_loc) + &
905+
dYs_ds(i)
906906
end do
907-
end if
907+
end if
908908
end do
909909
910910
! First-Order Temporal Derivatives of Primitive Variables
@@ -919,7 +919,7 @@ contains
919919
call s_compute_slip_wall_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds)
920920
else if ((cbc_loc == -1 .and. bc${XYZ}$b == BC_CHAR_NR_SUB_BUFFER) .or. &
921921
(cbc_loc == 1 .and. bc${XYZ}$e == BC_CHAR_NR_SUB_BUFFER)) then
922-
call s_compute_nonreflecting_subsonic_buffer_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds)
922+
call s_compute_nonreflecting_subsonic_buffer_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds, dYs_ds)
923923
else if ((cbc_loc == -1 .and. bc${XYZ}$b == BC_CHAR_NR_SUB_INFLOW) .or. &
924924
(cbc_loc == 1 .and. bc${XYZ}$e == BC_CHAR_NR_SUB_INFLOW)) then
925925
call s_compute_nonreflecting_subsonic_inflow_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds)
@@ -1051,19 +1051,19 @@ contains
10511051
sum_Enthalpies = 0._wp
10521052
!$acc loop seq
10531053
do i = 1, num_species
1054-
h_k(i) = h_k(i)*gas_constant/molecular_weights(i)*T
1055-
sum_Enthalpies = sum_Enthalpies+ (rho*h_k(i) - pres*Mw/molecular_weights(i)*Cp/R_gas)*dYs_dt(i)
1054+
h_k(i) = h_k(i)*gas_constant/molecular_weights(i)*T
1055+
sum_Enthalpies = sum_Enthalpies + (rho*h_k(i) - pres*Mw/molecular_weights(i)*Cp/R_gas)*dYs_dt(i)
10561056
end do
10571057
flux_rs${XYZ}$_vf_l(-1, k, r, E_idx) = flux_rs${XYZ}$_vf_l(0, k, r, E_idx) &
1058-
+ ds(0)*((E/rho + pres/rho)*drho_dt + rho*vel_dv_dt_sum + Cp*T*L(2)/(c*c) + sum_Enthalpies)
1058+
+ ds(0)*((E/rho + pres/rho)*drho_dt + rho*vel_dv_dt_sum + Cp*T*L(2)/(c*c) + sum_Enthalpies)
10591059
!$acc loop seq
1060-
do i=1,num_species
1061-
flux_rs${XYZ}$_vf_l(-1, k, r, i - 1 + chemxb) = flux_rs${XYZ}$_vf_l(0, k, r, chemxb + i - 1) &
1062-
+ ds(0)*(drho_dt*Ys(i) + rho*dYs_dt(i))
1060+
do i = 1, num_species
1061+
flux_rs${XYZ}$_vf_l(-1, k, r, i - 1 + chemxb) = flux_rs${XYZ}$_vf_l(0, k, r, chemxb + i - 1) &
1062+
+ ds(0)*(drho_dt*Ys(i) + rho*dYs_dt(i))
10631063
end do
10641064
else
10651065
flux_rs${XYZ}$_vf_l(-1, k, r, E_idx) = flux_rs${XYZ}$_vf_l(0, k, r, E_idx) &
1066-
+ ds(0)*(pres*dgamma_dt &
1066+
+ ds(0)*(pres*dgamma_dt &
10671067
+ gamma*dpres_dt &
10681068
+ dpi_inf_dt &
10691069
+ dqv_dt &

0 commit comments

Comments
 (0)