Skip to content

Commit 658f188

Browse files
Information Geometric Regularization (#900)
Co-authored-by: Spencer Bryngelson <[email protected]>
1 parent 0ed69c5 commit 658f188

File tree

107 files changed

+10468
-841
lines changed

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

107 files changed

+10468
-841
lines changed

CMakeLists.txt

Lines changed: 9 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -135,17 +135,17 @@ if (CMAKE_Fortran_COMPILER_ID STREQUAL "GNU")
135135
if (CMAKE_BUILD_TYPE STREQUAL "Debug")
136136
add_compile_options(
137137
-Wall
138-
-Wextra
138+
-Wextra
139139
-fcheck=all,no-array-temps
140140
-fbacktrace
141141
-fimplicit-none
142142
-fsignaling-nans
143143
-finit-real=snan
144144
-finit-integer=-99999999
145-
-Wintrinsic-shadow
146-
-Wunderflow
147-
-Wrealloc-lhs
148-
-Wsurprising
145+
-Wintrinsic-shadow
146+
-Wunderflow
147+
-Wrealloc-lhs
148+
-Wsurprising
149149
)
150150
endif()
151151

@@ -163,6 +163,7 @@ elseif (CMAKE_Fortran_COMPILER_ID STREQUAL "Cray")
163163
"SHELL:-h acc_model=auto_async_none"
164164
"SHELL: -h acc_model=no_fast_addr"
165165
"SHELL: -h list=adm"
166+
"SHELL: -munsafe-fp-atomics" # Not unsafe for operations we do
166167
)
167168

168169
add_link_options("SHELL:-hkeepfiles")
@@ -172,7 +173,6 @@ elseif (CMAKE_Fortran_COMPILER_ID STREQUAL "Cray")
172173
"SHELL:-h acc_model=auto_async_none"
173174
"SHELL: -h acc_model=no_fast_addr"
174175
"SHELL: -K trap=fp" "SHELL: -G2"
175-
176176
)
177177
add_link_options("SHELL: -K trap=fp" "SHELL: -G2")
178178
endif()
@@ -200,10 +200,10 @@ elseif ((CMAKE_Fortran_COMPILER_ID STREQUAL "NVHPC") OR (CMAKE_Fortran_COMPILER_
200200
if (CMAKE_BUILD_TYPE STREQUAL "Debug")
201201
add_compile_options(
202202
$<$<COMPILE_LANGUAGE:Fortran>:-O0>
203-
$<$<COMPILE_LANGUAGE:Fortran>:-C>
203+
$<$<COMPILE_LANGUAGE:Fortran>:-C>
204204
$<$<COMPILE_LANGUAGE:Fortran>:-g>
205-
$<$<COMPILE_LANGUAGE:Fortran>:-traceback>
206-
$<$<COMPILE_LANGUAGE:Fortran>:-Minform=inform>
205+
$<$<COMPILE_LANGUAGE:Fortran>:-traceback>
206+
$<$<COMPILE_LANGUAGE:Fortran>:-Minform=inform>
207207
$<$<COMPILE_LANGUAGE:Fortran>:-Mbounds>
208208
)
209209
endif()

README.md

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -155,6 +155,8 @@ They are organized below.
155155
* Runge-Kutta orders 1-3 (SSP TVD), adaptive time stepping
156156
* RK4-5 operator splitting for Euler-Lagrange modeling
157157
* Interface sharpening (THINC-like)
158+
* Information geometric regularization (IGR)
159+
* Shock capturing without WENO and Riemann solvers
158160

159161
### Large-scale and accelerated simulation
160162

benchmarks/5eq_rk3_weno3_hllc/case.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -188,8 +188,8 @@
188188
"cyl_coord": "F",
189189
"dt": dt,
190190
"t_step_start": 0,
191-
"t_step_stop": int(30 * (95 * size + 5)),
192-
"t_step_save": int(30 * (95 * size + 5)),
191+
"t_step_stop": int(20 * (5 * size + 5)),
192+
"t_step_save": int(20 * (5 * size + 5)),
193193
# Simulation Algorithm Parameters
194194
"num_patches": 3,
195195
"model_eqns": 2,

benchmarks/hypo_hll/case.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -41,8 +41,8 @@
4141
"p": Nz,
4242
"dt": 1e-8,
4343
"t_step_start": 0,
44-
"t_step_stop": int(30 * (95 * size + 5)),
45-
"t_step_save": int(30 * (95 * size + 5)),
44+
"t_step_stop": int(20 * (5 * size + 5)),
45+
"t_step_save": int(20 * (5 * size + 5)),
4646
# Simulation Algorithm Parameters
4747
"num_patches": 2,
4848
"model_eqns": 2,

benchmarks/ibm/case.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -45,8 +45,8 @@
4545
"p": Nz,
4646
"dt": mydt,
4747
"t_step_start": 0,
48-
"t_step_stop": int(20 * (95 * size + 5)),
49-
"t_step_save": int(20 * (95 * size + 5)),
48+
"t_step_stop": int(20 * (5 * size + 5)),
49+
"t_step_save": int(20 * (5 * size + 5)),
5050
# Simulation Algorithm Parameters
5151
"num_patches": 1,
5252
"model_eqns": 2,

benchmarks/igr/case.py

Lines changed: 113 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,113 @@
1+
#!/usr/bin/env python3
2+
# Benchmark igr_T_viscous_T
3+
# Additional Benchmarked Features
4+
# - igr : T
5+
# - viscous : T
6+
# - igr_order : 5
7+
8+
import json, math, argparse
9+
10+
parser = argparse.ArgumentParser(prog="Benchmarking Case 5", description="This MFC case was created for the purposes of benchmarking MFC.", formatter_class=argparse.ArgumentDefaultsHelpFormatter)
11+
12+
parser.add_argument("--mfc", type=json.loads, default="{}", metavar="DICT", help="MFC's toolchain's internal state.")
13+
parser.add_argument("--gbpp", type=int, metavar="MEM", default=16, help="Adjusts the problem size per rank to fit into [MEM] GB of GPU memory per GPU.")
14+
15+
ARGS = vars(parser.parse_args())
16+
DICT = ARGS["mfc"]
17+
18+
size = 1 if DICT["gpu"] else 0
19+
20+
ppg = 8000000 / 16.0
21+
procs = DICT["nodes"] * DICT["tasks_per_node"]
22+
ncells = math.floor(ppg * procs * ARGS["gbpp"])
23+
s = math.floor((ncells) ** (1 / 3))
24+
Nx, Ny, Nz = s, s, s
25+
26+
Re = 1600
27+
L = 1
28+
P0 = 101325
29+
rho0 = 1
30+
C0 = math.sqrt(1.4 * P0)
31+
V0 = 0.1 * C0
32+
mu = V0 * L / Re
33+
34+
cfl = 0.5
35+
dx = 2 * math.pi * L / (Nx + 1)
36+
37+
dt = cfl * dx / (C0)
38+
39+
tC = L / V0
40+
tEnd = 20 * tC
41+
42+
Nt = int(tEnd / dt)
43+
44+
# Configuring case dictionary
45+
print(
46+
json.dumps(
47+
{
48+
# Logistics
49+
"run_time_info": "T",
50+
# Computational Domain Parameters
51+
"x_domain%beg": -math.pi * L,
52+
"x_domain%end": math.pi * L,
53+
"y_domain%beg": -math.pi * L,
54+
"y_domain%end": math.pi * L,
55+
"z_domain%beg": -math.pi * L,
56+
"z_domain%end": math.pi * L,
57+
"m": Nx,
58+
"n": Ny,
59+
"p": Nz,
60+
"cyl_coord": "F",
61+
"dt": dt,
62+
"t_step_start": 0,
63+
"t_step_stop": int(20 * (5 * size + 5)),
64+
"t_step_save": int(20 * (5 * size + 5)),
65+
# Simulation Algorithm Parameters
66+
"num_patches": 1,
67+
"model_eqns": 2,
68+
"num_fluids": 1,
69+
"time_stepper": 3,
70+
"bc_x%beg": -1,
71+
"bc_x%end": -1,
72+
"bc_y%beg": -1,
73+
"bc_y%end": -1,
74+
"bc_z%beg": -1,
75+
"bc_z%end": -1,
76+
"igr": "T",
77+
"igr_order": 5,
78+
"igr_iter_solver": 1,
79+
"num_igr_iters": 3,
80+
"num_igr_warm_start_iters": 3,
81+
"alf_factor": 10,
82+
"viscous": "T",
83+
# Formatted Database Files Structure Parameters
84+
"format": 1,
85+
"precision": 2,
86+
"prim_vars_wrt": "T",
87+
"omega_wrt(1)": "T",
88+
"omega_wrt(2)": "T",
89+
"omega_wrt(3)": "T",
90+
"qm_wrt": "T",
91+
"fd_order": 4,
92+
"parallel_io": "T",
93+
# Patch 1: Background (AIR - 2)
94+
"patch_icpp(1)%geometry": 9,
95+
"patch_icpp(1)%x_centroid": 0,
96+
"patch_icpp(1)%y_centroid": 0,
97+
"patch_icpp(1)%z_centroid": 0,
98+
"patch_icpp(1)%length_x": 2 * math.pi * L,
99+
"patch_icpp(1)%length_y": 2 * math.pi * L,
100+
"patch_icpp(1)%length_z": 2 * math.pi * L,
101+
"patch_icpp(1)%vel(1)": f"{V0}*sin(x/{L})*cos(y/{L})*sin(z/{L})",
102+
"patch_icpp(1)%vel(2)": f"-{V0}*cos(x/{L})*sin(y/{L})*sin(z/{L})",
103+
"patch_icpp(1)%vel(3)": 0,
104+
"patch_icpp(1)%pres": f"{P0} + ({rho0}*{V0}**2/16)*(cos(2*x/{L}) + cos(2*y/{L}))*(cos(2*z/{L}) + 2)",
105+
"patch_icpp(1)%alpha_rho(1)": 1,
106+
"patch_icpp(1)%alpha(1)": 1,
107+
# Fluids Physical Parameters
108+
"fluid_pp(1)%gamma": 1.0e00 / (1.4 - 1),
109+
"fluid_pp(1)%pi_inf": 0,
110+
"fluid_pp(1)%Re(1)": 1 / mu,
111+
}
112+
)
113+
)

benchmarks/viscous_weno5_sgb_acoustic/case.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -107,8 +107,8 @@
107107
"p": Nz,
108108
"dt": dt,
109109
"t_step_start": 0,
110-
"t_step_stop": int(15 * (25 * size + 5)),
111-
"t_step_save": int(15 * (25 * size + 5)),
110+
"t_step_stop": int(15 * (5 * size + 5)),
111+
"t_step_save": int(15 * (5 * size + 5)),
112112
# Simulation Algorithm Parameters
113113
"num_patches": 2,
114114
"model_eqns": 2,

docs/documentation/case.md

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -409,6 +409,13 @@ Details of implementation of viscosity in MFC can be found in [Coralic (2015)](r
409409
| `surface_tension` | Logical | Activate surface tension |
410410
| `viscous` | Logical | Activate viscosity |
411411
| `hypoelasticity` | Logical | Activate hypoelasticity* |
412+
| `igr` | Logical | Enable solution via information geometric regularization (IGR) [Cao (2024)](references.md) |
413+
| `igr_order` | Integer | Order of reconstruction for IGR [3,5] |
414+
| `alf_factor` | Real | Alpha factor for IGR entropic pressure (default 10) |
415+
| `igr_pres_lim` | Logical | Limit IGR pressure to avoid negative values (default F) |
416+
| `igr_iter_solver` | Integer | Solution method for IGR elliptic solve [1] Jacobi [2] Gauss-Seidel |
417+
| `num_igr_iters` | Integer | Number of iterations for for the IGR elliptic solve (default 2) |
418+
| `num_igr_warm_start_iters` | Integer | Number of iterations for the IGR elliptic solve at the first time step (default 50) |
412419

413420
- \* Options that work only with `model_eqns = 2`.
414421
- † Options that work only with ``cyl_coord = 'F'``.
@@ -502,7 +509,7 @@ This option requires `weno_Re_flux` to be true because cell boundary values are
502509
| `type`* | Integer | The geometry of the patch. [1]: Line [2]: Circle [3]: Rectangle |
503510
| `x[y,z]_centroid`* | Real | Centroid of the boundary patch in the x[y,z]-direction |
504511
| `length_x[y,z]`* | Real | Length of the boundary patch in the x[y,z]-direction |
505-
| `radius`* | Real | Radius of the boundary patch |
512+
| `radiue`* | Real | Radius of the boundary patch |
506513
*: These parameters should be prepended with `patch_bc(j)%` where $j$ is the patch index.
507514

508515
Boundary condition patches can be used with the following boundary condition types:

docs/documentation/references.md

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,8 @@
1010

1111
- <a id="Bryngelson19">Bryngelson, S. H., Schmidmayer, K., Coralic, V., Meng, J. C., Maeda, K., and Colonius, T. (2019). Mfc: An open-source high-order multi-component, multi-phase, and multi-scale compressible flow solver. arXiv preprint arXiv:1907.10512.</a>
1212

13+
- <a id="Cao24">Cao A. and. Sch&auml;fer F. (2024). Information Geometric Regularization of the Barotropic Euler Equations. arXiv preprint arXiv:2308.14127.</a>
14+
1315
- <a id="Chen22">Chen, S. S., Li, J. P., Li, Z., Yuan, W., & Gao, Z. H. (2022). Anti-dissipation pressure correction under low Mach numbers for Godunov-type schemes. Journal of Computational Physics, 456, 111027. </a>
1416

1517
- <a id="Childs12">Childs, H., Brugger, E., Whitlock, B., Meredith, J., Ahern, S., Pugmire, D., Biagas, K., Miller, M., Harrison, C., Weber, G. H., Krishnan, H., Fogal, T., Sanderson, A., Garth, C., Bethel, E. W., Camp, D., R¨ubel, O., Durant, M., Favre, J. M., and Navr´atil, P. (2012). VisIt: An End-User Tool For Visualizing and Analyzing Very Large Data. In High Performance Visualization–Enabling Extreme-Scale Scientific Insight, pages 357–372.</a>
@@ -66,4 +68,4 @@
6668

6769
- <a id="Powell94">Powell, K. G. (1994). An approximate Riemann solver for magnetohydrodynamics: (That works in more than one dimension). In Upwind and high-resolution schemes (pp. 570-583). Springer.</a>
6870

69-
- <a id="Cao19">Cao, S., Zhang, Y., Liao, D., Zhong, P., and Wang, K. G. (2019). Shock-induced damage and dynamic fracture in cylindrical bodies submerged in liquid. International Journal of Solids and Structures, 169:55–71. Elsevier.</a>
71+
- <a id="Cao19">Cao, S., Zhang, Y., Liao, D., Zhong, P., and Wang, K. G. (2019). Shock-induced damage and dynamic fracture in cylindrical bodies submerged in liquid. International Journal of Solids and Structures, 169:55–71. Elsevier.</a>

examples/2D_IGR_2fluid/case.py

Lines changed: 105 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,105 @@
1+
2+
#!/usr/bin/env python3
3+
# This case file demonstrates the Laplace pressure jump of a water droplet in air. The laplace pressure jump
4+
# in 2D is given by delta = sigma / r where delta is the pressure jump, sigma is the surface tension coefficient,
5+
# and r is the radius of the droplet. The results of this simulation agree with theory to well within 1%
6+
# relative error.
7+
8+
import math
9+
import json
10+
11+
l = 1
12+
eps = 1e-6
13+
14+
# Numerical setup
15+
r0 = l / 4
16+
x0 = -l
17+
x1 = l
18+
y0 = -l
19+
y1 = l
20+
21+
Nx = 200
22+
Ny = 200
23+
24+
dt = 1e-4
25+
26+
# Configuration case dictionary
27+
data = {
28+
# Logistics
29+
"run_time_info": "T",
30+
# Computational Domain
31+
"x_domain%beg": x0,
32+
"x_domain%end": x1,
33+
"y_domain%beg": y0,
34+
"y_domain%end": y1,
35+
"m": Nx,
36+
"n": Ny,
37+
"p": 0,
38+
"cyl_coord": "F",
39+
"dt": dt,
40+
"t_step_start": 0,
41+
"t_step_stop": 5000,
42+
"t_step_save": 50,
43+
# Simulation Algorithm
44+
"model_eqns": 2,
45+
"alt_soundspeed": "F",
46+
"mixture_err": "T",
47+
"mpp_lim": "F",
48+
"time_stepper": 3,
49+
"bc_x%beg": -1,
50+
"bc_x%end": -1,
51+
"bc_y%beg": -1,
52+
"bc_y%end": -1,
53+
"num_patches": 2,
54+
"num_fluids": 2,
55+
"elliptic_smoothing": "T",
56+
"elliptic_smoothing_iters": 10,
57+
"igr": "T",
58+
"igr_order": 5,
59+
"igr_iter_solver": 1,
60+
"num_igr_iters": 3,
61+
"num_igr_warm_start_iters": 30,
62+
"alf_factor": 10,
63+
"viscous": "T",
64+
# Database Structure Parameters
65+
"format": 1,
66+
"precision": 2,
67+
"prim_vars_wrt": "T",
68+
"parallel_io": "T",
69+
# Fluid Parameters (Gas)
70+
"fluid_pp(1)%gamma": 1.0e00 / (1.4e00 - 1.0e00),
71+
"fluid_pp(1)%pi_inf": 0.0e00,
72+
"fluid_pp(1)%Re(1)": 1e5,
73+
"fluid_pp(2)%gamma": 1.0e00 / (1.4e00 - 1.0e00),
74+
"fluid_pp(2)%pi_inf": 0.0e00,
75+
"fluid_pp(2)%Re(1)": 1e5,
76+
# Ambient pressure
77+
"patch_icpp(1)%geometry": 3,
78+
"patch_icpp(1)%x_centroid": 0,
79+
"patch_icpp(1)%y_centroid": 0,
80+
"patch_icpp(1)%length_x": 2,
81+
"patch_icpp(1)%length_y": 2,
82+
"patch_icpp(1)%vel(1)": 0.0,
83+
"patch_icpp(1)%vel(2)": 0.0,
84+
"patch_icpp(1)%vel(3)": 0.0,
85+
"patch_icpp(1)%pres": 1,
86+
"patch_icpp(1)%alpha_rho(1)": 1 - eps,
87+
"patch_icpp(1)%alpha(1)": 1 - eps,
88+
"patch_icpp(1)%alpha_rho(2)": eps,
89+
"patch_icpp(1)%alpha(2)": eps,
90+
# High pressure
91+
"patch_icpp(2)%alter_patch(1)": "T",
92+
"patch_icpp(2)%geometry": 2,
93+
"patch_icpp(2)%x_centroid": 0,
94+
"patch_icpp(2)%y_centroid": 0,
95+
"patch_icpp(2)%radius": r0,
96+
"patch_icpp(2)%vel(1)": 0.0,
97+
"patch_icpp(2)%vel(2)": 0.0,
98+
"patch_icpp(2)%pres": 9.518,
99+
"patch_icpp(2)%alpha_rho(1)": eps,
100+
"patch_icpp(2)%alpha(1)": eps,
101+
"patch_icpp(2)%alpha_rho(2)": (1 - eps) * 5,
102+
"patch_icpp(2)%alpha(2)": 1 - eps,
103+
}
104+
105+
print(json.dumps(data))

0 commit comments

Comments
 (0)