From 381ef1ed30b2b6469fab5eb2f301f3dc78ed113a Mon Sep 17 00:00:00 2001 From: Spencer Bryngelson Date: Tue, 5 May 2026 17:44:59 -0400 Subject: [PATCH] Fix 1D flow perturbation applying random factor twice due to index aliasing In 1D, eqn_idx%mom%end == eqn_idx%mom%beg (same slot), so the old code wrote rand*v0 into the slot then read that modified value to compute (1+rand)*rand*v0, doubling the perturbation. Fix: save the original mom%beg value before any writes, scale mom%beg first, and only assign mom%end when num_vels > 1 (i.e., multi-D). Fixes #1369 --- src/pre_process/m_perturbation.fpp | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/pre_process/m_perturbation.fpp b/src/pre_process/m_perturbation.fpp index cc83378396..8472b47430 100644 --- a/src/pre_process/m_perturbation.fpp +++ b/src/pre_process/m_perturbation.fpp @@ -63,7 +63,7 @@ contains type(scalar_field), dimension(sys_size), intent(inout) :: q_prim_vf integer :: i, j, k - real(wp) :: rand_real + real(wp) :: rand_real, v_beg call random_seed() @@ -72,8 +72,9 @@ contains do i = 0, m call random_number(rand_real) rand_real = rand_real*perturb_flow_mag - q_prim_vf(eqn_idx%mom%end)%sf(i, j, k) = rand_real*q_prim_vf(eqn_idx%mom%beg)%sf(i, j, k) - q_prim_vf(eqn_idx%mom%beg)%sf(i, j, k) = (1._wp + rand_real)*q_prim_vf(eqn_idx%mom%beg)%sf(i, j, k) + v_beg = q_prim_vf(eqn_idx%mom%beg)%sf(i, j, k) + q_prim_vf(eqn_idx%mom%beg)%sf(i, j, k) = (1._wp + rand_real)*v_beg + if (num_vels > 1) q_prim_vf(eqn_idx%mom%end)%sf(i, j, k) = rand_real*v_beg if (bubbles_euler) then q_prim_vf(eqn_idx%alf)%sf(i, j, k) = (1._wp + rand_real)*q_prim_vf(eqn_idx%alf)%sf(i, j, k) end if