AMReX_InterpBndryData.H

Go to the documentation of this file.

 
#ifndef AMREX_INTERPBNDRYDATA_H_
#define AMREX_INTERPBNDRYDATA_H_
#include <AMReX_Config.H>
 
#include <AMReX_BCRec.H>
#include <AMReX_BndryData.H>
#include <AMReX_InterpBndryData_K.H>
 
namespace amrex {
 
template <typename MF>
class InterpBndryDataT
    :
    public BndryDataT<MF>
{
public:
 
    using value_type = typename MF::value_type;
 
    InterpBndryDataT () noexcept = default;
 
    InterpBndryDataT (const BoxArray& _grids,
                      const DistributionMapping& _dmap,
                      int             _ncomp,
                      const Geometry& geom);
 
    ~InterpBndryDataT () = default;
 
    InterpBndryDataT (const InterpBndryDataT<MF>& rhs) = delete;
    InterpBndryDataT (InterpBndryDataT<MF>&& rhs) = delete;
    InterpBndryDataT<MF>& operator= (const InterpBndryDataT<MF>& rhs) = delete;
    InterpBndryDataT<MF>& operator= (InterpBndryDataT<MF>&& rhs) = delete;
 
    void setPhysBndryValues (const MF& mf, int mf_start, int bnd_start, int num_comp);
 
    void setBndryValues (BndryRegisterT<MF> const& crse, int c_start, const MF& fine, int f_start,
                         int bnd_start, int num_comp, const IntVect& ratio,
                         int max_order = IBD_max_order_DEF, int max_width = 2);
 
    void updateBndryValues (BndryRegisterT<MF>& crse, int c_start, int bnd_start, int num_comp,
                            const IntVect& ratio, int max_order = IBD_max_order_DEF, int max_width = 2);
 
    void setHomogValues ();
 
    // For sliding parabolic interp in bdfuncs
    static constexpr int IBD_max_order_DEF = 3;
 
};
 
template <typename MF>
InterpBndryDataT<MF>::InterpBndryDataT (const BoxArray& _grids,
                                        const DistributionMapping& _dmap,
                                        int             _ncomp,
                                        const Geometry& _geom)
    :
    BndryDataT<MF>(_grids,_dmap,_ncomp,_geom)
{}
 
template <typename MF>
void
InterpBndryDataT<MF>::setPhysBndryValues (const MF& mf, int mf_start, int bnd_start, int num_comp)
{
    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(mf,MFItInfo().SetDynamic(true)); mfi.isValid(); ++mfi)
    {
        const Box& bx = mfi.validbox();
        for (OrientationIter fi; fi; ++fi) {
            const Orientation face = fi();
            if (bx[face] == fine_domain[face] && !(this->geom.isPeriodic(face.coordDir())))
            {
                // 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();
                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);
                });
            }
        }
    }
}
 
template <typename MF>
void
InterpBndryDataT<MF>::setBndryValues (BndryRegisterT<MF> const& crse, int c_start, const MF& fine,
                                      int f_start, int bnd_start, int num_comp, const IntVect& ratio,
                                      int max_order, int max_width)
{
    AMREX_ASSERT(this->grids == fine.boxArray());
 
    const Box& fine_domain = this->geom.Domain();
 
    // Order-3 interpolation reads the boundary mask at tangential offsets of
    // ±ratio.  The masks are built with a fixed tangential half-width
    // (BndryDataT::NTangHalfWidth == 5, i.e. ref_ratio + 1), so ratios above 4
    // would walk past the allocated mask halo.  Order-1 does not touch the mask.
    if (max_order == 3) {
        AMREX_ALWAYS_ASSERT_WITH_MESSAGE(ratio.allLE(4),
            "InterpBndryDataT<MF>::setBndryValues: max_order=3 requires refinement ratio <= 4");
    }
 
    if (max_order==3 || max_order==1)
    {
        MFItInfo info;
        if (Gpu::notInLaunchRegion()) { info.SetDynamic(true); }
#ifdef AMREX_USE_OMP
#pragma omp parallel if (Gpu::notInLaunchRegion())
#endif
        for (MFIter mfi(fine,info); mfi.isValid(); ++mfi)
        {
            const Box& fine_bx  = mfi.validbox();
            for (OrientationIter fi; fi; ++fi) {
                const Orientation face = fi();
                const int dir = face.coordDir();
                if (fine_bx[face] != fine_domain[face] || this->geom.isPeriodic(dir))
                { // coarse/fine boundary: interpolate from coarse data stored in BndryRegister crse
                    auto const& crse_array = crse[face].const_array(mfi);
                    auto const& bdry_array = this->bndry[face].array(mfi);
                    Box const& b = this->bndry[face][mfi].box();
                    const auto rr = ratio.dim3();
 
                    if (max_order == 1)
                    {
                        AMREX_HOST_DEVICE_FOR_4D ( b, num_comp, it, jt, kt, n,
                        {
                            interpbndrydata_o1(it,jt,kt,n,bdry_array,bnd_start,
                                               crse_array,c_start,rr);
                        });
                    }
                    else
                    {
                        auto const& mask_array = this->masks[face].const_array(mfi);
                        int is_not_covered = BndryData::not_covered;
                        switch (dir) {
                        case 0:
                        {
                            AMREX_HOST_DEVICE_FOR_4D ( b, num_comp, it, jt, kt, n,
                            {
                                interpbndrydata_x_o3(it,jt,kt,n,bdry_array,bnd_start,
                                                     crse_array,c_start,rr,
                                                     mask_array, is_not_covered, max_width);
                            });
                            break;
                        }
#if (AMREX_SPACEDIM >= 2)
                        case 1:
                        {
                            AMREX_HOST_DEVICE_FOR_4D ( b, num_comp, it, jt, kt, n,
                            {
                                interpbndrydata_y_o3(it,jt,kt,n,bdry_array,bnd_start,
                                                     crse_array,c_start,rr,
                                                     mask_array, is_not_covered, max_width);
                            });
                            break;
                        }
#if (AMREX_SPACEDIM == 3)
                        case 2:
                        {
                            AMREX_HOST_DEVICE_FOR_4D ( b, num_comp, it, jt, kt, n,
                            {
                                interpbndrydata_z_o3(it,jt,kt,n,bdry_array,bnd_start,
                                                     crse_array,c_start,rr,
                                                     mask_array, is_not_covered, max_width);
                            });
                            break;
                        }
#endif
#endif
                        default: {}
                        }
                    }
                }
                else if (fine.defined(mfi)) {
                    // Physical bndry
                    auto      & bnd_fab = this->bndry[face][mfi];
                    auto const& src_fab = fine[mfi];
                    auto const& bnd_array = bnd_fab.array();
                    auto const& src_array = src_fab.const_array();
                    const Box& b = bnd_fab.box() & src_fab.box();
                    AMREX_HOST_DEVICE_FOR_4D ( b, num_comp, ii, jj, kk, nn,
                    {
                        bnd_array(ii,jj,kk,nn+bnd_start) = src_array(ii,jj,kk,nn+f_start);
                    });
                }
            }
        }
    }
    else
    {
        amrex::Abort("InterpBndryDataT<MF>::setBndryValues supports only max_order=1 or 3");
    }
}
 
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 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, max_width);
}
 
template <typename MF>
void
InterpBndryDataT<MF>::setHomogValues ()
{
    this->setVal(0);
}
 
using InterpBndryData = InterpBndryDataT<MultiFab>;
using fInterpBndryData = InterpBndryDataT<fMultiFab>;
 
}
 
#endif /*AMREX_INTERPBNDRYDATA_H_*/