Skip to content

Commit d55b19e

Browse files
DimAdam-01DimAdam-01
committed
Improvements
1 parent 7e61509 commit d55b19e

File tree

5 files changed

+301
-89
lines changed

5 files changed

+301
-89
lines changed

‎docs/documentation/case.md‎

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -236,6 +236,15 @@ and use `patch_icpp(i)%%geometry = 7` and `patch_icpp(i)%%hcid = 200` in the inp
236236
Additional variables can be declared in `Hardcoded1[2,3]DVariables` and used in `hardcoded1[2,3]D`.
237237
As a convention, any hard coded patches that are part of the MFC master branch should be identified as 1[2,3]xx where the first digit indicates the number of dimensions.
238238

239+
Also Several custom hard-coded patches are already defined to support common workflows. These include:
240+
241+
case(100) in 1dHardcodedIC.fpp: used to initialize the domain with data from a previous 1D simulation or external source (e.g. Cantera), allowing reuse of existing 1D profiles.
242+
243+
case(270) in 2dHardcodedIC.fpp: used to create 2D fields by extruding a 1D profile along the second spatial dimensions, This allows the solution of a 1D case at a
244+
specific timestep to be extended uniformly in the transverse direction and used as a full 2D initial condition.
245+
246+
case(370) in 3dHardcodedIC.fpp: used to build 3D domains by extruding 2D data along the third dimension.
247+
239248
#### Parameter Descriptions
240249

241250
- `num_patches` defines the total number of patches defined in the domain.

‎examples/2D_Detonation/case.py‎

Lines changed: 170 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,170 @@
1+
import json, argparse
2+
import cantera as ct
3+
4+
parser = argparse.ArgumentParser(
5+
prog="2D_detonation",
6+
formatter_class=argparse.ArgumentDefaultsHelpFormatter)
7+
8+
parser.add_argument("--mfc", type=json.loads, default='{}', metavar="DICT",
9+
help="MFC's toolchain's internal state.")
10+
parser.add_argument("--no-chem", dest='chemistry', default=True, action="store_false",
11+
help="Disable chemistry.")
12+
parser.add_argument("--scale", type=float, default=1, help="Scale.")
13+
14+
args = parser.parse_args()
15+
16+
ctfile = 'h2o2.yaml'
17+
sol_L = ct.Solution(ctfile)
18+
sol_L.DPX = 0.072, 7173, 'H2:2,O2:1,AR:7'
19+
20+
sol_R = ct.Solution(ctfile)
21+
sol_R.DPX = 0.18075, 35594, 'H2:2,O2:1,AR:7'
22+
23+
u_l = 0
24+
u_r = -487.34
25+
26+
L = 0.12
27+
Nx = 800
28+
Ny = 200
29+
dx = L/Nx
30+
dy = (L/4)/Ny
31+
dt = min(dx,dy)/abs(u_r)*0.05*0.1*0.2
32+
Tend=830e-6
33+
34+
NT=int(Tend/dt)
35+
SAVE_COUNT=100
36+
NS=NT//SAVE_COUNT
37+
38+
case = {
39+
# Logistics ================================================================
40+
'run_time_info' : 'T',
41+
# ==========================================================================
42+
43+
# Computational Domain Parameters ==========================================
44+
'x_domain%beg' : 0,
45+
'x_domain%end' : L,
46+
'y_domain%beg' : 0,
47+
'y_domain%end' : L/4,
48+
'm' : Nx,
49+
'n' : Ny,
50+
'p' : 0,
51+
'dt' : float(dt),
52+
't_step_start' : 0,
53+
't_step_stop' : 1,
54+
't_step_save' : 1,
55+
't_step_print' : 1,
56+
'parallel_io' : 'F', #if args.mfc.get("mpi", True) else 'F',
57+
58+
# Simulation Algorithm Parameters ==========================================
59+
'model_eqns' : 2,
60+
'num_fluids' : 1,
61+
'num_patches' : 2,
62+
'mpp_lim' : 'F',
63+
'mixture_err' : 'F',
64+
'time_stepper' : 3,
65+
'weno_order' : 5,
66+
'weno_eps' : 1E-16,
67+
'weno_avg' : 'F',
68+
'mapped_weno' : 'T',
69+
'mp_weno' : 'T',
70+
'riemann_solver' : 2,
71+
'wave_speeds' : 1,
72+
'avg_state' : 2,
73+
'bc_x%beg' :-3,
74+
'bc_x%end' :-3,
75+
'bc_y%beg' :-1,
76+
'bc_y%end' :-1,
77+
78+
# Chemistry ================================================================
79+
'chemistry' : 'T' if not args.chemistry else 'T',
80+
'chem_params%diffusion' : 'F',
81+
'chem_params%reactions' : 'F',
82+
'cantera_file' : ctfile,
83+
# ==========================================================================
84+
85+
# Formatted Database Files Structure Parameters ============================
86+
'format' : 1,
87+
'precision' : 2,
88+
'prim_vars_wrt' : 'T',
89+
'chem_wrt_T' : 'T',
90+
'vel_wrt(1)' : 'T',
91+
'vel_wrt(2)' : 'T',
92+
'pres_wrt' : 'T',
93+
# ==========================================================================
94+
# 'rho_wrt' : 'T',
95+
96+
97+
'patch_icpp(1)%geometry' : 3,
98+
'patch_icpp(1)%x_centroid' : L/2,
99+
'patch_icpp(1)%y_centroid' : L/8,
100+
'patch_icpp(1)%length_x' : L,
101+
'patch_icpp(1)%length_y' : L/4,
102+
'patch_icpp(1)%vel(1)' : -487.34,
103+
'patch_icpp(1)%vel(2)' : 0.0,
104+
'patch_icpp(1)%pres' : sol_R.P,
105+
'patch_icpp(1)%alpha(1)' : 1,
106+
'patch_icpp(1)%alpha_rho(1)' : sol_R.density,
107+
# 'patch_icpp(1)%Y(1)' : 0.5,
108+
# 'patch_icpp(1)%Y(2)' : 0.5,
109+
# ==========================================================================
110+
# ==========================================================================
111+
# ==========================================================================
112+
113+
'patch_icpp(2)%geometry' : 7,
114+
'patch_icpp(2)%x_centroid' : L/4,
115+
'patch_icpp(2)%y_centroid' : L/8,
116+
'patch_icpp(2)%length_x' : L/2,
117+
'patch_icpp(2)%length_y' : L/4,
118+
'patch_icpp(2)%hcid' : 270,
119+
'patch_icpp(2)%vel(1)' : 0,
120+
'patch_icpp(2)%vel(2)' : 0,
121+
'patch_icpp(2)%pres' : sol_R.P,
122+
'patch_icpp(2)%alpha(1)' : 1,
123+
'patch_icpp(2)%alpha_rho(1)' : sol_R.density,
124+
# 'patch_icpp(1)%alter_patch(1)' : 'F',
125+
'patch_icpp(2)%alter_patch(1)' : 'T',
126+
# ==========================================================================
127+
128+
# ==========================================================================
129+
#'patch_icpp(3)%geometry' : 21,
130+
#'patch_icpp(3)%x_centroid' : 3*L/2,
131+
#'patch_icpp(3)%y_centroid' : L/4,
132+
#'patch_icpp(3)%length_x' : L,
133+
#'patch_icpp(3)%length_y' : L/2,
134+
#'patch_icpp(3)%vel(1)' : 0,
135+
#'patch_icpp(3)%vel(2)' : 100,
136+
#'patch_icpp(3)%model_scale(1)' : 1/20.0,
137+
#'patch_icpp(3)%model_scale(2)' : 1/2.0,
138+
#'patch_icpp(3)%model_scale(3)' : 1/10.0,
139+
#'patch_icpp(3)%model_translate(1)' : L/2+L/1.2,
140+
#'patch_icpp(3)%model_translate(2)' : L/4,
141+
#'patch_icpp(3)%model_translate(3)' : 0,
142+
#'patch_icpp(3)%model_threshold' : 0.5,
143+
#'patch_icpp(3)%model_spc' : 'mfc-small.stl',
144+
#'patch_icpp(3)%pres' : 1000*sol_R.P,
145+
#'patch_icpp(3)%alpha(1)' : 1,
146+
#'patch_icpp(3)%alpha_rho(1)' : sol_R.density,
147+
#'patch_icpp(3)%model_filepath' : 'mfc-small.stl',
148+
#'patch_icpp(3)%alter_patch(1)' : 'T',
149+
#'patch_icpp(3)%alter_patch(2)' : 'T',
150+
# ==========================================================================
151+
152+
# 'vel_wrt(1)' : 'T',
153+
#'vel_wrt(2)' : 'T',
154+
# 'pres_wrt' : 'T',
155+
156+
# Fluids Physical Parameters ===============================================
157+
'fluid_pp(1)%gamma' : 1.0E+00/(4.4E+00-1.0E+00),
158+
'fluid_pp(1)%pi_inf' : 0,
159+
# ==========================================================================
160+
}
161+
162+
if args.chemistry:
163+
for i in range(len(sol_L.Y)):
164+
case[f'chem_wrt_Y({1})'] = 'T'
165+
case[f'patch_icpp(1)%Y({i+1})'] = sol_L.Y[i]
166+
# case[f'patch_icpp(2)%Y({i+1})'] = sol_L.Y[i]
167+
# #case[f'patch_icpp(3)%Y({i+1})'] = sol_L.Y[i]
168+
169+
if __name__ == '__main__':
170+
print(json.dumps(case))

‎src/pre_process/include/1dHardcodedIC.fpp‎

Lines changed: 31 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -1,18 +1,38 @@
11
#:def Hardcoded1DVariables()
2+
3+
!> @file 1DHardcodedIC.fpp
4+
!> @brief Initialize a simulation using data from a prior 1D run (for example from Cantera).
5+
!>
6+
!> @details
7+
!> Reads a 1D primitive‐variable dataset and sets the
8+
!> initial conditions for your simulation based on those 1D results.
9+
!>
10+
!> It provides:
11+
!> - Import of files named
12+
!> `prim.<variable>.00.<timestep>.dat`
13+
!> (produced when parallel I/O is disabled in `case.py`).
14+
!> - Configurable parameters:
15+
!> - @c nFiles: total number of primitive‐variable files.
16+
!> - @c nRows: total number of grid points in the 1D profile.
17+
!> - Default file directory: `examples/Case_File/D`
18+
!>
19+
!> @param nFiles Total number of primitive‐variable files to read.
20+
!> @param nRows Total number of grid points in the imported 1D profile.
21+
222
! Place any declaration of intermediate variables here
3-
integer, parameter :: nFiles = 14 ! Can be changed to any number
4-
integer, parameter :: nRows = 512!
23+
integer, parameter :: nFiles = 14 ! Number of files (variables) that are being read
24+
integer, parameter :: nRows = 512 ! Number of grid points
525
integer :: f, iter, ios, unit, idx
626
real(wp) :: x_len, x_step
7-
integer :: global_offset
8-
real(wp) :: delta
9-
character(len=100), dimension(nFiles) :: fileNames
10-
! Arrays to store all data from files
11-
real(wp), dimension(nRows, nFiles) :: stored_values
27+
integer :: global_offset ! MPI subdomain offset
28+
real(wp) :: delta
29+
character(len=100), dimension(nFiles) :: fileNames ! Arrays to store all data from files
30+
character(len=200) :: errmsg
31+
real(wp), dimension(nRows, nFiles) :: stored_values ! Imported Data
1232
real(wp), dimension(nRows) :: x_coords
1333
logical :: files_loaded = .false.
1434
real(wp) :: domain_start, domain_end
15-
character(len=*), parameter :: init_dir = "/home/YourDirectory"
35+
character(len=*), parameter :: init_dir = "/home/YourDirectory" ! For example /home/MFC/examples/1D_Shock/D/
1636
character(len=20) :: file_num_str ! For storing the file number as a string
1737
character(len=20) :: zeros_part ! For the trailing zeros part
1838
character(len=6), parameter :: zeros_default = "000000" ! Default zeros (can be changed)
@@ -37,21 +57,19 @@
3757
case (100)
3858
! Put your variable assignments here
3959
if (.not. files_loaded) then
40-
! Print status message
41-
print *, "Loading all data files..."
4260
do f = 1, nFiles
4361
! Open the file for reading
4462
open (newunit=unit, file=trim(fileNames(f)), status='old', action='read', iostat=ios)
4563
if (ios /= 0) then
46-
print *, "Error opening file: ", trim(fileNames(f))
4764
cycle ! Skip this file on error
4865
end if
4966
! Read all rows at once into memory
5067
do iter = 1, nRows
5168
read (unit, *, iostat=ios) x_coords(iter), stored_values(iter, f)
5269
if (ios /= 0) then
53-
print *, "Error reading file ", trim(fileNames(f)), " at row ", iter
54-
exit ! Exit loop on error
70+
write(errmsg, '(A,A,A,I0,A)') ' Error reading "', trim(fileNames(f)), &
71+
'" at index (', iter, ')'
72+
call s_mpi_abort(trim(errmsg)) ! Exit loop on error
5573
end if
5674
end do
5775
close (unit)
@@ -63,7 +81,6 @@
6381

6482
delta = x_cc(0) - domain_start - x_step/2.0
6583
global_offset = nint(abs(delta)/x_step)
66-
print *, "All data files loaded. Domain range:", domain_start, "to", domain_end
6784
files_loaded = .true.
6885
end if
6986
! Simple mapping - find the closest index

0 commit comments

Comments
 (0)