From c2c067fef980af8e76fc260bc73d1a40a2f3c632 Mon Sep 17 00:00:00 2001
From: Weiqun Zhang <weiqunzhang@lbl.gov>
Date: Mon, 27 May 2024 07:51:04 -0700
Subject: [PATCH] Fix bug introduced in #3959

---
 Src/Boundary/AMReX_InterpBndryData.H       | 33 ++++++++--------------
 Src/LinearSolvers/MLMG/AMReX_MLCellLinOp.H | 10 +++++--
 2 files changed, 19 insertions(+), 24 deletions(-)

diff --git a/Src/Boundary/AMReX_InterpBndryData.H b/Src/Boundary/AMReX_InterpBndryData.H
index 0a561fae71b..6a549883464 100644
--- a/Src/Boundary/AMReX_InterpBndryData.H
+++ b/Src/Boundary/AMReX_InterpBndryData.H
@@ -104,7 +104,7 @@ public:
      * \param max_order interpolation order
      */
     void updateBndryValues (BndryRegisterT<MF>& crse, int c_start, int bnd_start, int num_comp,
-                            const IntVect& ratio, int max_order = IBD_max_order_DEF);
+                            const IntVect& ratio, int max_order = IBD_max_order_DEF, int max_width = 2);
 
     //! Set boundary values to zero
     void setHomogValues ();
@@ -127,15 +127,14 @@ template <typename MF>
 void
 InterpBndryDataT<MF>::setPhysBndryValues (const MF& mf, int mf_start, int bnd_start, int num_comp)
 {
-    AMREX_ASSERT(mf.empty() || this->grids == mf.boxArray());
+    AMREX_ASSERT(this->grids == mf.boxArray());
 
     const Box& fine_domain = this->geom.Domain();
 
 #ifdef AMREX_USE_OMP
 #pragma omp parallel if (Gpu::notInLaunchRegion())
 #endif
-    for (MFIter mfi(this->grids, this->DistributionMap(),
-                    MFItInfo().SetDynamic(true)); mfi.isValid(); ++mfi)
+    for (MFIter mfi(mf,MFItInfo().SetDynamic(true)); mfi.isValid(); ++mfi)
     {
         const Box& bx = mfi.validbox();
         for (OrientationIter fi; fi; ++fi) {
@@ -144,22 +143,14 @@ InterpBndryDataT<MF>::setPhysBndryValues (const MF& mf, int mf_start, int bnd_st
             {
                 // Physical bndry, copy from grid.
                 auto      & bnd_fab = this->bndry[face][mfi];
+                auto const& src_fab = mf[mfi];
                 auto const& bnd_array = bnd_fab.array();
-                if (mf.empty()) {
-                    const Box& b = bnd_fab.box();
-                    AMREX_HOST_DEVICE_PARALLEL_FOR_4D ( b, num_comp, i, j, k, n,
-                    {
-                        bnd_array(i,j,k,n+bnd_start) = value_type(0);
-                    });
-                } else {
-                    auto const& src_fab = mf[mfi];
-                    auto const& src_array = src_fab.const_array();
-                    const Box& b = src_fab.box() & bnd_fab.box();
-                    AMREX_HOST_DEVICE_PARALLEL_FOR_4D ( b, num_comp, i, j, k, n,
-                    {
-                        bnd_array(i,j,k,n+bnd_start) = src_array(i,j,k,n+mf_start);
-                    });
-                }
+                auto const& src_array = src_fab.const_array();
+                const Box& b = src_fab.box() & bnd_fab.box();
+                AMREX_HOST_DEVICE_PARALLEL_FOR_4D ( b, num_comp, i, j, k, n,
+                {
+                    bnd_array(i,j,k,n+bnd_start) = src_array(i,j,k,n+mf_start);
+                });
             }
         }
     }
@@ -270,10 +261,10 @@ InterpBndryDataT<MF>::setBndryValues (BndryRegisterT<MF> const& crse, int c_star
 template <typename MF>
 void
 InterpBndryDataT<MF>::updateBndryValues (BndryRegisterT<MF>& crse, int c_start, int bnd_start,
-                                         int num_comp, const IntVect& ratio, int max_order)
+                                         int num_comp, const IntVect& ratio, int max_order, int max_width)
 {
     MF foo(this->grids, this->bndry[0].DistributionMap(), 1, num_comp, MFInfo().SetAlloc(false));
-    setBndryValues(crse, c_start, foo, 0, bnd_start, num_comp, ratio, max_order);
+    setBndryValues(crse, c_start, foo, 0, bnd_start, num_comp, ratio, max_order, max_width);
 }
 
 template <typename MF>
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLCellLinOp.H b/Src/LinearSolvers/MLMG/AMReX_MLCellLinOp.H
index 2abc093cafd..536d4c82b04 100644
--- a/Src/LinearSolvers/MLMG/AMReX_MLCellLinOp.H
+++ b/Src/LinearSolvers/MLMG/AMReX_MLCellLinOp.H
@@ -505,11 +505,15 @@ MLCellLinOpT<MF>::setLevelBC (int amrlev, const MF* a_levelbcdata, const MF* rob
 
     const int ncomp = this->getNComp();
 
+    MF zero;
     IntVect ng(1);
     if (this->hasHiddenDimension()) { ng[this->hiddenDirection()] = 0; }
-    AMREX_ALWAYS_ASSERT(a_levelbcdata == nullptr || a_levelbcdata->nGrowVect().allGE(ng));
-
-    MF zero;
+    if (a_levelbcdata == nullptr) {
+        zero.define(this->m_grids[amrlev][0], this->m_dmap[amrlev][0], ncomp, ng);
+        zero.setVal(RT(0.0));
+    } else {
+        AMREX_ALWAYS_ASSERT(a_levelbcdata->nGrowVect().allGE(ng));
+    }
     const MF& bcdata = (a_levelbcdata == nullptr) ? zero : *a_levelbcdata;
 
     IntVect br_ref_ratio(-1);