Skip to content

Large-Scale Startup Optimizations #1686

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
danieljvickers wants to merge 20 commits into
base: master
Choose a base branch
from
Open

Large-Scale Startup Optimizations #1686

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
20 commits
Select commit Hold shift + click to select a range
04cb847
Swapped the particle bed generation to be neighborhood aware during i…
Jul 8, 2026
a22bea9
Fixed an out-of-bounds error caused by reading a bad length for the n…
Jul 9, 2026
a1bec2f
Cleaned up comments
Jul 9, 2026
ba5a7a5
file access stagger and cleanup from timings in full-scale run
Jul 27, 2026
efec23a
Missed a few prints
Jul 27, 2026
3a0d013
Fixing issue with rejection packing algorithm not having enough alloc…
Jul 30, 2026
f35421e
Merge branch 'master' into smarter-particle-bed-instantiation
danieljvickers Jul 30, 2026
2349818
Update to the tuo mako file for tme limits value in the mako
Jul 31, 2026
75c7c61
Addressing some review comments
danieljvickers Jul 31, 2026
92afde1
Addressing more review comments
danieljvickers Aug 4, 2026
7f2f90c
Merge branch 'master' into smarter-particle-bed-instantiation
danieljvickers Aug 4, 2026
1b3d86a
Accidentally modified a subroutine name
danieljvickers Aug 4, 2026
54f1c1b
Swapped iterator and array name
danieljvickers Aug 4, 2026
c5b7968
prohibit requires mpi_proxy module
danieljvickers Aug 4, 2026
52bd1f7
Merge branch 'master' into smarter-particle-bed-instantiation
sbryngelson Aug 9, 2026
48ca630
Merge branch 'master' into smarter-particle-bed-instantiation
sbryngelson Aug 11, 2026
1eaaadd
Merge branch 'master' into smarter-particle-bed-instantiation
danieljvickers Aug 11, 2026
1361dd9
Reordered prohibit checks
danieljvickers Aug 13, 2026
79c914d
Merge branch 'master' into smarter-particle-bed-instantiation
danieljvickers Aug 13, 2026
70ba5bf
Merge branch 'master' into smarter-particle-bed-instantiation
sbryngelson Aug 15, 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
4 changes: 2 additions & 2 deletions src/common/m_boundary_io.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -130,7 +130,7 @@ contains

call s_mpi_barrier()

call DelayFileAccess(proc_rank)
call s_delay_file_access(proc_rank)

write (proc_rank_str, '(I7.7)') proc_rank
file_path = trim(file_loc) // '/bc_' // trim(proc_rank_str) // '.dat'
Expand Down Expand Up @@ -233,7 +233,7 @@ contains

call s_mpi_barrier()

call DelayFileAccess(proc_rank)
call s_delay_file_access(proc_rank)

write (proc_rank_str, '(I7.7)') proc_rank
file_path = trim(file_loc) // '/bc_' // trim(proc_rank_str) // '.dat'
Expand Down
3 changes: 1 addition & 2 deletions src/common/m_constants.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -27,10 +27,9 @@ module m_constants
integer, parameter :: num_ib_airfoils_max = 5 !< Maximum number of ib_airfoil instances
integer, parameter :: num_stl_models_max = 10
!> Maximum number of immersed boundary patches (legacy, not used for patch_ib sizing)
integer, parameter :: num_ib_patches_max = 2050000
!> Fixed capacity of patch_ib (namelist patches + local particle bed subset after reduction)
integer, parameter :: num_ib_patches_max_namelist = 54000
integer, parameter :: num_local_ibs_max = 2000 !< Maximum number of immersed boundary patches (patch_ib)
integer, parameter :: num_ib_patches_max_namelist = 54000
integer, parameter :: num_particle_clouds_max = 10 !< Maximum number of particle bed patch specifications
integer, parameter :: num_bc_patches_max = 10 !< Maximum number of boundary condition patches
integer, parameter :: max_2d_fourier_modes = 10 !< Max Fourier mode index for 2D modal patch (geometry 13)
Expand Down
20 changes: 10 additions & 10 deletions src/common/m_delay_file_access.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,26 +11,26 @@ module m_delay_file_access

private

public :: DelayFileAccess
public :: s_delay_file_access
Comment thread
danieljvickers marked this conversation as resolved.

integer, private, parameter :: N_PROCESSES_FILE_ACCESS = 128, FILE_ACCESS_DELAY_UNIT = 10000

contains

!> Introduce a rank-dependent busy-wait delay to stagger parallel file access and reduce I/O contention.
impure subroutine DelayFileAccess(ProcessRank)
impure subroutine s_delay_file_access(process_rank)

integer, intent(in) :: ProcessRank
integer :: iDelay, nFileAccessDelayIterations
real(wp) :: Number, Dummy
integer, intent(in) :: process_rank
integer :: i, n_file_access_delay_iterations
real(wp) :: num, dummy

nFileAccessDelayIterations = (ProcessRank/N_PROCESSES_FILE_ACCESS)*FILE_ACCESS_DELAY_UNIT
n_file_access_delay_iterations = (process_rank/N_PROCESSES_FILE_ACCESS)*FILE_ACCESS_DELAY_UNIT

do iDelay = 1, nFileAccessDelayIterations
call random_number(Number)
Dummy = Number*Number
do i = 1, n_file_access_delay_iterations
call random_number(num)
dummy = num*num
end do

end subroutine DelayFileAccess
end subroutine s_delay_file_access

end module m_delay_file_access
2 changes: 1 addition & 1 deletion src/pre_process/m_data_output.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -454,7 +454,7 @@ contains
call s_create_directory(trim(file_loc))
end if
call s_mpi_barrier()
call DelayFileAccess(proc_rank)
call s_delay_file_access(proc_rank)

if (down_sample) then
call s_initialize_mpi_data_ds(m_ds, n_ds, p_ds)
Expand Down
62 changes: 8 additions & 54 deletions src/simulation/m_collisions.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -16,6 +16,7 @@ module m_collisions
use m_compute_levelset
use m_ib_patches
use m_model
use m_mpi_proxy

implicit none

Expand Down Expand Up @@ -65,8 +66,7 @@ contains

! get is distance used in the force calculation with each IB and each wall
call s_detect_wall_collisions()
! call s_detect_ib_collisions(ghost_points, ib_markers, num_gps, num_considered_collisions)
call s_detect_ib_collisions_n2(num_considered_collisions)
call s_detect_ib_collisions(ghost_points, ib_markers, num_gps, num_considered_collisions)

select case (collision_model)
case (1) ! soft sphere model
Expand Down Expand Up @@ -253,7 +253,7 @@ contains
integer :: num_raw, local_num_raw

num_raw = 0
z_bound = 0; if (num_dims == 3) z_bound = 1
z_bound = 0; if (num_dims == 3) z_bound = 2

$:GPU_PARALLEL_LOOP(private='[gp_idx, gp_patch_id, neighbor_patch_id, local_num_raw, i, j, k, ii, jj, kk]', &
& copy='[raw_pairs, num_raw]', copyin='[z_bound]')
Expand All @@ -264,8 +264,8 @@ contains
gp_patch_id = ib_markers%sf(i, j, k)

! search in a cube around the BG for Ib markers belonging to another patch
neighbor_search: do ii = i - 1, i + 1
do jj = j - 1, j + 1
neighbor_search: do ii = i - 2, i + 2
do jj = j - 2, j + 2
do kk = k - z_bound, k + z_bound
neighbor_patch_id = ib_markers%sf(ii, jj, kk)

Expand Down Expand Up @@ -324,6 +324,9 @@ contains
! and if it is not, append it to the list of pairs
if (.not. already_found) then
num_considered_collisions = num_considered_collisions + 1
@:PROHIBIT(num_considered_collisions > size(collision_lookup, 1) , &
& "More collisions detected than memory to hold them. Consider increasing the size of the collision_lookup array")

collision_lookup(num_considered_collisions, 1) = decoded_pairs(1)
collision_lookup(num_considered_collisions, 2) = decoded_pairs(2)
collision_lookup(num_considered_collisions, 3) = raw_pairs(pair_idx, 1)
Expand All @@ -334,55 +337,6 @@ contains

end subroutine s_detect_ib_collisions

subroutine s_detect_ib_collisions_n2(num_considered_collisions)

integer, intent(out) :: num_considered_collisions
integer :: pid1, pid2, encoded_pid2, current_collisions
integer :: xp_lower, xp_upper, yp_lower, yp_upper, zp_lower, zp_upper, xp, yp, zp
real(wp), dimension(3) :: centroid_1, centroid_2, distance_vec

num_considered_collisions = 0

call s_get_periodicities(xp_lower, xp_upper, yp_lower, yp_upper, zp_lower, zp_upper)

$:GPU_PARALLEL_LOOP(private='[pid1, pid2, encoded_pid2, centroid_1, centroid_2, xp, yp, zp, distance_vec, &
& current_collisions]', copyin='[xp_lower, xp_upper, yp_lower, yp_upper, zp_lower, zp_upper]', copy='[num_considered_collisions]')
do pid1 = 1, num_ibs - 1
centroid_1 = [patch_ib(pid1)%x_centroid, patch_ib(pid1)%y_centroid, 0._wp]
if (num_dims == 3) centroid_1(3) = patch_ib(pid1)%z_centroid
do pid2 = pid1 + 1, num_ibs
periodic_search: do xp = xp_lower, xp_upper
do yp = yp_lower, yp_upper
do zp = zp_lower, zp_upper
centroid_2(1) = patch_ib(pid2)%x_centroid + real(xp, wp)*(glb_bounds(1)%end - glb_bounds(1)%beg)
centroid_2(2) = patch_ib(pid2)%y_centroid + real(yp, wp)*(glb_bounds(2)%end - glb_bounds(2)%beg)
if (num_dims == 3) centroid_2(3) = patch_ib(pid2)%z_centroid + real(zp, &
& wp)*(glb_bounds(3)%end - glb_bounds(3)%beg)
distance_vec = centroid_2 - centroid_1

if (norm2(distance_vec) < patch_ib(pid1)%radius + patch_ib(pid2)%radius) then
$:GPU_ATOMIC(atomic='capture')
num_considered_collisions = num_considered_collisions + 1
current_collisions = num_considered_collisions
$:END_GPU_ATOMIC_CAPTURE()

call s_encode_patch_periodicity(patch_ib(pid2)%gbl_patch_id, xp, yp, zp, encoded_pid2)

collision_lookup(current_collisions, 1) = pid1
collision_lookup(current_collisions, 2) = pid2
collision_lookup(current_collisions, 3) = patch_ib(pid1)%gbl_patch_id
collision_lookup(current_collisions, 4) = encoded_pid2
exit periodic_search
end if
end do
end do
end do periodic_search
end do
end do
$:END_GPU_PARALLEL_LOOP()

end subroutine s_detect_ib_collisions_n2

!> @brief uses boundary conditions and particle locations to check for wall conditions
subroutine s_detect_wall_collisions()

Expand Down
8 changes: 6 additions & 2 deletions src/simulation/m_data_output.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -717,7 +717,7 @@ contains
call s_create_directory(trim(file_loc))
end if
call s_mpi_barrier()
call DelayFileAccess(proc_rank)
call s_delay_file_access(proc_rank)

call s_initialize_mpi_data(q_cons_vf, qbmm_pb=pb_ts(1), qbmm_mv=mv_ts(1))

Expand Down Expand Up @@ -902,6 +902,10 @@ contains

write (file_loc, '(A)') 'ib.dat'
file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // trim(file_loc)

call s_mpi_barrier()
call s_delay_file_access(proc_rank)

@sbryngelson sbryngelson Jul 31, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

Does this actually help?

MPI_FILE_OPEN(MPI_COMM_WORLD, ...) is collective -- staggering rank arrival cannot reduce contention on a single shared-file open; every rank just ends up waiting for the largest delay (rank 78000 => ~6.1M random_number calls). Every other s_delay_file_access call site precedes a per-rank file open, where the staggering does do something.

If you measured a real win here at full-system scale, please put the numbers in the PR -- otherwise I would drop the barrier and the delay.

@sbryngelson sbryngelson Aug 1, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

Still open, and I want to be precise about why this one differs from the other call sites, because the rename made them all look alike.

Every pre-existing s_delay_file_access call precedes a per-rank open:

  • m_data_output.fpp:716 -> MPI_FILE_OPEN(MPI_COMM_SELF, ...) on <t_step>_<rank>.dat
  • m_data_output.fpp:964 -> per-rank ib_state_<t_step>_<rank>.dat
  • m_boundary_io.fpp:133/:236 -> per-rank bc_<rank>.dat
  • pre_process/m_data_output.fpp:457 -> same pattern

There, staggering spreads 78k independent metadata operations across time, which is the thing that helps.

The new site is before MPI_FILE_OPEN(MPI_COMM_WORLD, file_loc, ...) on a single shared ib.dat (:905). That call is collective: no rank returns from it until all ranks have entered, so staggering arrival cannot spread anything -- it just idles the early ranks. The added wall time per write is the maximum delay over all ranks, (num_procs/128)*10000 random_number calls, which at 78k ranks is ~6.1M iterations on the last rank, paid on every ib.dat write. The s_mpi_barrier() on the preceding line makes it strictly worse: it forces every rank to wait for the slowest before any of them starts its delay.

If you measured a win at full scale I will take the numbers over my reasoning -- MPI-IO on Lustre does surprising things, and it is possible you were working around something in the file-creation path specifically. But as written I would expect this to cost ~seconds per output step at 78k ranks and buy nothing, so please either post the timings or drop both lines.


call MPI_FILE_OPEN(MPI_COMM_WORLD, file_loc, ior(MPI_MODE_WRONLY, MPI_MODE_CREATE), mpi_info_int, ifile, ierr)

var_MOK = int(sys_size + 1, MPI_OFFSET_KIND)
Expand Down Expand Up @@ -961,7 +965,7 @@ contains
call s_create_directory(trim(file_loc))
end if
call s_mpi_barrier()
call DelayFileAccess(proc_rank)
call s_delay_file_access(proc_rank)

write (file_loc, '(A,I0,A,i7.7,A)') 'ib_state_', t_step, '_', proc_rank, '.dat'
file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string) // '/' // trim(file_loc)
Expand Down
Loading