Fill all halo cells for Flux, Value/Gradient, and Open BCs

[ewquon]

Summary

• The _fill_*_halo! functions for Flux, Value/Gradient/Mixed, and Open BCs previously filled only the first halo cell (h=1). High-order schemes requiring large halos (e.g., WENO9 with halo=5) had outer halo cells left unfilled (stale or zero).
• This fix loops over all halo cells (1:grid.Hx/Hy/Hz) for each BC type:
  • Flux BCs: call the per-cell flux fill for each halo depth
  • Value/Gradient/Mixed BCs: linearly extrapolate from the boundary into all halo cells
  • Open BCs: set boundary value, then zero-gradient extrapolate into deeper cells

This bug was discovered while debugging NaN crashes with WENO9 (halo=5) on a neutral ABL simulation (see #5485). Diagnostic analysis confirmed that ρw halos were all-zero throughout the simulation. While fixing the halos alone did not resolve the NaN crash (which has a separate Float32 precision root cause), the unfilled halos affect correctness of stencil computations near bounded boundaries for any scheme with halo > 1.

Test plan

• Existing boundary condition tests pass (these use small halos, so behavior is unchanged for halo=1)
• Verify halo fill with halo=(5,5,5): all 5 halo cells should be properly filled for Flux, Value, Gradient, and Open BCs
• Run WENO9 simulation and confirm halo cells are non-zero near Bounded boundaries

🤖 Generated with Claude Code

[simone-silvestri]

I am not sure we can call this a bug. This was a deliberate choice we made long time ago (back in #2035).
In general halos > 1 should not be touched on a bounded direction, we can change the behavior and fill those halos if there are other benefits though.

[simone-silvestri]

I think this might be a nice idea if we want to detangle the "ordering" of the BCs, which now follow a specific order because the corners are filled only by the Periodic and Distributed BCs.

[glwagner]

Right just to clarify, the description is wrong -- high-order schemes do not touch halo regions, because they are "stepped down" to lower order schemes approaching the boundary. We call these "buffer schemes" in the struct for Centered, Upwind, WENO.

I am wary about this because of the cost of filling points that are not used. I don't totally understand the point about boundary condition ordering.

[ewquon]

Okay let me verify that the halo “fix” has no effect due to the reduction in order at the boundaries then close this issue