Skip to content

Commit 6a155cc

Browse files
committed
hi
1 parent d1b8028 commit 6a155cc

File tree

10 files changed

+12336
-13
lines changed

10 files changed

+12336
-13
lines changed

examples/1D_reactive_shocktube/case.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -30,7 +30,7 @@
3030
u_r = -487.34
3131

3232
L = 0.12
33-
Nx = 400 * args.scale
33+
Nx = 2 * 2 * 400 * args.scale
3434
dx = L/Nx
3535
dt = dx/abs(u_r)*0.02
3636
Tend=230e-6
Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,21 @@
1+
import mfc.viz
2+
from tqdm import tqdm
3+
4+
case = mfc.viz.Case(".")
5+
lines = []
6+
7+
for var in tqdm(range(1, 15+1), desc="Loading and Generating Arrays"):
8+
case.load_variable(f"{var}", f"prim.{var}")
9+
elems = [str(x) for x in case.get_data()[59760][str(var)]]
10+
lines.append(f"real(kind(0d0)) :: var{var}(0:400) = [ &")
11+
for i, element in enumerate(elems):
12+
if i == len(elems) - 1:
13+
lines.append(f"{element} &")
14+
else:
15+
lines.append(f"{element}, &")
16+
lines.append("]")
17+
18+
for var in tqdm(range(1, 15+1), desc="Loading Conserved Variables"):
19+
lines.append(f"q_prim_vf({var})%sf(i, j, 0) = var{var}(i)")
20+
21+
print("\n".join(lines))
Lines changed: 128 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,128 @@
1+
#!/usr/bin/env python3
2+
3+
# References:
4+
# + https://doi.org/10.1016/j.ijhydene.2023.03.190: Verification of numerical method
5+
# + https://doi.org/10.1016/j.compfluid.2013.10.014: 4.7. Multi-species reactive shock tube
6+
7+
import json, argparse
8+
import cantera as ct
9+
10+
parser = argparse.ArgumentParser(
11+
prog="1D_reactive_shocktube",
12+
formatter_class=argparse.ArgumentDefaultsHelpFormatter)
13+
14+
parser.add_argument("--mfc", type=json.loads, default='{}', metavar="DICT",
15+
help="MFC's toolchain's internal state.")
16+
parser.add_argument("--no-chem", dest='chemistry', default=True, action="store_false",
17+
help="Disable chemistry.")
18+
parser.add_argument("--scale", type=float, default=1, help="Scale.")
19+
20+
args = parser.parse_args()
21+
22+
ctfile = 'h2o2.yaml'
23+
sol_L = ct.Solution(ctfile)
24+
sol_L.DPX = 0.072, 7173, 'H2:2,O2:1,AR:7'
25+
26+
sol_R = ct.Solution(ctfile)
27+
sol_R.DPX = 0.18075, 35594, 'H2:2,O2:1,AR:7'
28+
29+
u_l = 0
30+
u_r = -487.34
31+
32+
L = 0.12
33+
NScale = 2
34+
Nx = 1 * NScale * 400 * args.scale
35+
Ny = 1 * NScale * 200 * args.scale
36+
dx = L/Nx
37+
dt = dx/abs(u_r)*0.05*0.1*0.2
38+
Tend=230e-6/2
39+
40+
NT=int(Tend/dt)
41+
SAVE_COUNT=100
42+
NS=NT//SAVE_COUNT
43+
44+
case = {
45+
# Logistics ================================================================
46+
'run_time_info' : 'T',
47+
# ==========================================================================
48+
49+
# Computational Domain Parameters ==========================================
50+
'x_domain%beg' : 0,
51+
'x_domain%end' : L,
52+
'y_domain%beg' : 0,
53+
'y_domain%end' : L/2,
54+
'm' : Nx,
55+
'n' : Ny,
56+
'p' : 0,
57+
'dt' : float(dt),
58+
't_step_start' : 0,
59+
't_step_stop' : NT,
60+
't_step_save' : NS,
61+
't_step_print' : 10,
62+
'parallel_io' : 'F' if args.mfc.get("mpi", True) else 'F',
63+
64+
# Simulation Algorithm Parameters ==========================================
65+
'model_eqns' : 2,
66+
'num_fluids' : 1,
67+
'num_patches' : 1,
68+
'mpp_lim' : 'F',
69+
'mixture_err' : 'F',
70+
'time_stepper' : 3,
71+
'weno_order' : 5,
72+
'weno_eps' : 1E-16,
73+
'weno_avg' : 'F',
74+
'mapped_weno' : 'T',
75+
'mp_weno' : 'T',
76+
'riemann_solver' : 2,
77+
'wave_speeds' : 1,
78+
'avg_state' : 2,
79+
'bc_x%beg' :-3,
80+
'bc_x%end' :-3,
81+
'bc_y%beg' :-1,
82+
'bc_y%end' :-1,
83+
84+
# Chemistry ================================================================
85+
'chemistry' : 'F' if not args.chemistry else 'T',
86+
'chem_params%diffusion' : 'F',
87+
'chem_params%reactions' : 'T',
88+
'cantera_file' : ctfile,
89+
# ==========================================================================
90+
91+
# Formatted Database Files Structure Parameters ============================
92+
'format' : 1,
93+
'precision' : 2,
94+
#'prim_vars_wrt' : 'T',
95+
'chem_wrt_T' : 'F',
96+
# ==========================================================================
97+
'rho_wrt' : 'T',
98+
99+
# ==========================================================================
100+
'patch_icpp(1)%geometry' : 7,
101+
'patch_icpp(1)%x_centroid' : L/2,
102+
'patch_icpp(1)%y_centroid' : L/4,
103+
'patch_icpp(1)%length_x' : L,
104+
'patch_icpp(1)%length_y' : L/2,
105+
'patch_icpp(1)%hcid' : 666,
106+
'patch_icpp(1)%vel(1)' : 0,
107+
'patch_icpp(1)%vel(2)' : 0,
108+
'patch_icpp(1)%pres' : sol_L.P,
109+
'patch_icpp(1)%alpha(1)' : 1,
110+
'patch_icpp(1)%alpha_rho(1)' : sol_L.density,
111+
# ==========================================================================
112+
113+
# Fluids Physical Parameters ===============================================
114+
'fluid_pp(1)%gamma' : 1.0E+00/(4.4E+00-1.0E+00),
115+
'fluid_pp(1)%pi_inf' : 0,
116+
# ==========================================================================
117+
}
118+
119+
if args.chemistry:
120+
for i in range(len(sol_L.Y)):
121+
#if sol_L.species_name(i) in ['H2', 'O2', 'OH', 'H2O']:
122+
# case[f'chem_wrt_Y({i + 1})'] = 'T'
123+
#else:
124+
# case[f'chem_wrt_Y({i + 1})'] = 'F'
125+
case[f'patch_icpp(1)%Y({i+1})'] = sol_L.Y[i]
126+
127+
if __name__ == '__main__':
128+
print(json.dumps(case))

examples/nD_perfect_reactor/case.py

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -26,10 +26,10 @@
2626

2727
sol.TPX = 1_600, ct.one_atm, 'H2:0.04, O2:0.02, AR:0.94'
2828

29-
Nx = 25 * args.scale
29+
Nx = int(25 * args.scale)
3030
Tend = 1e-4
3131
s = 1e-2
32-
dt = 1e-7
32+
dt = 1e-8
3333

3434
NT = int(Tend / dt)
3535
SAVE_COUNT = 20
@@ -52,9 +52,9 @@
5252
'p' : Nx,
5353
'dt' : float(dt),
5454
't_step_start' : 0,
55-
't_step_stop' : NT,
56-
't_step_save' : NS,
57-
't_step_print' : NS,
55+
't_step_stop' : 100,
56+
't_step_save' : 100,
57+
't_step_print' : 1,
5858
'parallel_io' : 'T' if args.ndim > 1 and args.mfc.get("mpi", True) else 'F',
5959

6060
# Simulation Algorithm Parameters ==========================================

src/common/m_variables_conversion.fpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -23,7 +23,7 @@ module m_variables_conversion
2323
use m_helper
2424

2525
use m_thermochem, only: &
26-
num_species, get_temperature, get_pressure, &
26+
num_species, get_temperature, get_pressure, gas_constant, &
2727
get_mixture_molecular_weight, get_mixture_energy_mass
2828

2929
! ==========================================================================
@@ -1114,7 +1114,7 @@ contains
11141114
end do
11151115

11161116
call get_mixture_molecular_weight(Ys, mix_mol_weight)
1117-
T = q_prim_vf(T_idx)%sf(j, k, l)
1117+
T = q_prim_vf(E_idx)%sf(j, k, l)*mix_mol_weight/(gas_constant*rho)
11181118
call get_mixture_energy_mass(T, Ys, e_mix)
11191119

11201120
q_cons_vf(E_idx)%sf(j, k, l) = &

src/post_process/m_start_up.f90

Lines changed: 10 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -191,12 +191,15 @@ subroutine s_save_data(t_step, varname, pres, c, H)
191191

192192
integer :: i, j, k, l
193193

194+
print*, proc_rank, "open"
194195
! Opening a new formatted database file
195196
call s_open_formatted_database_file(t_step)
196197

197198
! Adding the grid to the formatted database file
199+
print*, proc_rank, "wrt 1"
198200
call s_write_grid_to_formatted_database_file(t_step)
199201

202+
print*, proc_rank, "wrt 2"
200203
! Computing centered finite-difference coefficients in x-direction
201204
if (omega_wrt(2) .or. omega_wrt(3) .or. qm_wrt .or. schlieren_wrt) then
202205
call s_compute_finite_difference_coefficients(m, x_cc, &
@@ -296,6 +299,7 @@ subroutine s_save_data(t_step, varname, pres, c, H)
296299

297300
! Adding the species' concentrations to the formatted database file ----
298301
if (chemistry) then
302+
print*, proc_rank, "hi 4"
299303
do i = 1, num_species
300304
if (chem_wrt_Y(i) .or. prim_vars_wrt) then
301305
q_sf = q_prim_vf(chemxb + i - 1)%sf(-offset_x%beg:m + offset_x%end, &
@@ -309,7 +313,7 @@ subroutine s_save_data(t_step, varname, pres, c, H)
309313

310314
end if
311315
end do
312-
316+
print*, proc_rank, "hi 3"
313317
if (chem_wrt_T) then
314318
q_sf = q_prim_vf(T_idx)%sf(-offset_x%beg:m + offset_x%end, &
315319
-offset_y%beg:n + offset_y%end, &
@@ -382,7 +386,7 @@ subroutine s_save_data(t_step, varname, pres, c, H)
382386

383387
end if
384388
! ----------------------------------------------------------------------
385-
389+
print*, proc_rank, "hi 2"
386390
! Adding the volume fraction(s) to the formatted database file ---------
387391
if (((model_eqns == 2) .and. (bubbles .neqv. .true.)) &
388392
.or. (model_eqns == 3) &
@@ -483,6 +487,8 @@ subroutine s_save_data(t_step, varname, pres, c, H)
483487
end if
484488
! ----------------------------------------------------------------------
485489

490+
print*, proc_rank, "hi"
491+
486492
! Adding the sound speed to the formatted database file ----------------
487493
if (c_wrt) then
488494
do k = -offset_z%beg, p + offset_z%end
@@ -671,6 +677,8 @@ subroutine s_save_data(t_step, varname, pres, c, H)
671677
end if
672678
end if
673679

680+
681+
print*, proc_rank, "clos 1"
674682
! Closing the formatted database file
675683
call s_close_formatted_database_file()
676684
end subroutine s_save_data

src/post_process/p_main.fpp

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -44,9 +44,16 @@ program p_main
4444

4545
! Time-Marching Loop =======================================================
4646
do
47+
! If not all time-steps are ready to be post-processed, if one process
48+
! is faster than another, a slower processor processing the last available
49+
! step might be killed if a faster one failing to load the next. Therefore,
50+
! we force synchronization here.
51+
call s_mpi_barrier()
52+
4753
call s_perform_time_step(t_step)
4854

4955
call s_save_data(t_step, varname, pres, c, H)
56+
print*, proc_rank, "c"
5057

5158
if (cfl_dt) then
5259
if (t_step == n_save - 1) then
@@ -61,6 +68,7 @@ program p_main
6168
if ((t_step_stop - t_step) < t_step_save .and. t_step_stop /= t_step) then
6269
t_step = t_step_stop - t_step_save
6370
elseif (t_step == t_step_stop) then
71+
print*, proc_rank, "is don"
6472
exit
6573
end if
6674
end if

0 commit comments

Comments
 (0)