Skip to content

Fix 6 low-risk pre-process bugs (batch) #1241

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

Draft
sbryngelson wants to merge 17 commits into
base: master
Choose a base branch
from
+26 −18
Draft

Fix 6 low-risk pre-process bugs (batch) #1241

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
17 commits
Select commit Hold shift + click to select a range
32cb2da
Fix moncon_cutoff declared as integer, truncating 1e-8 to 0
sbryngelson Feb 21, 2026
c34851c
Fix grid stretching using wrong array bounds for y_cc and z_cc
sbryngelson Feb 21, 2026
f7ab2da
Remove duplicate R3bar accumulation line that doubles bubble perturba…
sbryngelson Feb 21, 2026
630bb63
Fix y-velocity perturbation using already-modified x-velocity
sbryngelson Feb 21, 2026
e4552bc
Fix hyper_cleaning psi initialization only covering k=0 plane
sbryngelson Feb 21, 2026
38f3be9
Guard z_cc access in hyper_cleaning init for 2D case (p=0)
sbryngelson Feb 21, 2026
b985e53
Fix mixed-precision literals in hyper_cleaning psi initialization
sbryngelson Feb 21, 2026
2ae93a2
Fix IB patch vel/angular_vel/angles broadcast inside wrong loop
sbryngelson Feb 21, 2026
1c1d8ac
Simplify hyper_cleaning init and use size() for IB broadcast
sbryngelson Feb 23, 2026
61f39b1
Guard y_cc and z_cc access in hyper_cleaning init for 1D/2D safety
sbryngelson Feb 23, 2026
3b87205
Add ASSERT for psi_idx and clarifying comment for IB loop scope
sbryngelson Feb 24, 2026
7cb0512
Prohibit hyper_cleaning in 2D simulations (p == 0)
sbryngelson Feb 24, 2026
44986bf
Fix hyper_cleaning 2D prohibition double-triggering for 1D cases
sbryngelson Feb 24, 2026
0ed89d3
Remove incorrect 2D hyper_cleaning prohibition
sbryngelson Feb 25, 2026
3d15838
Fix typos in m_assign_variables.fpp and m_perturbation.fpp
sbryngelson Feb 26, 2026
732146c
Add missing macros.fpp include for ASSERT macro in m_start_up.fpp
sbryngelson Feb 26, 2026
1cdd751
Merge branch 'master' into fix/low-risk-bugfixes-batch
sbryngelson Mar 2, 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
2 changes: 1 addition & 1 deletion src/common/m_constants.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,7 @@ module m_constants
! Interface Compression
real(wp), parameter :: dflt_ic_eps = 1e-4_wp !< Ensure compression is only applied to surface cells in THINC
real(wp), parameter :: dflt_ic_beta = 1.6_wp !< Sharpness parameter's default value used in THINC
integer, parameter :: moncon_cutoff = 1e-8_wp !< Monotonicity constraint's limiter to prevent extremas in THINC
real(wp), parameter :: moncon_cutoff = 1e-8_wp !< Monotonicity constraint's limiter to prevent extremas in THINC

! Chemistry
real(wp), parameter :: dflt_T_guess = 1200._wp ! Default guess for temperature (when a previous value is not available)
Expand Down
7 changes: 3 additions & 4 deletions src/pre_process/m_assign_variables.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -38,9 +38,9 @@ module m_assign_variables
!> Skeleton of s_assign_patch_mixture_primitive_variables
!! and s_assign_patch_species_primitive_variables
!! @param patch_id is the patch identifier
!! @param j (x) cell index in which the mixture or species primitive variables from the indicated patch areassigned
!! @param k (y,th) cell index in which the mixture or species primitive variables from the indicated patch areassigned
!! @param l (z) cell index in which the mixture or species primitive variables from the indicated patch areassigned
!! @param j (x) cell index in which the mixture or species primitive variables from the indicated patch are assigned
!! @param k (y,th) cell index in which the mixture or species primitive variables from the indicated patch are assigned
!! @param l (z) cell index in which the mixture or species primitive variables from the indicated patch are assigned
!! @param eta pseudo volume fraction
!! @param q_prim_vf Primitive variables
!! @param patch_id_fp Array to track patch ids
Expand Down Expand Up @@ -233,7 +233,6 @@ contains
if (qbmm) then
do i = 1, nb
R3bar = R3bar + weight(i)*0.5_wp*(q_prim_vf(bubxb + 1 + (i - 1)*nmom)%sf(j, k, l))**3._wp
R3bar = R3bar + weight(i)*0.5_wp*(q_prim_vf(bubxb + 1 + (i - 1)*nmom)%sf(j, k, l))**3._wp
end do
else
do i = 1, nb
Expand Down
4 changes: 2 additions & 2 deletions src/pre_process/m_grid.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -131,7 +131,7 @@ impure subroutine s_generate_serial_grid
end do

y_cb = y_cb*length
y_cc(0:m) = (y_cb(0:n) + y_cb(-1:n - 1))/2._wp
y_cc(0:n) = (y_cb(0:n) + y_cb(-1:n - 1))/2._wp

dy = minval(y_cb(0:n) - y_cb(-1:n - 1))

Expand Down Expand Up @@ -168,7 +168,7 @@ impure subroutine s_generate_serial_grid
end do

z_cb = z_cb*length
z_cc(0:m) = (z_cb(0:p) + z_cb(-1:p - 1))/2._wp
z_cc(0:p) = (z_cb(0:p) + z_cb(-1:p - 1))/2._wp

dz = minval(z_cb(0:p) - z_cb(-1:p - 1))

Expand Down
10 changes: 5 additions & 5 deletions src/pre_process/m_mpi_proxy.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -76,10 +76,6 @@ contains
call MPI_BCAST(patch_bc(i)%${VAR}$, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
#:endfor

#:for VAR in ['vel', 'angular_vel', 'angles']
call MPI_BCAST(patch_ib(i)%${VAR}$, 3, mpi_p, 0, MPI_COMM_WORLD, ierr)
#:endfor

call MPI_BCAST(patch_bc(i)%radius, 1, mpi_p, 0, MPI_COMM_WORLD, ierr)

#:for VAR in ['centroid', 'length']
Expand Down Expand Up @@ -120,7 +116,11 @@ contains
if (chemistry) then
call MPI_BCAST(patch_icpp(i)%Y, size(patch_icpp(i)%Y), mpi_p, 0, MPI_COMM_WORLD, ierr)
end if
! Broadcast IB variables
! Broadcast IB variables: patch_ib is indexed 1:num_patches_max,
! not 1:num_bc_patches_max, so these must live in the num_patches_max loop.
#:for VAR in ['vel', 'angular_vel', 'angles']
Copy link
Member

@danieljvickers danieljvickers Feb 24, 2026

Choose a reason for hiding this comment

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

This change is pointless. It just copies the values one at a time instead copying 3 at once. In fact, it is probably slightly slower.

Copy link
Contributor

@wilfonba wilfonba Feb 25, 2026

Choose a reason for hiding this comment

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

I'm pretty sure this still broadcasts all three elements of each variable at the same time unless I'm missing something

call MPI_BCAST(patch_ib(i)%${VAR}$, size(patch_ib(i)%${VAR}$), mpi_p, 0, MPI_COMM_WORLD, ierr)
#:endfor
call MPI_BCAST(patch_ib(i)%geometry, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
call MPI_BCAST(patch_ib(i)%model_filepath, len(patch_ib(i)%model_filepath), MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr)
call MPI_BCAST(patch_ib(i)%model_threshold, 1, mpi_p, 0, MPI_COMM_WORLD, ierr)
Expand Down
4 changes: 2 additions & 2 deletions src/pre_process/m_perturbation.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -83,8 +83,8 @@ contains
perturb_alpha = q_prim_vf(E_idx + perturb_flow_fluid)%sf(i, j, k)
call random_number(rand_real)
rand_real = rand_real*perturb_flow_mag
q_prim_vf(mom_idx%beg)%sf(i, j, k) = (1._wp + rand_real)*q_prim_vf(mom_idx%beg)%sf(i, j, k)
q_prim_vf(mom_idx%end)%sf(i, j, k) = rand_real*q_prim_vf(mom_idx%beg)%sf(i, j, k)
q_prim_vf(mom_idx%beg)%sf(i, j, k) = (1._wp + rand_real)*q_prim_vf(mom_idx%beg)%sf(i, j, k)
if (bubbles_euler) then
q_prim_vf(alf_idx)%sf(i, j, k) = (1._wp + rand_real)*q_prim_vf(alf_idx)%sf(i, j, k)
end if
Expand Down Expand Up @@ -245,7 +245,7 @@ contains

!> This subroutine computes velocity perturbations for a temporal mixing
!! layer with a hyperbolic tangent mean streamwise velocity
!! profile, using an inverter version of the spectrum-based
!! profile, using an inverted version of the spectrum-based
!! synthetic turbulence generation method proposed by
!! Guo et al. (2023, JFM).
subroutine s_perturb_mixlayer(q_prim_vf)
Expand Down
17 changes: 13 additions & 4 deletions src/pre_process/m_start_up.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,8 @@
!! @file
!! @brief Contains module m_start_up

#:include 'macros.fpp'

!> @brief Reads and validates user inputs, loads existing grid/IC data, and initializes pre-process modules
module m_start_up

Expand Down Expand Up @@ -767,7 +769,8 @@ contains

real(wp), intent(inout) :: start, finish

integer :: j, k
integer :: j, k, l
real(wp) :: r2

! Setting up the grid and the initial condition. If the grid is read in from
! preexisting grid data files, it is checked for consistency. If the grid is
Expand All @@ -787,10 +790,16 @@ contains

! hard-coded psi
if (hyper_cleaning) then
do j = 0, m
@:ASSERT(psi_idx > 0, "hyper_cleaning requires psi_idx to be set")
do l = 0, p
do k = 0, n
q_cons_vf(psi_idx)%sf(j, k, 0) = 1d-2*exp(-(x_cc(j)**2 + y_cc(k)**2)/(2.0*0.05**2))
q_prim_vf(psi_idx)%sf(j, k, 0) = q_cons_vf(psi_idx)%sf(j, k, 0)
do j = 0, m
r2 = x_cc(j)**2
if (n > 0) r2 = r2 + y_cc(k)**2
if (p > 0) r2 = r2 + z_cc(l)**2
q_cons_vf(psi_idx)%sf(j, k, l) = 1.0e-2_wp*exp(-r2/(2.0_wp*0.05_wp**2))
q_prim_vf(psi_idx)%sf(j, k, l) = q_cons_vf(psi_idx)%sf(j, k, l)
Copy link
Member

@danieljvickers danieljvickers Feb 24, 2026

Choose a reason for hiding this comment

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

This is a good catch, i think. May be correct to ask @ChrisZYJ since I though the point of giving psi initial values was just to add some stability at the start, so I am surprised we went this long without finding this.

Copy link
Contributor

@ChrisZYJ ChrisZYJ Feb 25, 2026

Choose a reason for hiding this comment

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

Yeah this fix makes sense. It was my oversight for only initializing the 2D plane as I was initially only testing 2D cases. By construction, psi encodes the div B error so it propagates and dampens, and an initial mismatch just goes away, which is probably why it went unidentified. Now it's consistent and more stable.

end do
end do
end do
end if
Expand Down
Loading