Skip to content

Commit 08d12f8

Browse files
danieljvickerssbryngelsonDaniel VickersDaniel VickersDaniel Vickers
authored
Local Aware IBM (#1378)
Co-authored-by: Spencer Bryngelson <sbryngelson@gmail.com> Co-authored-by: Daniel Vickers <danieljvickers@login05.frontier.olcf.ornl.gov> Co-authored-by: Daniel Vickers <danieljvickers@login07.frontier.olcf.ornl.gov> Co-authored-by: Daniel Vickers <danieljvickers@login09.frontier.olcf.ornl.gov> Co-authored-by: Daniel Vickers <danieljvickers@login08.frontier.olcf.ornl.gov> Co-authored-by: Daniel Vickers <danieljvickers@frontier10305.frontier.olcf.ornl.gov>
1 parent 40dde5e commit 08d12f8

29 files changed

Lines changed: 1716 additions & 381 deletions

docs/documentation/case.md

Lines changed: 6 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -312,6 +312,8 @@ This is enabled by adding ``'elliptic_smoothing': "T",`` and ``'elliptic_smoothi
312312
| Parameter | Type | Description |
313313
| ---: | :----: | :--- |
314314
| `num_ibs` | Integer | Number of immersed boundary patches |
315+
| `num_particle_beds` | Integer | Number of particle bed specifications to generate immersed boundary patches from |
316+
| `ib_neighborhood_radius` | Integer | Parameter that controls the neighborhood size for IB detection. |
315317
| `geometry` | Integer | Geometry configuration of the patch.|
316318
| `x[y,z]_centroid` | Real | Centroid of the applied geometry in the [x,y,z]-direction. |
317319
| `length_x[y,z]` | Real | Length, if applicable, in the [x,y,z]-direction. |
@@ -373,6 +375,8 @@ Additional details on this specification can be found in [NACA airfoil](https://
373375

374376
- `ib_coefficient_of_friction` is the coefficient of friction used in IB collisions.
375377

378+
- `ib_neighborhood_radius` controls the size of the neighborhood size. This value defaults to 1, which indicates that any given rank is aware of IB's up to 1 ranks away. This parameter is required to strong-scale a case when IB's eventually grow to be larger than one full processor domain wide.
379+
376380
### 5. Fluid Material's {#sec-fluid-materials}
377381

378382
| Parameter | Type | Description |
@@ -643,7 +647,7 @@ To restart the simulation from $k$-th time step, see @ref running "Restarting Ca
643647
| `alpha_wrt(i)` | Logical | Add the volume fraction of fluid $i$ to the database |
644648
| `gamma_wrt` | Logical | Add the specific heat ratio function to the database |
645649
| `heat_ratio_wrt` | Logical | Add the specific heat ratio to the database |
646-
| `ib_state_wrt` | Logical | Write IB state and loads to a datafile at each time step |
650+
| `ib_state_wrt` | Logical | Parameter to handle writing IB state on saves and outputting the state as a point mesh to SILO files. |
647651
| `pi_inf_wrt` | Logical | Add the liquid stiffness function to the database |
648652
| `pres_inf_wrt` | Logical | Add the liquid stiffness to the formatted database |
649653
| `c_wrt` | Logical | Add the sound speed to the database |
@@ -711,7 +715,7 @@ If `file_per_process` is true, then pre_process, simulation, and post_process mu
711715

712716
- `probe_wrt` activates the output of state variables at coordinates specified by `probe(i)%[x;y,z]`.
713717

714-
- `ib_state_wrt` activates the output of data specified by patch_ib(i)%force(:) (and torque, vel, angular_vel, angles, [x,y,z]_centroid) into a single binary datafile for all IBs at all timesteps. During post_processing, this file is converted into separate time histories for each IB.
718+
- `ib_state_wrt` is used to trigger post-processing of the IB state to be written out as a point mesh in the SILO files. When no IBs are moving, it also triggers force and torque calculation so that those values may be written to the output state files.
715719

716720
- `output_partial_domain` activates the output of part of the domain specified by `[x,y,z]_output%%beg` and `[x,y,z]_output%%end`.
717721
This is useful for large domains where only a portion of the domain is of interest.

docs/module_categories.json

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -38,6 +38,7 @@
3838
"m_compute_cbc",
3939
"m_boundary_common",
4040
"m_ibm",
41+
"m_particle_bed",
4142
"m_igr",
4243
"m_ib_patches",
4344
"m_compute_levelset"
Lines changed: 136 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,136 @@
1+
import json
2+
import math
3+
4+
# 2D shock wave interacting with a bed of 20 free-floating circular particles.
5+
6+
gam_a = 1.4
7+
8+
# Shock parameters (Mach 1.5)
9+
mach_number = 1.5
10+
pre_shock_pressure = 1
11+
pre_shock_density = 1.4
12+
pre_shock_speed = 0.0
13+
post_shock_pressure = 2.4583
14+
post_shock_density = 2.6069
15+
post_shock_speed = 0.6944
16+
17+
domain_size = 4.0
18+
wave_front = -1.5
19+
20+
total_time = 1.5
21+
num_time_steps = 2000
22+
dt = float(total_time / num_time_steps)
23+
num_saves = 100
24+
steps_to_save = int(num_time_steps / num_saves)
25+
26+
# Soft-sphere collision parameters (from 3D_mibm_sphere_head_on_collision)
27+
collision_time = 20.0 * dt
28+
29+
# Particle bed parameters
30+
bed_x = 0.5
31+
bed_y = 0.0
32+
bed_lx = 2.0
33+
bed_ly = 3.5
34+
particle_radius = 0.15
35+
particle_mass = 0.25
36+
particle_min_spacing = 0.05
37+
38+
print(
39+
json.dumps(
40+
{
41+
# Logistics
42+
"run_time_info": "T",
43+
# Computational Domain Parameters
44+
"x_domain%beg": -domain_size * 0.5,
45+
"x_domain%end": domain_size * 0.5,
46+
"y_domain%beg": -domain_size * 0.5,
47+
"y_domain%end": domain_size * 0.5,
48+
"cyl_coord": "F",
49+
"m": 256,
50+
"n": 256,
51+
"p": 0,
52+
"dt": dt,
53+
"t_step_start": 0,
54+
"t_step_stop": num_time_steps,
55+
"t_step_save": steps_to_save,
56+
# Simulation Algorithm Parameters
57+
"num_patches": 2,
58+
"model_eqns": 2,
59+
"alt_soundspeed": "F",
60+
"num_fluids": 1,
61+
"mpp_lim": "F",
62+
"mixture_err": "T",
63+
"time_stepper": 3,
64+
"weno_order": 5,
65+
"weno_eps": 1.0e-16,
66+
"weno_Re_flux": "T",
67+
"weno_avg": "T",
68+
"avg_state": 2,
69+
"mapped_weno": "T",
70+
"null_weights": "F",
71+
"mp_weno": "T",
72+
"riemann_solver": 2,
73+
"wave_speeds": 1,
74+
"bc_x%beg": -17,
75+
"bc_x%end": -8,
76+
"bc_y%beg": -15,
77+
"bc_y%end": -15,
78+
# Immersed boundaries — all circles come from the particle bed
79+
"ib": "T",
80+
"num_ibs": 0,
81+
"viscous": "T",
82+
# Collision model (soft-sphere, from 3D_mibm_sphere_head_on_collision)
83+
"collision_model": 1,
84+
"coefficient_of_restitution": 0.9,
85+
"collision_time": collision_time,
86+
"ib_coefficient_of_friction": 0.1,
87+
# Particle bed: 20 free-floating circles placed randomly in region
88+
"num_particle_beds": 1,
89+
"particle_bed(1)%x_centroid": bed_x,
90+
"particle_bed(1)%y_centroid": bed_y,
91+
"particle_bed(1)%z_centroid": 0.0,
92+
"particle_bed(1)%length_x": bed_lx,
93+
"particle_bed(1)%length_y": bed_ly,
94+
"particle_bed(1)%length_z": 0.0,
95+
"particle_bed(1)%num_particles": 20,
96+
"particle_bed(1)%radius": particle_radius,
97+
"particle_bed(1)%mass": particle_mass,
98+
"particle_bed(1)%min_spacing": particle_min_spacing,
99+
"particle_bed(1)%moving_ibm": 2,
100+
"particle_bed(1)%seed": 42,
101+
# Output
102+
"format": 1,
103+
"precision": 2,
104+
"prim_vars_wrt": "T",
105+
"E_wrt": "T",
106+
"ib_state_wrt": "F",
107+
"parallel_io": "T",
108+
# IC Patch 1: post-shock region (left of wave front)
109+
"patch_icpp(1)%geometry": 3,
110+
"patch_icpp(1)%x_centroid": 0.5 * wave_front - 0.25 * domain_size,
111+
"patch_icpp(1)%y_centroid": 0.0,
112+
"patch_icpp(1)%length_x": 0.5 * domain_size + wave_front,
113+
"patch_icpp(1)%length_y": domain_size,
114+
"patch_icpp(1)%vel(1)": post_shock_speed,
115+
"patch_icpp(1)%vel(2)": 0.0,
116+
"patch_icpp(1)%pres": post_shock_pressure,
117+
"patch_icpp(1)%alpha_rho(1)": post_shock_density,
118+
"patch_icpp(1)%alpha(1)": 1.0,
119+
# IC Patch 2: pre-shock region (right of wave front)
120+
"patch_icpp(2)%geometry": 3,
121+
"patch_icpp(2)%x_centroid": 0.5 * wave_front + 0.25 * domain_size,
122+
"patch_icpp(2)%y_centroid": 0.0,
123+
"patch_icpp(2)%length_x": 0.5 * domain_size - wave_front,
124+
"patch_icpp(2)%length_y": domain_size,
125+
"patch_icpp(2)%vel(1)": pre_shock_speed,
126+
"patch_icpp(2)%vel(2)": 0.0,
127+
"patch_icpp(2)%pres": pre_shock_pressure,
128+
"patch_icpp(2)%alpha_rho(1)": pre_shock_density,
129+
"patch_icpp(2)%alpha(1)": 1.0,
130+
# Fluid properties: air
131+
"fluid_pp(1)%gamma": 1.0 / (gam_a - 1.0),
132+
"fluid_pp(1)%pi_inf": 0,
133+
"fluid_pp(1)%Re(1)": 2500000,
134+
}
135+
)
136+
)

examples/2D_mibm_shock_cylinder/case.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -87,7 +87,7 @@
8787
"precision": 2,
8888
"prim_vars_wrt": "T",
8989
"E_wrt": "T",
90-
"ib_state_wrt": "T",
90+
"ib_state_wrt": "F",
9191
"parallel_io": "T",
9292
# Patch: Constant Tube filled with air
9393
# Specify the cylindrical air tube grid geometry

src/common/include/case.fpp

Lines changed: 0 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -4,10 +4,8 @@
44

55
! For pre-process.
66
#:def analytical()
7-
87
#:enddef
98

109
! For moving immersed boundaries in simulation
1110
#:def mib_analytical()
12-
1311
#:enddef

src/common/m_constants.fpp

Lines changed: 14 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -16,15 +16,20 @@ module m_constants
1616
real(wp), parameter :: verysmall = 1.e-12_wp !< Very small number
1717
!> Radius cutoff to avoid division by zero for 3D spherical harmonic patch (geometry 14)
1818
real(wp), parameter :: small_radius = 1.e-32_wp
19-
integer, parameter :: num_stcls_min = 5 !< Minimum # of stencils
20-
integer, parameter :: path_len = 400 !< Maximum path length
21-
integer, parameter :: name_len = 50 !< Maximum name length
22-
integer, parameter :: dflt_int = -100 !< Default integer value
23-
integer, parameter :: fourier_rings = 5 !< Fourier filter ring limit
24-
integer, parameter :: num_fluids_max = 10 !< Maximum number of fluids in the simulation
25-
integer, parameter :: num_probes_max = 10 !< Maximum number of flow probes in the simulation
26-
integer, parameter :: num_patches_max = 10 !< Maximum number of IC patches
27-
integer, parameter :: num_ib_patches_max = 50000 !< Maximum number of immersed boundary patches (patch_ib)
19+
integer, parameter :: num_stcls_min = 5 !< Minimum # of stencils
20+
integer, parameter :: path_len = 400 !< Maximum path length
21+
integer, parameter :: name_len = 50 !< Maximum name length
22+
integer, parameter :: dflt_int = -100 !< Default integer value
23+
integer, parameter :: fourier_rings = 5 !< Fourier filter ring limit
24+
integer, parameter :: num_fluids_max = 10 !< Maximum number of fluids in the simulation
25+
integer, parameter :: num_probes_max = 10 !< Maximum number of flow probes in the simulation
26+
integer, parameter :: num_patches_max = 10 !< Maximum number of IC patches
27+
!> Maximum number of immersed boundary patches (legacy, not used for patch_ib sizing)
28+
integer, parameter :: num_ib_patches_max = 2050000
29+
!> Fixed capacity of patch_ib (namelist patches + local particle bed subset after reduction)
30+
integer, parameter :: num_ib_patches_max_namelist = 54000
31+
integer, parameter :: num_local_ibs_max = 2000 !< Maximum number of immersed boundary patches (patch_ib)
32+
integer, parameter :: num_particle_beds_max = 10 !< Maximum number of particle bed patch specifications
2833
integer, parameter :: num_bc_patches_max = 10 !< Maximum number of boundary condition patches
2934
integer, parameter :: max_2d_fourier_modes = 10 !< Max Fourier mode index for 2D modal patch (geometry 13)
3035
integer, parameter :: max_sph_harm_degree = 5 !< Max degree L for 3D spherical harmonic patch (geometry 14)

src/common/m_derived_types.fpp

Lines changed: 17 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -202,12 +202,10 @@ module m_derived_types
202202

203203
type :: t_model_array
204204
! Original CPU-side fields (unchanged)
205-
type(t_model), allocatable :: model !< STL/OBJ geometry model
206-
real(wp), allocatable, dimension(:,:,:) :: boundary_v !< Boundary vertices
207-
real(wp), allocatable, dimension(:,:) :: interpolated_boundary_v !< Interpolated boundary vertices
208-
integer :: boundary_edge_count !< Number of boundary edges
209-
integer :: total_vertices !< Total vertex count
210-
integer :: interpolate !< Interpolation flag
205+
type(t_model), allocatable :: model !< STL/OBJ geometry model
206+
real(wp), allocatable, dimension(:,:,:) :: boundary_v !< Boundary vertices
207+
integer :: boundary_edge_count !< Number of boundary edges
208+
integer :: total_vertices !< Total vertex count
211209

212210
! GPU-friendly flattened arrays
213211
integer :: ntrs !< Copy of model%ntrs
@@ -272,9 +270,10 @@ module m_derived_types
272270
end type ic_patch_parameters
273271

274272
type ib_patch_parameters
275-
276273
integer :: geometry !< Type of geometry for the patch
274+
integer :: gbl_patch_id
277275
real(wp) :: x_centroid, y_centroid, z_centroid !< Geometric center coordinates of the patch
276+
278277
!> Centroid locations of intermediate steps in the time_stepper module
279278
real(wp) :: step_x_centroid, step_y_centroid, step_z_centroid
280279
real(wp), dimension(1:3) :: centroid_offset !< offset of center of mass from computed cell center for odd-shaped IBs
@@ -307,6 +306,17 @@ module m_derived_types
307306
real(wp), dimension(1:3) :: step_angular_vel !< velocity array used to store intermediate steps in the time_stepper module
308307
end type ib_patch_parameters
309308

309+
type particle_bed_parameters
310+
real(wp) :: x_centroid, y_centroid, z_centroid !< Center of the particle bed region
311+
real(wp) :: length_x, length_y, length_z !< Dimensions of the particle bed region
312+
integer :: num_particles !< Number of particles to generate
313+
real(wp) :: radius !< Particle radius
314+
real(wp) :: mass !< Particle mass
315+
real(wp) :: min_spacing !< Minimum surface-to-surface gap (particle centers are 2*radius + min_spacing apart)
316+
integer :: moving_ibm !< Motion flag: 0=static, 1=moving (forces), 2=forced path
317+
integer :: seed !< Random seed for reproducible placement
318+
end type particle_bed_parameters
319+
310320
!> Derived type annexing the physical parameters (PP) of the fluids. These include the specific heat ratio function and liquid
311321
!! stiffness function.
312322
type physical_parameters

src/common/m_model.fpp

Lines changed: 8 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -983,12 +983,18 @@ contains
983983
dx_local = minval(dx); dy_local = minval(dy)
984984
if (p /= 0) dz_local = minval(dz)
985985

986+
! AMD HSA rejects zero-size GPU allocations; only create GPU mapping when there are actual models.
987+
if (num_ibs > 0) then
988+
@:ALLOCATE(models(num_ibs))
989+
else
990+
allocate (models(0))
991+
end if
986992
allocate (stl_bounding_boxes(num_ibs,1:3,1:3))
987993

988994
do patch_id = 1, num_ibs
989995
if (patch_ib(patch_id)%geometry == 5 .or. patch_ib(patch_id)%geometry == 12) then
990996
allocate (models(patch_id)%model)
991-
print *, " * Reading model: " // trim(patch_ib(patch_id)%model_filepath)
997+
if (proc_rank == 0) print *, " * Reading model: " // trim(patch_ib(patch_id)%model_filepath)
992998

993999
model = f_model_read(patch_ib(patch_id)%model_filepath)
9941000
params%scale(:) = patch_ib(patch_id)%model_scale(:)
@@ -1001,9 +1007,7 @@ contains
10011007
params%scale(:) = 1._wp
10021008
end if
10031009

1004-
if (proc_rank == 0) then
1005-
print *, " * Transforming model."
1006-
end if
1010+
if (proc_rank == 0) print *, " * Transforming model."
10071011

10081012
! Get the model center before transforming the model
10091013
bbox_old = f_create_bbox(model)

src/common/m_mpi_common.fpp

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -94,6 +94,8 @@ contains
9494
proc_rank = 0
9595
#endif
9696

97+
$:GPU_UPDATE(device='[num_procs, proc_rank]')
98+
9799
end subroutine s_mpi_initialize
98100

99101
!> Set up MPI I/O data views and variable pointers for parallel file output.

0 commit comments

Comments
 (0)