How to distinguish interior vs boundary cells in ParallelFor

[Afraid31]

I'm working with cell-centered MultiFabs in AMReX and would like to distinguish between:

Interior valid cells, which are far from domain boundaries or coarse/fine interfaces;

Boundary valid cells, which are near physical or AMR level boundaries.

Here is a simplified version of the code I'm using:

MultiFab& macro = M_old[lev];
for (MFIter mfi(macro); mfi.isValid(); ++mfi)
{
    const Box& bx = mfi.validbox();
    const Array4<Real>& macro_ = macro.array(mfi);

    amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j)
    {
        if (/* cell (i,j) is interior */)
            macro_(i, j) = myfun1().;
        else
            macro_(i, j) = myfun2().;
    });
}

My question:
What is the recommended way to determine, inside ParallelFor, whether a valid cell (i,j) is part of the interior of the full level or belongs to a boundary?

I want to avoid treating shared edges between adjacent Boxes as "boundary";

The boundary definition should reflect true physical and coarse/fine AMR boundaries.

Additional context:
One possible idea is to treat boundary cells as ghost cells and keep only interior cells as valid. This would make the classification inside ParallelFor much easier.

However, in my application I need to average values from fine grid boundary cells into coarse grid cells, e.g., by calling amrex::average_down(...). But it seems that average_down only works with valid cells, and will not include ghost cells in the averaging.

Is there any recommended way to:

Clearly distinguish boundary and interior cells;

Possibly include ghost cell values in an average-down operation, or

Modify the valid/ghost cell layout to support this while keeping AMReX functionality?

Thanks in advance for any suggestions or examples!

[WeiqunZhang]

You can build a mask. For example,

iMultiFab maskg1(macro.boxArray(), macro.DistributionMap(), 1, 1);
maskg1.BuildMask(geom.Domain(), geom.perioicity(), 1, 0, 0, 1);

// You can probably cache the mask                                                                                                                                           
iMultiFab mask(maskg1.boxArray(), maskg1.DistributionMap(), 1, 0);
auto const& ma = mask.arrays();
auto const& ma1 = maskg1.const_arrays();
ParallelFor(mask, [=] AMREX_GPU_DEVICE (int b, int i, int j, int k)
{
    auto const& a1 = ma1[b];
    auto const& a = ma[b];
    // assuming this is 2d                                                                                                                                                   
    if (a.contains(i-1,j-1,0)) {
        a(i,j,k) = 1;
    } else {
        bool interior = true;
        for (int joff = -1; joff <= 1 && interior; ++joff) {
            for (int ioff = -1; ioff <= 1 && interior; ++ioff) {
                interior = interior && a1(i+ioff,j+off,k);
            }
        }
        a(i,j,k) = interior ? 1 : 0;
    }
});

As for average down, you can write one yourself. Something like

// assuming IntVect rr is refinement ratio and IntVect ng is the number of                                                                                                   
// fine ghost cells involved in average_down.                                                                                                                                
MultiFab cfine(amrex::coarsen(fine.boxArray(),rr), fine.DistributionMap(),
               ncomp, ng/rr); // ng must multiple of rr                                                                                                                      
auto const& cma = cfine.arrays();
auto const& fma = fine.const_arrays();
ParallelFor(cma, cfine.nGrowVect(), ncomp,
            [=] AMREX_GPU_DEVICE (int b, int i, int j, int k)
            {
                amrex_avgdown(i,j,k,n,cma[b],fma[b],0,0,rr);
            });
coarse.ParallelCopy(cfine, 0, 0, ncomp, cfine.nGrowVect(), IntVect(0), geom.perioicity());