Skip to content

Periodic Immersed Boundary Wrapping for Circular and Spherical Geometries #1117

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
conraddelgado wants to merge 8 commits into
base: master
Choose a base branch
from
+1,076 −70
Open

Periodic Immersed Boundary Wrapping for Circular and Spherical Geometries #1117

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
8906810
periodic ib wrapping for spheres
conraddelgado Jan 28, 2026
837b18b
periodic ib wrap for circle geometry + tests
conraddelgado Jan 28, 2026
3441ef4
Merge branch 'MFlowCode:master' into master
conraddelgado Jan 28, 2026
6338eea
formatting
conraddelgado Jan 28, 2026
1aa4b08
Merge branch 'master' of github.com:conraddelgado/periodic_ibs
conraddelgado Jan 28, 2026
8571caf
regenerate test data, gpu update
conraddelgado Jan 30, 2026
85903b3
formatting
conraddelgado Jan 30, 2026
067a337
r2 on GPU
conraddelgado Jan 30, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
141 changes: 141 additions & 0 deletions examples/2D_ibm_circle_periodic/case.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,141 @@
import json
import math
import numpy as np

gam_a = 1.4

# sphere diameter
D = 0.1

# circle density
rho_c = 200

# domain length
L = 10 * D

# mach number
M = 1.2

# reynolds number
Re = 400.0

# pressure
P = 101325

# density
rho = 1.0

# fluid velocity
v1 = M * np.sqrt(gam_a * P / rho)

# dynamic viscosity
mu = rho * v1 * D / Re

dt = 1.0e-06
Nt = 5000
t_save = 20

# to fully resolve requires 40-60 cells across sphere diameter
Nx = 511
Ny = Nx
Comment on lines +38 to +40
Copy link
Contributor

@coderabbitai coderabbitai bot Jan 30, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🟡 Minor

Minor doc fix: “sphere” → “circle.”

This is a 2D circle case, so the resolution comment should reference a circle.

✏️ Suggested edit 
-# to fully resolve requires 40-60 cells across sphere diameter
+# to fully resolve requires 40-60 cells across circle diameter
📝 Committable suggestion

‼️ IMPORTANT
Carefully review the code before committing. Ensure that it accurately replaces the highlighted code, contains no missing lines, and has no issues with indentation. Thoroughly test & benchmark the code to ensure it meets the requirements.

Suggested change
# to fully resolve requires 40-60 cells across sphere diameter
Nx = 511
Ny = Nx
# to fully resolve requires 40-60 cells across circle diameter
Nx = 511
Ny = Nx
🤖 Prompt for AI Agents 
In `@examples/2D_ibm_circle_periodic/case.py` around lines 38 - 40, The inline
comment above the grid size (near Nx and Ny variables) incorrectly refers to a
"sphere"; update that comment to say "circle" since this is a 2D case (e.g.,
change the comment string "# to fully resolve requires 40-60 cells across sphere
diameter" to use "circle" instead). Ensure the wording remains accurate and
brief and that it sits immediately above or next to the Nx/Ny declarations
(symbols: Nx, Ny).

Nz = 0

# immersed boundary dictionary
ib_dict = {}
ib_dict.update(
{
f"patch_ib({1})%geometry": 2,
f"patch_ib({1})%x_centroid": 0.5 * L,
f"patch_ib({1})%y_centroid": 0.5 * L,
f"patch_ib({1})%radius": D / 2,
f"patch_ib({1})%slip": "F",
f"patch_ib({1})%moving_ibm": 0,
f"patch_ib({1})%mass": rho_c * np.pi * (D / 2.0) ** 2,
}
)

# Configuring case dictionary
case_dict = {
# Logistics
"run_time_info": "T",
# Computational Domain Parameters
# x direction
"x_domain%beg": -5.0 * D,
"x_domain%end": 5.0 * D,
# y direction
"y_domain%beg": -5.0 * D,
"y_domain%end": 5.0 * D,
"cyl_coord": "F",
"m": Nx,
"n": Ny,
"p": Nz,
"dt": dt,
"t_step_start": 0,
"t_step_stop": Nt,
"t_step_save": t_save,
# Simulation Algorithm Parameters
# Only one patch is necessary for one fluid
"num_patches": 1,
# Use the 5 equation model
"model_eqns": 2,
# 6 equations model does not need the K \div(u) term
"alt_soundspeed": "F",
# One fluids: air
"num_fluids": 1,
# time step
"mpp_lim": "F",
# Correct errors when computing speed of sound
"mixture_err": "T",
# Use TVD RK3 for time marching
"time_stepper": 3,
# Reconstruct the primitive variables to minimize spurious
# Use WENO5
"weno_order": 5,
"weno_eps": 1.0e-14,
"weno_Re_flux": "T",
"weno_avg": "T",
"avg_state": 2,
"mapped_weno": "T",
"null_weights": "F",
"mp_weno": "T",
"riemann_solver": 2,
"wave_speeds": 1,
# periodic bc
"bc_x%beg": -1,
"bc_x%end": -1,
"bc_y%beg": -1,
"bc_y%end": -1,
# Set IB to True and add 1 patch for every sphere
"ib": "T",
"num_ibs": 1,
"viscous": "T",
# Formatted Database Files Structure Parameters
"format": 1,
"precision": 2,
"prim_vars_wrt": "T",
"E_wrt": "T",
"parallel_io": "T",
# Patch: square filled with air
"patch_icpp(1)%geometry": 3,
# Uniform properties, centroid is at the center of the domain
"patch_icpp(1)%x_centroid": 0.0,
"patch_icpp(1)%y_centroid": 0.0,
"patch_icpp(1)%length_x": 10 * D,
"patch_icpp(1)%length_y": 10 * D,
# Specify the patch primitive variables
"patch_icpp(1)%vel(1)": v1,
"patch_icpp(1)%vel(2)": 0.0,
"patch_icpp(1)%pres": P,
"patch_icpp(1)%alpha_rho(1)": rho,
"patch_icpp(1)%alpha(1)": 1.0e00,
# Fluids Physical Parameters
"fluid_pp(1)%gamma": 1.0e00 / (gam_a - 1.0e00), # 2.50(Not 1.40)
"fluid_pp(1)%pi_inf": 0,
"fluid_pp(1)%Re(1)": 1.0 / mu,
# ibs wrap around domain
"periodic_ibs": "T",
}

case_dict.update(ib_dict)

print(json.dumps(case_dict))
150 changes: 150 additions & 0 deletions examples/3D_ibm_sphere_periodic/case.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,150 @@
import json
import math
import numpy as np

gam_a = 1.4

# sphere diameter
D = 0.1

# sphere density
rho_s = 200

# domain length
L = 10 * D

# mach number
M = 1.2

# reynolds number
Re = 400.0

# pressure
P = 101325

# density
rho = 1.0

# fluid velocity
v1 = M * np.sqrt(gam_a * P / rho)

# dynamic viscosity
mu = rho * v1 * D / Re

dt = 2.0e-06
Nt = 2000
t_save = Nt // 250

# to fully resolve requires 40-60 cells across sphere diameter
Nx = 127
Ny = Nx
Nz = Ny
Comment on lines +38 to +41
Copy link
Contributor

@coderabbitai coderabbitai bot Jan 30, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🟡 Minor

Align the resolution comment with the chosen grid.

With L=1 and Nx=127, the sphere diameter spans ~13 cells, not 40–60. Consider updating the comment to avoid confusion.

💡 Suggested comment tweak 
-# to fully resolve requires 40-60 cells across sphere diameter
+# NOTE: fully resolving needs 40–60 cells across diameter; this example uses ~13 for faster runs
🤖 Prompt for AI Agents 
In `@examples/3D_ibm_sphere_periodic/case.py` around lines 38 - 41, The comment
about required resolution is inconsistent with the grid settings: with L = 1 and
Nx = 127 (Ny = Nz = Nx), the sphere diameter covers roughly 13 cells, not 40–60;
update the comment above the Nx/Ny/Nz definitions (referencing Nx, Ny, Nz and L)
to either state the actual cells-per-diameter for this grid (~13) or note that
40–60 cells would require increasing Nx accordingly, so the comment accurately
reflects the chosen resolution.


# immersed boundary dictionary
ib_dict = {}
ib_dict.update(
{
f"patch_ib({1})%geometry": 8,
f"patch_ib({1})%x_centroid": 0.5 * L,
f"patch_ib({1})%y_centroid": 0.5 * L,
f"patch_ib({1})%z_centroid": 0.5 * L,
f"patch_ib({1})%radius": D / 2,
f"patch_ib({1})%slip": "F",
f"patch_ib({1})%moving_ibm": 0,
f"patch_ib({1})%mass": rho_s * 4.0 / 3.0 * np.pi * (D / 2.0) ** 3,
}
)

# Configuring case dictionary
case_dict = {
# Logistics
"run_time_info": "T",
# Computational Domain Parameters
# x direction
"x_domain%beg": -5.0 * D,
"x_domain%end": 5.0 * D,
# y direction
"y_domain%beg": -5.0 * D,
"y_domain%end": 5.0 * D,
# z direction
"z_domain%beg": -5.0 * D,
"z_domain%end": 5.0 * D,
"cyl_coord": "F",
"m": Nx,
"n": Ny,
"p": Nz,
"dt": dt,
"t_step_start": 0,
"t_step_stop": Nt,
"t_step_save": t_save,
# Simulation Algorithm Parameters
# Only one patch is necessary for one fluid
"num_patches": 1,
# Use the 5 equation model
"model_eqns": 2,
# 6 equations model does not need the K \div(u) term
"alt_soundspeed": "F",
# One fluids: air
"num_fluids": 1,
# time step
"mpp_lim": "F",
# Correct errors when computing speed of sound
"mixture_err": "T",
# Use TVD RK3 for time marching
"time_stepper": 3,
# Reconstruct the primitive variables to minimize spurious
# Use WENO5
"weno_order": 5,
"weno_eps": 1.0e-14,
"weno_Re_flux": "T",
"weno_avg": "T",
"avg_state": 2,
"mapped_weno": "T",
"null_weights": "F",
"mp_weno": "T",
"riemann_solver": 2,
"wave_speeds": 1,
# periodic bc
"bc_x%beg": -1,
"bc_x%end": -1,
"bc_y%beg": -1,
"bc_y%end": -1,
"bc_z%beg": -1,
"bc_z%end": -1,
# Set IB to True and add 1 patch for every sphere
"ib": "T",
"num_ibs": 1,
"viscous": "T",
# Formatted Database Files Structure Parameters
"format": 1,
"precision": 2,
"prim_vars_wrt": "T",
"E_wrt": "T",
"parallel_io": "T",
# Patch: cube filled with air
"patch_icpp(1)%geometry": 9,
# Uniform properties, centroid is at the center of the domain
"patch_icpp(1)%x_centroid": 0.0,
"patch_icpp(1)%y_centroid": 0.0,
"patch_icpp(1)%z_centroid": 0.0,
"patch_icpp(1)%length_x": 10 * D,
"patch_icpp(1)%length_y": 10 * D,
"patch_icpp(1)%length_z": 10 * D,
# Specify the patch primitive variables
"patch_icpp(1)%vel(1)": v1,
"patch_icpp(1)%vel(2)": 0.0e00,
"patch_icpp(1)%vel(3)": 0.0e00,
"patch_icpp(1)%pres": P,
"patch_icpp(1)%alpha_rho(1)": rho,
"patch_icpp(1)%alpha(1)": 1.0e00,
# Fluids Physical Parameters
"fluid_pp(1)%gamma": 1.0e00 / (gam_a - 1.0e00), # 2.50(Not 1.40)
"fluid_pp(1)%pi_inf": 0,
"fluid_pp(1)%Re(1)": 1.0 / mu,
# ibs wrap around domain
"periodic_ibs": "T",
}

case_dict.update(ib_dict)

print(json.dumps(case_dict))
Loading
Loading