From 4c65007db5785ec04ff36488b5aa3cd13e8cec22 Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Fri, 20 Feb 2026 22:17:01 -0500 Subject: [PATCH 01/10] Fix z-boundary backward FD stencil using fields(2) instead of fields(3) The z-direction upper-boundary backward difference at iz_s%end uses fields(2) (y-component) instead of fields(3) (z-component) in the third term, corrupting the divergence in all 3D simulations using this finite difference routine. Co-Authored-By: Claude Opus 4.6 --- src/common/m_finite_differences.fpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/common/m_finite_differences.fpp b/src/common/m_finite_differences.fpp index 09dba0a903..e44b6905c0 100644 --- a/src/common/m_finite_differences.fpp +++ b/src/common/m_finite_differences.fpp @@ -50,7 +50,7 @@ contains if (z == iz_s%beg) then divergence = divergence + (-3._wp*fields(3)%sf(x, y, z) + 4._wp*fields(3)%sf(x, y, z + 1) - fields(3)%sf(x, y, z + 2))/(z_cc(z + 2) - z_cc(z)) else if (z == iz_s%end) then - divergence = divergence + (+3._wp*fields(3)%sf(x, y, z) - 4._wp*fields(3)%sf(x, y, z - 1) + fields(2)%sf(x, y, z - 2))/(z_cc(z) - z_cc(z - 2)) + divergence = divergence + (+3._wp*fields(3)%sf(x, y, z) - 4._wp*fields(3)%sf(x, y, z - 1) + fields(3)%sf(x, y, z - 2))/(z_cc(z) - z_cc(z - 2)) else divergence = divergence + (fields(3)%sf(x, y, z + 1) - fields(3)%sf(x, y, z - 1))/(z_cc(z + 1) - z_cc(z - 1)) end if From 84577113f45f2ad3594bb69267e36c2b1406dcf9 Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Fri, 20 Feb 2026 22:17:33 -0500 Subject: [PATCH 02/10] Fix bc_x%ve3 never MPI-broadcast (duplicate ve2 instead) The bc_x velocity-end broadcast list has bc_x%ve2 duplicated where bc_x%ve3 should be. bc_y and bc_z rows are correct. Non-root ranks get uninitialized bc_x%ve3 in multi-rank 3D runs with x-boundary velocity BCs. Co-Authored-By: Claude Opus 4.6 --- src/simulation/m_mpi_proxy.fpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/simulation/m_mpi_proxy.fpp b/src/simulation/m_mpi_proxy.fpp index be737782be..95f61fc7d7 100644 --- a/src/simulation/m_mpi_proxy.fpp +++ b/src/simulation/m_mpi_proxy.fpp @@ -145,7 +145,7 @@ contains #:for VAR in [ 'dt','weno_eps','teno_CT','pref','rhoref','R0ref','Web','Ca', 'sigma', & & 'Re_inv', 'poly_sigma', 'palpha_eps', 'ptgalpha_eps', 'pi_fac', & - & 'bc_x%vb1','bc_x%vb2','bc_x%vb3','bc_x%ve1','bc_x%ve2','bc_x%ve2', & + & 'bc_x%vb1','bc_x%vb2','bc_x%vb3','bc_x%ve1','bc_x%ve2','bc_x%ve3', & & 'bc_y%vb1','bc_y%vb2','bc_y%vb3','bc_y%ve1','bc_y%ve2','bc_y%ve3', & & 'bc_z%vb1','bc_z%vb2','bc_z%vb3','bc_z%ve1','bc_z%ve2','bc_z%ve3', & & 'bc_x%pres_in','bc_x%pres_out','bc_y%pres_in','bc_y%pres_out', 'bc_z%pres_in','bc_z%pres_out', & From 25e4ac6c8be25035c49310bf0ebc7ea824646a83 Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Fri, 20 Feb 2026 22:18:19 -0500 Subject: [PATCH 03/10] Fix fluid_rho broadcast using MPI_LOGICAL instead of mpi_p fluid_rho is a real(wp) array but is broadcast with MPI_LOGICAL type, silently corrupting reference densities via bit reinterpretation on non-root ranks. Co-Authored-By: Claude Opus 4.6 --- src/pre_process/m_mpi_proxy.fpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pre_process/m_mpi_proxy.fpp b/src/pre_process/m_mpi_proxy.fpp index 30ef061689..cbfac0571b 100644 --- a/src/pre_process/m_mpi_proxy.fpp +++ b/src/pre_process/m_mpi_proxy.fpp @@ -58,7 +58,7 @@ contains & 'igr', 'down_sample', 'simplex_perturb','fft_wrt', 'hyper_cleaning' ] call MPI_BCAST(${VAR}$, 1, MPI_LOGICAL, 0, MPI_COMM_WORLD, ierr) #:endfor - call MPI_BCAST(fluid_rho(1), num_fluids_max, MPI_LOGICAL, 0, MPI_COMM_WORLD, ierr) + call MPI_BCAST(fluid_rho(1), num_fluids_max, mpi_p, 0, MPI_COMM_WORLD, ierr) #:for VAR in [ 'x_domain%beg', 'x_domain%end', 'y_domain%beg', & & 'y_domain%end', 'z_domain%beg', 'z_domain%end', 'a_x', 'a_y', & From 988ee4d55aeb9c5fc7b4f8932aea394a3c545cf4 Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Fri, 20 Feb 2026 22:20:18 -0500 Subject: [PATCH 04/10] Fix loc_violations used uninitialized in MPI_Allreduce loc_violations is never set to 0 before the conditional that may or may not assign it. Non-violating ranks sum garbage in the reduction. Co-Authored-By: Claude Opus 4.6 --- src/pre_process/m_data_output.fpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pre_process/m_data_output.fpp b/src/pre_process/m_data_output.fpp index bac4fdc038..e057234fbd 100644 --- a/src/pre_process/m_data_output.fpp +++ b/src/pre_process/m_data_output.fpp @@ -476,7 +476,7 @@ contains ! Generic loop iterators integer :: i, j, k, l - real(wp) :: loc_violations, glb_violations + real(wp) :: loc_violations = 0._wp, glb_violations ! Downsample variables integer :: m_ds, n_ds, p_ds From 83d261864e06529e3d0b9bcab43c93f988dff9a8 Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Fri, 20 Feb 2026 23:48:21 -0500 Subject: [PATCH 05/10] Fix loc_violations gaining implicit SAVE attribute Initializing a local variable in its declaration gives it the SAVE attribute in Fortran, meaning it would not reset to zero on subsequent calls. Move the initialization to an executable assignment so the variable is properly zeroed each time the subroutine is entered. Co-Authored-By: Claude Sonnet 4.6 --- src/pre_process/m_data_output.fpp | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/pre_process/m_data_output.fpp b/src/pre_process/m_data_output.fpp index e057234fbd..1311a105b7 100644 --- a/src/pre_process/m_data_output.fpp +++ b/src/pre_process/m_data_output.fpp @@ -476,13 +476,15 @@ contains ! Generic loop iterators integer :: i, j, k, l - real(wp) :: loc_violations = 0._wp, glb_violations + real(wp) :: loc_violations, glb_violations ! Downsample variables integer :: m_ds, n_ds, p_ds integer :: m_glb_ds, n_glb_ds, p_glb_ds integer :: m_glb_save, n_glb_save, p_glb_save ! Size of array being saved + loc_violations = 0._wp + if (down_sample) then if ((mod(m + 1, 3) > 0) .or. (mod(n + 1, 3) > 0) .or. (mod(p + 1, 3) > 0)) then loc_violations = 1._wp From 52e10102c40b47b8955ae9a184fa68e1cd92bb73 Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Sat, 21 Feb 2026 17:56:10 -0500 Subject: [PATCH 06/10] Fix GRCBC subsonic inflow using wrong L index In the GRCBC subsonic inflow loop (do i = 2, momxb), L(2) was hardcoded instead of L(i), causing only the second wave amplitude to be updated rather than each wave amplitude in the loop. Co-Authored-By: Claude Opus 4.6 --- src/simulation/m_cbc.fpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/simulation/m_cbc.fpp b/src/simulation/m_cbc.fpp index f9d3f3a929..f5f9cb4c30 100644 --- a/src/simulation/m_cbc.fpp +++ b/src/simulation/m_cbc.fpp @@ -930,7 +930,7 @@ contains if (bc_${XYZ}$%grcbc_in) then $:GPU_LOOP(parallelism='[seq]') do i = 2, momxb - L(2) = c**3._wp*Ma*(alpha_rho(i - 1) - alpha_rho_in(i - 1, ${CBC_DIR}$))/Del_in(${CBC_DIR}$) - c*Ma*(pres - pres_in(${CBC_DIR}$))/Del_in(${CBC_DIR}$) + L(i) = c**3._wp*Ma*(alpha_rho(i - 1) - alpha_rho_in(i - 1, ${CBC_DIR}$))/Del_in(${CBC_DIR}$) - c*Ma*(pres - pres_in(${CBC_DIR}$))/Del_in(${CBC_DIR}$) end do if (n > 0) then L(momxb + 1) = c*Ma*(vel(dir_idx(2)) - vel_in(${CBC_DIR}$, dir_idx(2)))/Del_in(${CBC_DIR}$) From 1b91f13467327a8475766fa5aaa5460ce7713259 Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Sat, 21 Feb 2026 17:57:11 -0500 Subject: [PATCH 07/10] Fix G_K exponential degradation in damage model The damage factor was applied inside the stress component loop, causing G_K (and G) to be multiplied by the damage factor on every iteration. With N stress components, the effective shear modulus was reduced by damage^N instead of damage^1. Move the damage application before the loop so it is applied exactly once per cell. Co-Authored-By: Claude Opus 4.6 --- src/common/m_variables_conversion.fpp | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/src/common/m_variables_conversion.fpp b/src/common/m_variables_conversion.fpp index 4f4aca7f05..9c4e72258f 100644 --- a/src/common/m_variables_conversion.fpp +++ b/src/common/m_variables_conversion.fpp @@ -852,11 +852,11 @@ contains end if if (hypoelasticity) then + if (cont_damage) G_K = G_K*max((1._wp - qK_cons_vf(damage_idx)%sf(j, k, l)), 0._wp) $:GPU_LOOP(parallelism='[seq]') do i = strxb, strxe ! subtracting elastic contribution for pressure calculation if (G_K > verysmall) then - if (cont_damage) G_K = G_K*max((1._wp - qK_cons_vf(damage_idx)%sf(j, k, l)), 0._wp) qK_prim_vf(E_idx)%sf(j, k, l) = qK_prim_vf(E_idx)%sf(j, k, l) - & ((qK_prim_vf(i)%sf(j, k, l)**2._wp)/(4._wp*G_K))/gamma_K ! Double for shear stresses @@ -1123,11 +1123,10 @@ contains end if if (hypoelasticity) then + if (cont_damage) G = G*max((1._wp - q_prim_vf(damage_idx)%sf(j, k, l)), 0._wp) do i = strxb, strxe ! adding elastic contribution if (G > verysmall) then - if (cont_damage) G = G*max((1._wp - q_prim_vf(damage_idx)%sf(j, k, l)), 0._wp) - q_cons_vf(E_idx)%sf(j, k, l) = q_cons_vf(E_idx)%sf(j, k, l) + & (q_prim_vf(i)%sf(j, k, l)**2._wp)/(4._wp*G) ! Double for shear stresses From 445b5b7157e733159c7d1f89562919477f76c065 Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Sat, 21 Feb 2026 17:57:29 -0500 Subject: [PATCH 08/10] Fix domain decomposition overwriting muscl_order The IGR conditional block unconditionally reset recon_order to weno_order in its else branch, overwriting the muscl_order that was correctly set by the recon_type check above. Remove the else branch so the original recon_order is preserved when IGR is inactive. Co-Authored-By: Claude Opus 4.6 --- src/common/m_mpi_common.fpp | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/common/m_mpi_common.fpp b/src/common/m_mpi_common.fpp index de3ed8150a..197b5b912c 100644 --- a/src/common/m_mpi_common.fpp +++ b/src/common/m_mpi_common.fpp @@ -1153,8 +1153,6 @@ contains if (igr) then recon_order = igr_order - else - recon_order = weno_order end if ! 3D Cartesian Processor Topology From f904239ccdc3b32787e7a497c54f29147915bcac Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Sat, 21 Feb 2026 17:57:51 -0500 Subject: [PATCH 09/10] Fix NaN check reading wrong index in MPI unpack The NaN diagnostic check used q_comm(i)%sf(j, k, l) but the value was unpacked into q_comm(i)%sf(j + unpack_offset, k, l). This meant the check was reading a stale or unrelated cell instead of the just- received value. Co-Authored-By: Claude Opus 4.6 --- src/common/m_mpi_common.fpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/common/m_mpi_common.fpp b/src/common/m_mpi_common.fpp index 197b5b912c..000203386a 100644 --- a/src/common/m_mpi_common.fpp +++ b/src/common/m_mpi_common.fpp @@ -936,7 +936,7 @@ contains (j + buff_size*((k + 1) + (n + 1)*l)) q_comm(i)%sf(j + unpack_offset, k, l) = real(buff_recv(r), kind=stp) #if defined(__INTEL_COMPILER) - if (ieee_is_nan(q_comm(i)%sf(j, k, l))) then + if (ieee_is_nan(q_comm(i)%sf(j + unpack_offset, k, l))) then print *, "Error", j, k, l, i error stop "NaN(s) in recv" end if @@ -991,7 +991,7 @@ contains ((k + buff_size) + buff_size*l)) q_comm(i)%sf(j, k + unpack_offset, l) = real(buff_recv(r), kind=stp) #if defined(__INTEL_COMPILER) - if (ieee_is_nan(q_comm(i)%sf(j, k, l))) then + if (ieee_is_nan(q_comm(i)%sf(j, k + unpack_offset, l))) then print *, "Error", j, k, l, i error stop "NaN(s) in recv" end if @@ -1050,7 +1050,7 @@ contains (l + buff_size))) q_comm(i)%sf(j, k, l + unpack_offset) = real(buff_recv(r), kind=stp) #if defined(__INTEL_COMPILER) - if (ieee_is_nan(q_comm(i)%sf(j, k, l))) then + if (ieee_is_nan(q_comm(i)%sf(j, k, l + unpack_offset))) then print *, "Error", j, k, l, i error stop "NaN(s) in recv" end if From 0289d4f1e966babfd5d68591278abf25f1c45b4a Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Mon, 23 Feb 2026 21:33:19 -0500 Subject: [PATCH 10/10] Replace error stop with s_mpi_abort in NaN diagnostic CLAUDE.md Critical Rules forbid error stop. Use s_mpi_abort for proper MPI-aware termination in the Intel NaN detection blocks. Co-Authored-By: Claude Opus 4.6 --- src/common/m_mpi_common.fpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/common/m_mpi_common.fpp b/src/common/m_mpi_common.fpp index 000203386a..8768298306 100644 --- a/src/common/m_mpi_common.fpp +++ b/src/common/m_mpi_common.fpp @@ -938,7 +938,7 @@ contains #if defined(__INTEL_COMPILER) if (ieee_is_nan(q_comm(i)%sf(j + unpack_offset, k, l))) then print *, "Error", j, k, l, i - error stop "NaN(s) in recv" + call s_mpi_abort("NaN(s) in recv") end if #endif end do @@ -993,7 +993,7 @@ contains #if defined(__INTEL_COMPILER) if (ieee_is_nan(q_comm(i)%sf(j, k + unpack_offset, l))) then print *, "Error", j, k, l, i - error stop "NaN(s) in recv" + call s_mpi_abort("NaN(s) in recv") end if #endif end do @@ -1052,7 +1052,7 @@ contains #if defined(__INTEL_COMPILER) if (ieee_is_nan(q_comm(i)%sf(j, k, l + unpack_offset))) then print *, "Error", j, k, l, i - error stop "NaN(s) in recv" + call s_mpi_abort("NaN(s) in recv") end if #endif end do