MLMG: Robin a-coeff terms and cross-type getFluxes write past single-component temporaries

[WeiqunZhang]

Severity: high/medium · Category: memory-ub · Subsystem: LinearSolvers (blast radius 4/5)

Locations: Src/LinearSolvers/MLMG/AMReX_MLEBABecLap.cpp:1555, Src/LinearSolvers/MLMG/AMReX_MLMG.H:1153

Based on commit 248fcfb7a5 (branch ai_audit; line numbers refer to that tree).

Two unrelated sites in the MLMG solvers size a buffer for one component and then hand it to a writer that emits more: detail::applyRobinBCTermsCoeffs adds per-component Robin terms into the deliberately single-component m_a_coeffs, and the type-conversion branches of MLMGT::getFluxes(Vector<AMF*>,Location) define ncomp-wide cell-centered/nodal temporaries for an AMREX_SPACEDIM-wide result.

The defect

Src/LinearSolvers/MLMG/AMReX_MLEBABecLap.cpp:1555 — With ncomp>1 and any Robin domain BC, detail::applyRobinBCTermsCoeffs writes afab(i,j,k,icomp) for icomp>=1 into m_a_coeffs, which is defined with only 1 component, an out-of-bounds heap write.

Src/LinearSolvers/MLMG/AMReX_MLMG.H:1153 — In the AMF!=MF (type-conversion) branch of getFluxes(Vector<AMF*>,Location), the temporary cell-centered flux MFs are defined with ncomp (=1) components, but the inner call writes AMREX_SPACEDIM components via (EB_)average_face_to_cellcenter (which asserts cc.nComp()>=dcomp+AMREX_SPACEDIM); same undersized define at line 1218 for the nodal branch (MLNodeLaplacian::getFluxes asserts nComp>=AMREX_SPACEDIM), and the copy-backs at 1157/1222 return only 1 of AMREX_SPACEDIM components.

Why it matters

F049: Any build: MLEBABecLap(...,a_ncomp=2), setDomainBC with LinOpBCType::Robin on any face for any component, setLevelBC(lev, bc, robin_a, robin_b, robin_f), MLMG::solve -> prepareForSolve -> applyRobinBCTermsCoeffs. m_a_coeffs[amrlev][0] is allocated with 1 component (AMReX_MLEBABecLap.cpp:115-117) but afab(...,icomp) with icomp=1 writes past the FAB's data -> silent arena/heap corruption in release (Array4 bound-check abort only in debug). Multi-component Robin is supported infrastructure (m_robin_bcval sized ncomp*3, per-icomp m_lobc_orig), and MLABecLaplacian shares the same defect.

F228: 2D/3D, single-precision solver with double user data (API exercised by Tests/LinearSolvers/ABecLap_SP): MLPoissonT + MLMGT; solve() with MultiFab sol/rhs; then getFluxes(Vector<MultiFab*>{&flux}, Location::CellCenter) with flux having AMREX_SPACEDIM comps. Debug: AMREX_ASSERT(nComp>=AMREX_SPACEDIM) aborts. Release: amrex_avg_fc_to_cc writes components 1..AMREX_SPACEDIM-1 past the 1-component fab allocation - heap corruption.

How to reach it

• F049: config: Default CMake (any SPACEDIM, no EB needed; MLEBABecLap variant needs AMReX_EB=ON). Built by every linux/macos/cuda CI job.
• F049: API path: MLABecLaplacian op(geom,grids,dmap,LPInfo(),{},2) (or MLEBABecLap(...,factory,2)); op.setDomainBC({Robin,..},{..}); op.setLevelBC(0,&sol,&ra,&rb,&rf); MLMG(op).solve(...) -> prepareForSolve -> applyRobinBCTermsCoeffs. Public API; no in-tree caller uses Robin with ncomp>1 (grep 'Robin' in Tests/ is empty).
• F049: test coverage: none — no Tests/ program uses Robin BC at all; Tests/LinearSolvers/ABecLaplacian_C is ncomp=1 Dirichlet/Neumann.
• F049: downstream use: incflo uses Robin only with 1-component ops: src/diffusion/DiffusionScalarOp.cpp:38,79 (default ncomp=1, per-comp setDomainBC at :375,:402,:661,:717); get_diffuse_tensor_bc (src/diffusion/incflo_diffusion.cpp:130-178) aborts on BC::mixed so MLTensorOp never receives Robin. AMReX-Hydro, ERF, REMORA, WarpX: none found.
• F228: config: -DAMReX_SPACEDIM=2 or 3, CPU, default double (MultiFab and fMultiFab are always distinct types). gcc.yml has Debug 2D/3D and AMReX_ASSERTIONS=ON jobs (assert path) and Release jobs (OOB path); no test exercises it.
• F228: API path: MLPoissonT lp({geom},{ba},{dm}); MLMGT mlmg(lp); MultiFab sol(ba,dm,1,1), rhs(ba,dm,1,0); mlmg.solve(Vector<MultiFab*>{&sol}, Vector<MultiFab const*>{&rhs}, tol, 0); MultiFab flux(ba,dm,AMREX_SPACEDIM,0); mlmg.getFluxes(Vector<MultiFab*>{&flux}, MLMG::Location::CellCenter); // assert / heap overwrite. Public overload with NO in-tree or downstream caller (Tests use same-type Vector<MF*> or the Array<AMF*> variant, which sizes per-direction correctly).
• F228: test coverage: none — Tests/LinearSolvers/ABecLap_SP exercises cross-type solve() only; EBflux_grad uses the Array<fMultiFab*> getGradSolution variant; Nodal_Projection_EB uses same-type getFluxes.

Suggested fix

Robin terms (F049). m_a_coeffs is single-component on purpose: MLABecLaplacian::setACoeffs asserts "alpha is supposed to be single component" and every kernel indexes a(i,j,k) with no component. So afab(...,icomp) for icomp>=1 in the shared detail::applyRobinBCTermsCoeffs (AMReX_MLABecLaplacian.H) has nowhere legitimate to land. The attached diff takes the conservative route and aborts there for ncomp != 1, which covers both MLEBABecLap::applyRobinBCTermsCoeffs and MLABecLaplacian (and the tensor ops through inheritance). Maintainer decision: (a) keep this late abort; (b) fail earlier and more informatively by making supportRobinBC() in both ops return getNComp()==1, so MLLinOp::setDomainBC rejects the BC where it is set; or (c) implement the feature, meaning an ncomp-component m_a_coeffs and a(i,j,k,n) in every ABecLap/EBABecLap kernel plus relaxing the setACoeffs assert. Open PRs #4930 and #4922 touch this file and may be related (unmatched by the verifier). getFluxes conversion branches (F228). The same-type entry asserts nComp(*a_flux[0]) >= AMREX_SPACEDIM, average_face_to_cellcenter writes exactly AMREX_SPACEDIM components, and MLNodeLaplacian::getFluxes does the same, so the temporaries at MLMG.H:1153 and :1218 and the LocalCopy calls at :1157 and :1222 should use AMREX_SPACEDIM, not ncomp (and not AMREX_SPACEDIM*ncomp). Consider also asserting the user's AMF width in these branches. PR #5484 touches this file but not these lines. The two fixes are independent and can land in either order.

For Src/LinearSolvers/MLMG/AMReX_MLEBABecLap.cpp:1555 (F049):

--- a/Src/LinearSolvers/MLMG/AMReX_MLABecLaplacian.H
+++ b/Src/LinearSolvers/MLMG/AMReX_MLABecLaplacian.H
@@ -550,6 +550,8 @@
     using RT = typename LP::RT;
 
     const int ncomp = linop.getNComp();
+    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(ncomp == 1,
+                                     "Robin BC is supported only for a single component");
     bool reset_alpha = false;
     if (linop.m_a_scalar == RT(0.0)) {
         linop.m_a_scalar = RT(1.0);

_{Guards the shared detail::applyRobinBCTermsCoeffs (root cause), so both MLEBABecLap::applyRobinBCTermsCoeffs and MLABecLaplacian are covered; a 1-component alpha cannot hold per-component Robin terms. Maintainer decision: keep this abort, fail earlier via supportRobinBC() returning getNComp()==1 in both ops, or the real feature (ncomp-component m_a_coeffs plus a(i,j,k,n) in every kernel). Open PRs #4930/#4922 touch this area.}

For Src/LinearSolvers/MLMG/AMReX_MLMG.H:1153 (F228):

--- a/Src/LinearSolvers/MLMG/AMReX_MLMG.H
+++ b/Src/LinearSolvers/MLMG/AMReX_MLMG.H
@@ -1150,9 +1150,9 @@
         Vector<MF> fluxes(namrlevs);
         for (int ilev = 0; ilev < namrlevs; ++ilev) {
             auto const& amf = *a_flux[ilev];
-            fluxes[ilev].define(boxArray(amf), DistributionMap(amf), ncomp, 0);
+            fluxes[ilev].define(boxArray(amf), DistributionMap(amf), AMREX_SPACEDIM, 0);
         }
         getFluxes(GetVecOfPtrs(fluxes), GetVecOfPtrs(sol), a_loc);
         for (int ilev = 0; ilev < namrlevs; ++ilev) {
-            LocalCopy(*a_flux[ilev], fluxes[ilev], 0, 0, ncomp, IntVect(0));
+            LocalCopy(*a_flux[ilev], fluxes[ilev], 0, 0, AMREX_SPACEDIM, IntVect(0));
         }
@@ -1215,9 +1215,9 @@
             Vector<MF> fluxes(namrlevs);
             for (int ilev = 0; ilev < namrlevs; ++ilev) {
                 auto const& amf = *a_flux[ilev];
-                fluxes[ilev].define(boxArray(amf), DistributionMap(amf), ncomp, 0);
+                fluxes[ilev].define(boxArray(amf), DistributionMap(amf), AMREX_SPACEDIM, 0);
             }
             linop.getFluxes(GetVecOfPtrs(fluxes), GetVecOfPtrs(sol));
             for (int ilev = 0; ilev < namrlevs; ++ilev) {
-                LocalCopy(*a_flux[ilev], fluxes[ilev], 0, 0, ncomp, IntVect(0));
+                LocalCopy(*a_flux[ilev], fluxes[ilev], 0, 0, AMREX_SPACEDIM, IntVect(0));
             }

_{Uses AMREX_SPACEDIM (not AMREX_SPACEDIM*ncomp): average_face_to_cellcenter asserts fc nComp==1 and writes exactly AMREX_SPACEDIM comps, and MLNodeLaplacian::getFluxes zeroes/writes AMREX_SPACEDIM comps, matching the nComp(*a_flux[0]) >= AMREX_SPACEDIM assert at the same-type entry. Maintainer may also want that assert on the user's AMF array in these conversion branches. PR #5484 does not touch these lines.}

_{Diff(s) are against 248fcfb7a5, 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

F049 — confirmed (two independent verifier lenses)

Lens 1 (refutation attempt): MLEBABecLap.cpp:115-117 and MLABecLaplacian.H:350-352 define m_a_coeffs with 1 component. Shared detail::applyRobinBCTermsCoeffs (MLABecLaplacian.H:596-650) loops for icomp<ncomp and writes afab(i+1,j,k,icomp) += ... (line 606 and 5 siblings). Array4 index = nnstride with nstride=numPts, so icomp>=1 lands past the FAB allocation (debug: index_assert on n>=ncomp, Array4.H:111; release: silent OOB write). No ncomp guard: setDomainBC only aborts if !supportRobinBC(), and MLABecLaplacian/MLEBABecLap (and by inheritance MLTensorOp/MLEBTensorOp, ncomp=SPACEDIM) return true. m_robin_bcval is sized ncomp3 per icomp, so multi-comp Robin is accepted as input. No dev commit or PR fixes it.

Lens 2 (reachability/intent): m_a_coeffs is defined with 1 component (AMReX_MLEBABecLap.cpp:115-117; AMReX_MLABecLaplacian.H:350-352), while detail::applyRobinBCTermsCoeffs (AMReX_MLABecLaplacian.H:594-650) loops icomp<getNComp() and does afab(i+1,j,k,icomp) += .... Array4 offset n*nstride (nstride=numPts) puts icomp>=1 past the FAB; only AMREX_DEBUG catches it (AMReX_Array4.H:20-27). No ncomp guard: setDomainBC aborts only on !supportRobinBC() (MLLinOp.H:1535), which MLABecLaplacian/MLEBABecLap (and MLTensorOp/MLEBTensorOp via MLABecLaplacian::prepareForSolve, MLTensorOp.cpp:161) return true. Per-icomp robin storage (MLCellLinOp.H:894-935) shows multi-comp intent, not a documented limitation. Introduced by PR #1848 (87b182a), kept through #3035/#3617/#3788; PR bodies say nothing about ncomp. Not changed on development (git log ai_audit..development empty; development copies identical). Consequence: heap corruption plus lost Robin modification for comps>=1.

Possibly related upstream: dedup flagged open PRs #4930 #4922 as touching this file (verifier did not match it).

F228 — confirmed (one verifier lens)

Lens 1 (refutation attempt): MLMG.H:1153 fluxes[ilev].define(boxArray(amf), DistributionMap(amf), ncomp, 0) (ncomp=linop.getNComp()=1 for Poisson/ABecLap) then :1155 calls getFluxes(Vector<MF*>,Vector<MF*>,Location) whose :1175 AMREX_ASSERT(nComp(*a_flux[0]) >= AMREX_SPACEDIM) and :1201/1203 (EB_)average_face_to_cellcenter(*a_flux[alev],0,…) (MultiFabUtil.H:996 asserts cc.nComp()>=dcomp+AMREX_SPACEDIM; kernel MultiFabUtil_3D_C.H:46-48 writes cc(i,j,k,1..2) -> OOB on a 1-comp fab in Release). Nodal branch :1218 same define; MLNodeLaplacian_misc.cpp:936 asserts a_flux[0]->nComp()>=AMREX_SPACEDIM (reachable with MF=MultiFab, AMF=fMultiFab; MLMGT nodal would Abort in MLLinOp::getFluxes since MLNodeLinOp is MultiFab-only). Copy-backs :1157/:1222 copy only ncomp comps. Introduced with #3035/#3681; untouched by PR #5484 diff; not on development.

Possibly related upstream: dedup flagged open PR #5484 as touching this file (verifier did not match it).

---

_{Based on commit 248fcfb7a5, the tree the audit verified against. From an automated audit of AMReX Src/. Audit finding ids: F049, F228. Reviewer unit(s): LS/EBABecLap-1, LS/MLMG-driver. Nothing here was compiled or run except where the evidence says so — the failure scenarios are code reasoning, so the reaching configuration above is the cheapest way to confirm or refute it.}