ApplyPBLHSmoothing: give the stencil the halo it reads, and cover it with a tiling parity test

[asalmgren]

Closes #3978.

The bug

ApplyPBLHSmoothing applies a 5-point stencil to the diagnosed PBL height, so every cell it writes depends on a column one further out. It looped over the whole planar work box and clamped only at the domain edge, so at every other edge of that box it read one column past the end of the FArrayBox.

The array it reads is allocated on PerpendicularBox<ZDir>(mfi.growntilebox(IntVect(1,1,0))), and growntilebox() is the trap: it does not grow at an interior tile edge, because the grown tile boxes have to keep partitioning the grown box. A tile in the middle of a box therefore has no halo at all, and the stencil reads outside the array on every one of its edges. With AMReX_BOUND_CHECK the run aborts; without it, whatever is next in the arena is mixed into the PBL height, which then feeds the entire K-profile.

On the ABL_MRF_Tiling deck with erf.enable_pblh_smoothing=true, tiled vs untiled:

 theta          0.5136862253    0.001666147717
 Kmv             4.420542564       0.2409735396
 Khv             5.196222304       0.2489480896
 Lturb           24.51525397       0.1206691032

Lturb (the PBL height) 24.5 m apart, 12% relative; Kmv/Khv ~24%; and it propagates into theta and the velocities.

The fix

The caller now names the region it wants smoothed and hands over a work array carrying passes columns of halo, grown from mfi.tilebox() explicitly rather than through growntilebox(), with the PBLH passes filling that halo. Pass p writes the region grown by passes-1-p and reads the one grown by passes-p, so each pass consumes one column of halo and the last lands exactly on the tile box. Every read is in bounds and the result no longer depends on how the domain is divided into boxes or tiles. An AMREX_ALWAYS_ASSERT states the halo contract — it is what caught the growntilebox() trap in the first place.

The stencil also now takes the domain periodicity. It clamped (Neumann) at every domain edge, but across a periodic edge the halo holds the wrapped columns and should be read directly; clamping folded an edge into a domain that has none. Both tiling decks are periodic in x and y.

Growing the work box duplicates a little PBLH work in the tile overlaps, so it is done only when smoothing is on. With erf.enable_pblh_smoothing=false — the default — the work box is the gbx the kernels already used, and the answers are unchanged.

Multi-pass smoothing now errors instead of lying

erf.pblh_smoothing_passes > 1 cannot be supported this way: the state and surface-layer arrays carry one ghost column, so there is no second column to consume. It now aborts naming the limit:

erf.pblh_smoothing_passes = 3 needs 3 halo columns, but the state and
surface-layer arrays carry only 1; reduce erf.pblh_smoothing_passes to at most 1

instead of silently returning a decomposition-dependent answer. At 3 passes with no tiling at all, one box and four boxes currently differ by 7.9 m in Lturb and 7.6% in Kmv.

Nothing in the tree sets this option — no deck, test or doc mentions enable_pblh_smoothing. Supporting more passes exactly needs the PBLH diagnosis split into its own MFIter loop writing a ghosted planar MultiFab, a FillBoundary, then a second loop for the K-profile. That is a much larger change and is left out here; say the word if it is wanted.

Also

The include guard's #endif in ERF_PBLModels.H sat before ApplyPBLHSmoothing, leaving the function outside the guard. It survives only because no translation unit includes the header twice. Moved to the end of the file.

Tests

ABL_MRF_Tiling_Smooth and ABL_YSUNew_Tiling_Smooth run the MRF deck with smoothing on through the tiling parity harness added in #3972/#3991. MRF and YSUNew size and fill the halo in separate code, so both are registered rather than just one. Both fail on the unfixed code with the differences quoted above, and pass here bit-identical.

GNU 13.3, Release, ERF_ENABLE_MPI=ON, 2 ranks:

• ctest -L regression: 110/110, no gold file changes
• one box vs four boxes with smoothing on, no tiling: PLOTFILE AGREE, every norm exactly 0
• all five tiling parity tests pass

🤖 Generated with Claude Code