SyncRegister::InitRHS: 3D bndry_mask threshold 26.5 never masks; outflow_dirs overflows

[WeiqunZhang]

Severity: high/medium · Category: correctness, memory-ub · Fix order: 7 of 21 — fix this 7th.

_{Filenames are numbered in reverse fix order: 001 = fix last, 021 = fix first. This file is 015.}

Locations: Source/SyncRegister.cpp:265, Source/SyncRegister.cpp:73

Based on commit 9bf664bf (line numbers refer to that tree).

Two independent defects in SyncRegister::InitRHS, both introduced by the Fortran-to-C++/MLMG port: the threshold used to convert the accumulated bndry_mask counts to 0/1, and the size of the outflow_dirs stack array.

The defect

Source/SyncRegister.cpp:265 — The 3D bndry_mask threshold is AMREX_D_TERM(AMREX_SPACEDIM,*AMREX_SPACEDIM,*AMREX_SPACEDIM)-0.5 = 26.5, but the max node count is 2^AMREX_SPACEDIM = 8 (the deleted Fortran convertmask used 7.5), so in 3D the mask is never 0 and nodes interior to the fine-grid union are never masked out of the sync-projection RHS. _{Reported independently at this same line by F058; each reviewer's own wording and evidence is under Verification evidence below.}

Source/SyncRegister.cpp:73 — outflow_dirs is a stack array of size AMREX_SPACEDIM-1, but the loop below can record up to AMREX_SPACEDIM outflow directions, writing one element past the end of the array. _{Reported independently at this same line by F059; each reviewer's own wording and evidence is under Verification evidence below.}

Why it matters

F007: Any 3D run with >=2 AMR levels, do_sync_proj=1 (default), fine level chopped into multiple grids (normal max_grid_size). InitRHS leaves spurious one-sided residual contributions at fine-fine grid interface nodes in rhnd for MLsyncProject, so the sync correction to pressure/velocity is wrong and results change with grid decomposition. 2D (3.5) is correct; only 3D is broken.

F035: A multilevel run (InitRHS is reached via MLsyncProject) with outflow physical BCs in every coordinate direction, e.g. 2D with xlo/xhi outflow and yhi outflow: nOutflow becomes 2 and outflow_dirs[1] is written into a 1-element array -> stack out-of-bounds write (UB, possible corruption of adjacent locals). Nothing elsewhere forbids outflow in all directions (OutFlowBC supports 2*SPACEDIM outflow faces).

How to reach it

• Exec/run3d/regtest.3d.rayleightaylor: amr.max_level=2, 64^3 cells, default max_grid_size=32 chops fine levels into multiple grids; ns.do_sync_proj defaults to 1, so post_timestep -> level_sync -> MLsyncProject -> InitRHS executes every coarse step.

Suggested fix

Both are one-token corrections that restore the pre-port kernels.

The sum at lines 194-197 counts the 2^AMREX_SPACEDIM cells touching a node, not SPACEDIM^SPACEDIM, and the physical-boundary doubling at lines 206-249 is calibrated so that a node interior to the coarsened fine-grid union reaches exactly 8 in 3D (4x2 on a face, 2x4 on an edge, 1x8 at a corner). The deleted Fortran convertmask confirms the intent: 7.5D0 in SYNCREG_3D.F90, 3.5D0 in SYNCREG_2D.F90 (a3f0406^). AMREX_D_TERM(2.,*2.,*2.)-0.5 reproduces both. Maintainer call: 2D is bit-unchanged, but 3D multilevel sync-projection answers will move, so 3D regression baselines must be regenerated.

outflow_dirs needs one entry per axis, since the loop appends for any axis with an outflow face on either side. Note also that ={-1} initializes only element 0; harmless while nOutflow bounds the reads, but filling with -1 is clearer.

For Source/SyncRegister.cpp:265 (F007):

--- a/Source/SyncRegister.cpp
+++ b/Source/SyncRegister.cpp
@@ -262,7 +262,7 @@
             FArrayBox& fab   = fs[fsi];
             const Box& bx    = fab.box();
             auto const& mask = fab.array();
-            const Real maxcount = AMREX_D_TERM(AMREX_SPACEDIM,*AMREX_SPACEDIM,*AMREX_SPACEDIM) - 0.5;
+            const Real maxcount = AMREX_D_TERM(2.,*2.,*2.) - 0.5;
             amrex::ParallelFor(bx, [mask,maxcount]
             AMREX_GPU_DEVICE (int i, int j, int k) noexcept
             {

_{Threshold must be 2^AMREX_SPACEDIM - 0.5 (3.5 in 2D, 7.5 in 3D) to match the 2^SPACEDIM-cell sum built at lines 194-197. Verified against the deleted Fortran twins: SYNCREG_2D.F90 convertmask used 3.5D0, SYNCREG_3D.F90 used 7.5D0. 2D behavior is unchanged; 3D masking becomes active.}

For Source/SyncRegister.cpp:73 (F035):

--- a/Source/SyncRegister.cpp
+++ b/Source/SyncRegister.cpp
@@ -64,7 +64,7 @@
     const int* phys_lo = phys_bc.lo();
     const int* phys_hi = phys_bc.hi();
 
-    int outflow_dirs[AMREX_SPACEDIM-1]={-1};
+    int outflow_dirs[AMREX_SPACEDIM]={-1};
     int nOutflow = 0;
     for (int dir = 0; dir < AMREX_SPACEDIM; dir++)
     {

_{The loop below can record one direction per dimension, so the array needs AMREX_SPACEDIM slots; it was one short (and zero-length in 1D). Pure sizing fix, no behavior change unless outflow exists in every direction, where it removes the stack out-of-bounds write.}

_{Diff(s) are against 9bf664bf, written from the current source and verified only with git apply --check — never compiled, never run, never applied to the tree. Treat them as precise intent, not tested patches.}

Verification evidence

F007 — confirmed (two independent verifier lenses)

Lens 1 (refutation attempt): Line 265: maxcount = AMREX_D_TERM(AMREX_SPACEDIM,*AMREX_SPACEDIM,*AMREX_SPACEDIM) - 0.5 = 26.5 in 3D, but the mask sum at lines 194-197 is over 2^3=8 cells and boundary doubling (lines 206-249) caps it at 8, so mask > 26.5 is never true and 3D masking is a no-op. Deleted Fortran convertmask (a3f0406^:Source/Src_3d/SYNCREG_3D.F90) used identical 8-cell counting with threshold 7.5D0. Reachable: Projection.cpp:507 InitRHS in MLsyncProject, do_sync_proj default 1.

Lens 2 (reachability/intent): SyncRegister.cpp:265 threshold is 26.5 in 3D while max mask count is 8 (boundary doubling at lines 206-249 caps all nodes at 8), so mask≡1 and the multiply at lines 276-284 is a no-op. Deleted Fortran SYNCREG_3D.F90::convertmask (removed in 4a3d3b3) used 7.5; port commit 7f67172 introduced D_TERM(SPACEDIM,*SPACEDIM,*SPACEDIM), correct only in 2D. rhnd is consumed unmasked by MLNodeLaplacian::compRHS (rhs[ilev]+=rhnd, no covered-node mask).

F035 — confirmed (one verifier lens)

Lens 1 (refutation attempt): Line 67: int outflow_dirs[AMREX_SPACEDIM-1]={-1};; loop over all AMREX_SPACEDIM dirs writes outflow_dirs[nOutflow]=dir (line 73) per outflow direction. Outflow in every direction gives nOutflow==AMREX_SPACEDIM, writing one past the array (2D: index 1 into 1-int array). No code forbids it: NavierStokes.cpp BC parsing accepts pressure_outflow on any face; MacProj.cpp:863/Projection.cpp:1742 size outflow arrays as 2*AMREX_SPACEDIM.

F059 — confirmed (one verifier lens)

Reported as: outflow_dirs is declared with AMREX_SPACEDIM-1 elements but the loop below can record one direction per coordinate axis (up to AMREX_SPACEDIM entries), so line 73 writes past the end of the stack array when outflow BCs exist in every coordinate direction (2D: any run with outflow faces in both x and y).

Failure scenario: 2D multilevel run with do_sync_proj (default) and outflow in both directions, e.g. ns.lo_bc="1 2", ns.hi_bc="2 2" (open-domain jet): SyncRegister::InitRHS writes outflow_dirs[1] into a 1-int array -- stack corruption/UB (may clobber nOutflow, crash with stack protector, or silently mis-zero outflow nodes). 3D triggers with outflow in all three directions.

Lens 1 (refutation attempt): Duplicate of F035; same defect confirmed in current source. outflow_dirs[AMREX_SPACEDIM-1] (line 67) overflows at line 73 when outflow BCs exist in all coordinate directions (e.g. 2D ns.lo_bc/hi_bc with outflow in both x and y). Reached via Projection.cpp:507 (MLsyncProject -> InitRHS) in any multilevel run with do_sync_proj=1 (default, NavierStokesBase.cpp:133).

F058 — confirmed (one verifier lens)

Reported as: The bndry_mask conversion threshold AMREX_D_TERM(AMREX_SPACEDIM,*AMREX_SPACEDIM,*AMREX_SPACEDIM)-0.5 evaluates to 26.5 in 3D, but the maximum possible mask count is 8 (2^3, including physical-boundary doubling), so in 3D the mask never zeroes sync-register nodes interior to the union of coarsened fine grids.

Failure scenario: Any 3D multilevel run with do_sync_proj (default) whose fine level has >=2 grids: InitRHS leaves un-masked residual sums at nodes on faces shared by adjacent fine grids, injecting spurious (level-solve-residual-sized) values into the sync-projection RHS; grows if proj_tol is loosened. The deleted Fortran convertmask (SYNCREG_3D.F90) used 7.5; 2D coincidentally matches (3.5).

Lens 1 (refutation attempt): Duplicate of F007; confirmed. Max attainable mask value is 8 (interior node 8 cells; boundary nodes: 4x2 face, 2x4 edge, 1x8 corner), threshold 26.5 in 3D never triggers, so bndry_mask stays 1.0 everywhere and lines 276-284 multiply rhs by 1. Fortran history proves intent was 2^SPACEDIM-0.5 (7.5 in SYNCREG_3D.F90, 3.5 in SYNCREG_2D.F90); 2D matches only because 2*2==2^2.

---

_{Based on commit 9bf664bf, which is also the tree the audit verified against. From an automated audit of Source/, Tutorials/ and Util/. Audit finding ids: F007, F035, F059, F058. Reviewer unit(s): NavierStokes-2 + SyncRegister, theme:sync-reflux-subcycling. Nothing here was compiled or run — the failure scenarios are code reasoning, so the reaching configuration above is the cheapest way to confirm or refute it.}