Skip to content

Commit

Permalink
Option to not fill ghost cells outside the domain
Browse files Browse the repository at this point in the history
  • Loading branch information
WeiqunZhang committed Sep 17, 2024
1 parent c01611a commit 66a84ee
Show file tree
Hide file tree
Showing 2 changed files with 23 additions and 5 deletions.
21 changes: 19 additions & 2 deletions Src/AmrCore/AMReX_FillPatchUtil.H
Original file line number Diff line number Diff line change
Expand Up @@ -702,11 +702,28 @@ namespace amrex
* this function. It's assumed that the fine level MultiFab mf's
* BoxArray is coarsenable by the refinement ratio. There is no support
* for EB.
*
* \tparam MF the MultiFab/FabArray type
* \tparam Interp spatial interpolater
*
* \param mf destination MF on the fine level
* \param nghost number of ghost cells of mf inside domain needed to be filled
* \param nghost_outside_domain number of ghost cells of mf outside domain needed to be filled
* \param cmf source MF on the coarse level
* \param scomp starting component of the source MF
* \param dcomp starting component of the destination MF
* \param ncomp number of components
* \param cgeom Geometry for the coarse level
* \param fgeom Geometry for the fine level
* \param ratio refinement ratio
* \param mapper spatial interpolater
* \param bcs boundar types for each component
*/
template <typename MF, typename Interp>
std::enable_if_t<IsFabArray<MF>::value>
InterpFromCoarseLevel (MF& mf, IntVect const& nghost, const MF& cmf,
int scomp, int dcomp, int ncomp,
InterpFromCoarseLevel (MF& mf, IntVect const& nghost,
IntVect const& nghost_outside_domain,
const MF& cmf, int scomp, int dcomp, int ncomp,
const Geometry& cgeom, const Geometry& fgeom,
const IntVect& ratio, Interp* mapper,
const Vector<BCRec>& bcs);
Expand Down
7 changes: 4 additions & 3 deletions Src/AmrCore/AMReX_FillPatchUtil_I.H
Original file line number Diff line number Diff line change
Expand Up @@ -1212,8 +1212,9 @@ InterpFromCoarseLevel (Array<MF*, AMREX_SPACEDIM> const& mf, IntVect const& ngho

template <typename MF, typename Interp>
std::enable_if_t<IsFabArray<MF>::value>
InterpFromCoarseLevel (MF& mf, IntVect const& nghost, const MF& cmf,
int scomp, int dcomp, int ncomp,
InterpFromCoarseLevel (MF& mf, IntVect const& nghost,
IntVect const& nghost_outside_domain,
const MF& cmf, int scomp, int dcomp, int ncomp,
const Geometry& cgeom, const Geometry& fgeom,
const IntVect& ratio, Interp* mapper,
const Vector<BCRec>& bcs)
Expand All @@ -1234,7 +1235,7 @@ InterpFromCoarseLevel (MF& mf, IntVect const& nghost, const MF& cmf,
cgeom.periodicity());

Box fdomain_g = amrex::convert(fgeom.Domain(),typ);
fdomain_g.grow(nghost);
fdomain_g.grow(nghost_outside_domain);
FillPatchInterp(mf, dcomp, mf_crse_patch, 0, ncomp, nghost, cgeom, fgeom,
fdomain_g, ratio, mapper, bcs, 0);
}
Expand Down

0 comments on commit 66a84ee

Please sign in to comment.