From b3e0a62ba4d8c66b7cc40ab439b94835a5f4247c Mon Sep 17 00:00:00 2001 From: Weiqun Zhang Date: Wed, 26 Oct 2022 15:02:13 -0700 Subject: [PATCH] Pre- and Post-interpolation hook interface (#2991) Support both Fab and MultiFab versions of pre- and post-interpolation hooks. Because the pre-interp hook might modify the data, we need to make a copy to avoid modifying cached coarse data. Close #2989. --- Src/AmrCore/AMReX_FillPatchUtil_I.H | 69 +++++++++++++---------------- Src/AmrCore/AMReX_FillPatcher.H | 31 +++++++------ 2 files changed, 45 insertions(+), 55 deletions(-) diff --git a/Src/AmrCore/AMReX_FillPatchUtil_I.H b/Src/AmrCore/AMReX_FillPatchUtil_I.H index 8d8f210a0fe..3e94abfad27 100644 --- a/Src/AmrCore/AMReX_FillPatchUtil_I.H +++ b/Src/AmrCore/AMReX_FillPatchUtil_I.H @@ -4,6 +4,31 @@ namespace amrex { +namespace detail { + +template +auto call_interp_hook (F const& f, MF& mf, int icomp, int ncomp) + -> decltype(f(mf[0],Box(),icomp,ncomp)) +{ +#ifdef AMREX_USE_OMP +#pragma omp parallel if (Gpu::notInLaunchRegion()) +#endif + for (MFIter mfi(mf); mfi.isValid(); ++mfi) { + auto& dfab = mf[mfi]; + const Box& dbx = dfab.box(); + f(dfab, dbx, icomp, ncomp); + } +} + +template +auto call_interp_hook (F const& f, MF& mf, int icomp, int ncomp) + -> decltype(f(mf,icomp,ncomp)) +{ + f(mf, icomp, ncomp); +} + +} + template bool ProperlyNested (const IntVect& ratio, const IntVect& blocking_factor, int ngrow, const IndexType& boxType, Interp* mapper) @@ -459,9 +484,6 @@ namespace { if ( ! fpc.ba_crse_patch.empty()) { - - using FAB = typename MF::FABType::value_type; - MF mf_crse_patch = make_mf_crse_patch (fpc, ncomp, mf.boxArray().ixType()); // Must make sure fine exists under needed coarse faces. // It stores values for the final (interior) interpolation, @@ -491,20 +513,12 @@ namespace { solve_mask.setVal(1); // Values to solve. solve_mask.setVal(0, mask_cpc, 0, 1); // Known values. - for (MFIter mfi(mf_refined_patch); mfi.isValid(); ++mfi) - { - FAB& sfab = mf_crse_patch[mfi]; - pre_interp(sfab, sfab.box(), 0, ncomp); - } + detail::call_interp_hook(pre_interp, mf_crse_patch, 0, ncomp); InterpFace(mapper, mf_crse_patch, 0, mf_refined_patch, 0, ncomp, ratio, solve_mask, cgeom, fgeom, bcscomp, RunOn::Gpu, bcs); - for (MFIter mfi(mf_refined_patch); mfi.isValid(); ++mfi) - { - FAB& dfab = mf_refined_patch[mfi]; - post_interp(dfab, dfab.box(), 0, ncomp); - } + detail::call_interp_hook(post_interp, mf_refined_patch, 0, ncomp); bool aliasing = false; for (auto const& fmf_a : fmf) { @@ -538,30 +552,14 @@ namespace { MF mf_fine_patch = make_mf_fine_patch(fpc, ncomp); -#ifdef AMREX_USE_OMP -#pragma omp parallel if (Gpu::notInLaunchRegion()) -#endif - for (MFIter mfi(mf_crse_patch); mfi.isValid(); ++mfi) - { - auto& sfab = mf_crse_patch[mfi]; - const Box& sbx = sfab.box(); - pre_interp(sfab, sbx, 0, ncomp); - } + detail::call_interp_hook(pre_interp, mf_crse_patch, 0, ncomp); FillPatchInterp(mf_fine_patch, 0, mf_crse_patch, 0, ncomp, IntVect(0), cgeom, fgeom, amrex::grow(amrex::convert(fgeom.Domain(),mf.ixType()),nghost), ratio, mapper, bcs, bcscomp); -#ifdef AMREX_USE_OMP -#pragma omp parallel if (Gpu::notInLaunchRegion()) -#endif - for (MFIter mfi(mf_fine_patch); mfi.isValid(); ++mfi) - { - auto& dfab = mf_fine_patch[mfi]; - const Box& dbx = dfab.box(); - post_interp(dfab, dbx, 0, ncomp); - } + detail::call_interp_hook(post_interp, mf_fine_patch, 0, ncomp); mf.ParallelCopy(mf_fine_patch, 0, dcomp, ncomp, IntVect{0}, nghost); } @@ -1024,14 +1022,7 @@ InterpFromCoarseLevel (MF& mf, IntVect const& nghost, Real time, cbc(mf_crse_patch, 0, ncomp, mf_crse_patch.nGrowVect(), time, cbccomp); -#ifdef AMREX_USE_OMP -#pragma omp parallel if (Gpu::notInLaunchRegion()) -#endif - for (MFIter mfi(mf_crse_patch); mfi.isValid(); ++mfi) - { - FAB& sfab = mf_crse_patch[mfi]; - pre_interp(sfab, sfab.box(), 0, ncomp); - } + detail::call_interp_hook(pre_interp, mf_crse_patch, 0, ncomp); FillPatchInterp(mf, dcomp, mf_crse_patch, 0, ncomp, nghost, cgeom, fgeom, fdomain_g, ratio, mapper, bcs, bcscomp); diff --git a/Src/AmrCore/AMReX_FillPatcher.H b/Src/AmrCore/AMReX_FillPatcher.H index 22b14d35c0d..d0e775416ee 100644 --- a/Src/AmrCore/AMReX_FillPatcher.H +++ b/Src/AmrCore/AMReX_FillPatcher.H @@ -330,26 +330,25 @@ FillPatcher::fillCoarseFineBoundary (MF& mf, IntVect const& nghost, Real tim } } - MF mf_crse_patch; + if (m_cf_crse_data_tmp == nullptr) { + m_cf_crse_data_tmp = std::make_unique + (make_mf_crse_patch(fpc, m_ncomp)); + } + if (m_cf_crse_data.size() > 0 && amrex::almostEqual(time, m_cf_crse_data[0].first,5)) { - mf_crse_patch = MF(*m_cf_crse_data[0].second, amrex::make_alias, - scomp, ncomp); + amrex::Copy(*m_cf_crse_data_tmp, *m_cf_crse_data[0].second, + scomp, 0, ncomp, 0); } else if (m_cf_crse_data.size() > 1 && amrex::almostEqual(time, m_cf_crse_data[1].first,5)) { - mf_crse_patch = MF(*m_cf_crse_data[1].second, amrex::make_alias, - scomp, ncomp); + amrex::Copy(*m_cf_crse_data_tmp, *m_cf_crse_data[1].second, + scomp, 0, ncomp, 0); } else if (m_cf_crse_data.size() == 2) { - if (m_cf_crse_data_tmp == nullptr) { - m_cf_crse_data_tmp = std::make_unique - (make_mf_crse_patch(fpc, m_ncomp)); - } - mf_crse_patch = MF(*m_cf_crse_data_tmp, amrex::make_alias, scomp, ncomp); int const ng_space_interp = 8; // Need to be big enough Box domain = m_cgeom.growPeriodicDomain(ng_space_interp); domain.convert(mf.ixType()); @@ -358,10 +357,10 @@ FillPatcher::fillCoarseFineBoundary (MF& mf, IntVect const& nghost, Real tim Real alpha = (t1-time)/(t1-t0); Real beta = (time-t0)/(t1-t0); AMREX_ASSERT(alpha >= 0._rt && beta >= 0._rt); - auto const& a = mf_crse_patch.arrays(); + auto const& a = m_cf_crse_data_tmp->arrays(); auto const& a0 = m_cf_crse_data[0].second->const_arrays(); auto const& a1 = m_cf_crse_data[1].second->const_arrays(); - amrex::ParallelFor(mf_crse_patch, IntVect(0), ncomp, + amrex::ParallelFor(*m_cf_crse_data_tmp, IntVect(0), ncomp, [=] AMREX_GPU_DEVICE (int bi, int i, int j, int k, int n) noexcept { if (domain.contains(i,j,k)) { @@ -377,17 +376,17 @@ FillPatcher::fillCoarseFineBoundary (MF& mf, IntVect const& nghost, Real tim amrex::Abort("FillPatcher: High order interpolation in time not supported. Or FillPatcher was not properly deleted."); } - cbc(mf_crse_patch, 0, ncomp, nghost, time, cbccomp); + cbc(*m_cf_crse_data_tmp, 0, ncomp, nghost, time, cbccomp); - pre_interp(mf_crse_patch, 0, ncomp); + detail::call_interp_hook(pre_interp, *m_cf_crse_data_tmp, 0, ncomp); - FillPatchInterp(*m_cf_fine_data, scomp, mf_crse_patch, 0, + FillPatchInterp(*m_cf_fine_data, scomp, *m_cf_crse_data_tmp, 0, ncomp, IntVect(0), m_cgeom, m_fgeom, amrex::grow(amrex::convert(m_fgeom.Domain(), mf.ixType()),nghost), m_ratio, m_interp, bcs, bcscomp); - post_interp(*m_cf_fine_data, scomp, ncomp); + detail::call_interp_hook(post_interp, *m_cf_fine_data, scomp, ncomp); mf.ParallelCopy(*m_cf_fine_data, scomp, dcomp, ncomp, IntVect{0}, nghost); }