Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix EB interpolation from cell centers to faces #1443

Merged
merged 1 commit into from
Oct 8, 2020
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
143 changes: 83 additions & 60 deletions Src/EB/AMReX_EBMultiFabUtil.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -122,10 +122,10 @@ EB_set_covered_faces (const Array<MultiFab*,AMREX_SPACEDIM>& umac, Real val)
{
AMREX_D_TERM(const Box& xbx = mfi.tilebox(IntVect::TheDimensionVector(0));,
const Box& ybx = mfi.tilebox(IntVect::TheDimensionVector(1));,
const Box& zbx = mfi.tilebox(IntVect::TheDimensionVector(2)););
const Box& zbx = mfi.tilebox(IntVect::TheDimensionVector(2)));
AMREX_D_TERM(Array4<Real> const& u = umac[0]->array(mfi);,
Array4<Real> const& v = umac[1]->array(mfi);,
Array4<Real> const& w = umac[2]->array(mfi););
Array4<Real> const& w = umac[2]->array(mfi));

auto fabtyp = flags[mfi].getType();
if (fabtyp == FabType::covered)
Expand Down Expand Up @@ -173,7 +173,7 @@ EB_set_covered_faces (const Array<MultiFab*,AMREX_SPACEDIM>& umac, Real val)
{
AMREX_D_TERM(Array4<Real const> const& ax = area[0]->const_array(mfi);,
Array4<Real const> const& ay = area[1]->const_array(mfi);,
Array4<Real const> const& az = area[2]->const_array(mfi););
Array4<Real const> const& az = area[2]->const_array(mfi));
AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM (
xbx, txbx,
{
Expand Down Expand Up @@ -238,10 +238,10 @@ EB_set_covered_faces (const Array<MultiFab*,AMREX_SPACEDIM>& umac, const int sco
{
AMREX_D_TERM(const Box& xbx = mfi.tilebox(IntVect::TheDimensionVector(0));,
const Box& ybx = mfi.tilebox(IntVect::TheDimensionVector(1));,
const Box& zbx = mfi.tilebox(IntVect::TheDimensionVector(2)););
const Box& zbx = mfi.tilebox(IntVect::TheDimensionVector(2)));
AMREX_D_TERM(Array4<Real> const& u = umac[0]->array(mfi);,
Array4<Real> const& v = umac[1]->array(mfi);,
Array4<Real> const& w = umac[2]->array(mfi););
Array4<Real> const& w = umac[2]->array(mfi));

auto fabtyp = flags[mfi].getType();
if (fabtyp == FabType::covered)
Expand Down Expand Up @@ -289,7 +289,7 @@ EB_set_covered_faces (const Array<MultiFab*,AMREX_SPACEDIM>& umac, const int sco
{
AMREX_D_TERM(Array4<Real const> const& ax = area[0]->const_array(mfi);,
Array4<Real const> const& ay = area[1]->const_array(mfi);,
Array4<Real const> const& az = area[2]->const_array(mfi););
Array4<Real const> const& az = area[2]->const_array(mfi));
AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM (
xbx, txbx,
{
Expand Down Expand Up @@ -718,7 +718,7 @@ void EB_computeDivergence (MultiFab& divu, const Array<MultiFab const*,AMREX_SPA
Array4<Real> const& divuarr = divu.array(mfi);
AMREX_D_TERM(Array4<Real const> const& uarr = umac[0]->const_array(mfi);,
Array4<Real const> const& varr = umac[1]->const_array(mfi);,
Array4<Real const> const& warr = umac[2]->const_array(mfi););
Array4<Real const> const& warr = umac[2]->const_array(mfi));

const auto fabtyp = flagfab.getType(bx);
if (fabtyp == FabType::covered) {
Expand All @@ -737,10 +737,10 @@ void EB_computeDivergence (MultiFab& divu, const Array<MultiFab const*,AMREX_SPA
Array4<Real const> const& vol = vfrac.const_array(mfi);
AMREX_D_TERM(Array4<Real const> const& apx = area[0]->const_array(mfi);,
Array4<Real const> const& apy = area[1]->const_array(mfi);,
Array4<Real const> const& apz = area[2]->const_array(mfi););
Array4<Real const> const& apz = area[2]->const_array(mfi));
AMREX_D_TERM(Array4<Real const> const& fcx = fcent[0]->const_array(mfi);,
Array4<Real const> const& fcy = fcent[1]->const_array(mfi);,
Array4<Real const> const& fcz = fcent[2]->const_array(mfi););
Array4<Real const> const& fcz = fcent[2]->const_array(mfi));
Array4<EBCellFlag const> const& flagarr = flagfab.const_array();
AMREX_HOST_DEVICE_FOR_4D(bx,divu.nComp(),i,j,k,n,
{
Expand Down Expand Up @@ -781,7 +781,7 @@ EB_average_face_to_cellcenter (MultiFab& ccmf, int dcomp,
Array4<Real> const& ccfab = ccmf.array(mfi);
AMREX_D_TERM(Array4<Real const> const& xfab = fmf[0]->const_array(mfi);,
Array4<Real const> const& yfab = fmf[1]->const_array(mfi);,
Array4<Real const> const& zfab = fmf[2]->const_array(mfi););
Array4<Real const> const& zfab = fmf[2]->const_array(mfi));
const auto fabtyp = flagfab.getType(bx);
if (fabtyp == FabType::covered) {
AMREX_HOST_DEVICE_FOR_3D(bx, i, j, k,
Expand All @@ -796,7 +796,7 @@ EB_average_face_to_cellcenter (MultiFab& ccmf, int dcomp,
} else {
AMREX_D_TERM(Array4<Real const> const& apx = area[0]->const_array(mfi);,
Array4<Real const> const& apy = area[1]->const_array(mfi);,
Array4<Real const> const& apz = area[2]->const_array(mfi););
Array4<Real const> const& apz = area[2]->const_array(mfi));
Array4<EBCellFlag const> const& flagarr = flagfab.const_array();
AMREX_HOST_DEVICE_FOR_3D(bx,i,j,k,
{
Expand Down Expand Up @@ -905,49 +905,61 @@ EB_interp_CC_to_FaceCentroid (const MultiFab& cc,
#endif
for (MFIter mfi(cc, mfi_info); mfi.isValid(); ++mfi)
{

const Box& vbx = mfi.tilebox();
AMREX_D_TERM(const Box& xbx = mfi.nodaltilebox(0);,
const Box& ybx = mfi.nodaltilebox(1);,
const Box& zbx = mfi.nodaltilebox(2));
const auto& ccfab = cc.array(mfi,scomp);
const auto fabtyp = flags[mfi].getType(amrex::grow(vbx,0));
const auto fabtyp_ghost = flags[mfi].getType(amrex::grow(vbx,nghost));

if (fabtyp != FabType::covered)
{

AMREX_D_TERM(Array4<Real> const& edg_x = fc_x.array(mfi,dcomp);,
Array4<Real> const& edg_y = fc_y.array(mfi,dcomp);,
Array4<Real> const& edg_z = fc_z.array(mfi,dcomp););
Array4<Real> const& edg_z = fc_z.array(mfi,dcomp));

if (fabtyp_ghost == FabType::regular )
{

AMREX_LAUNCH_HOST_DEVICE_LAMBDA ( vbx, thread_box,
{
eb_interp_cc2face(thread_box, ccfab,
AMREX_D_DECL(edg_x,edg_y,edg_z),
ncomp,
domain, d_bcs);
});

AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM
(xbx, txbx,
{
eb_interp_cc2face_x(txbx, ccfab, edg_x, ncomp, domain, d_bcs);
},
ybx, tybx,
{
eb_interp_cc2face_y(tybx, ccfab, edg_y, ncomp, domain, d_bcs);
},
zbx, tzbx,
{
eb_interp_cc2face_z(tzbx, ccfab, edg_z, ncomp, domain, d_bcs);
});
}
else
{
AMREX_D_TERM(Array4<Real const> const& apxfab = area[0]->const_array(mfi);,
Array4<Real const> const& apyfab = area[1]->const_array(mfi);,
Array4<Real const> const& apzfab = area[2]->const_array(mfi););
Array4<Real const> const& apzfab = area[2]->const_array(mfi));
AMREX_D_TERM(Array4<Real const> const& fcx = fcent[0]->const_array(mfi);,
Array4<Real const> const& fcy = fcent[1]->const_array(mfi);,
Array4<Real const> const& fcz = fcent[2]->const_array(mfi););

AMREX_LAUNCH_HOST_DEVICE_LAMBDA ( vbx, thread_box,
{
eb_interp_cc2facecent(thread_box, ccfab,
AMREX_D_DECL(apxfab,apyfab,apzfab),
AMREX_D_DECL(fcx,fcy,fcz),
AMREX_D_DECL(edg_x,edg_y,edg_z),
ncomp,
domain, d_bcs);
});
Array4<Real const> const& fcz = fcent[2]->const_array(mfi));

AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM
(xbx, txbx,
{
eb_interp_cc2facecent_x(txbx, ccfab, apxfab, fcx, edg_x,
ncomp, domain, d_bcs);
},
ybx, tybx,
{
eb_interp_cc2facecent_y(tybx, ccfab, apyfab, fcy, edg_y,
ncomp, domain, d_bcs);
},
zbx, tzbx,
{
eb_interp_cc2facecent_z(tzbx, ccfab, apzfab, fcz, edg_z,
ncomp, domain, d_bcs);
});
}
}
}
Expand Down Expand Up @@ -1032,8 +1044,10 @@ EB_interp_CellCentroid_to_FaceCentroid (const MultiFab& phi_centroid,
#endif
for (MFIter mfi(phi_centroid, mfi_info); mfi.isValid(); ++mfi)
{

const Box& vbx = mfi.tilebox();
AMREX_D_TERM(const Box& xbx = mfi.nodaltilebox(0);,
const Box& ybx = mfi.nodaltilebox(1);,
const Box& zbx = mfi.nodaltilebox(2));
const auto& ccfab = phi_centroid.array(mfi,scomp);
const auto fabtyp = flags[mfi].getType(amrex::grow(vbx,0));
const auto fabtyp_ghost = flags[mfi].getType(amrex::grow(vbx,nghost));
Expand All @@ -1042,44 +1056,53 @@ EB_interp_CellCentroid_to_FaceCentroid (const MultiFab& phi_centroid,
{
AMREX_D_TERM(Array4<Real> const& phi_x = phi_xface.array(mfi,dcomp);,
Array4<Real> const& phi_y = phi_yface.array(mfi,dcomp);,
Array4<Real> const& phi_z = phi_zface.array(mfi,dcomp););
Array4<Real> const& phi_z = phi_zface.array(mfi,dcomp));

if (fabtyp_ghost == FabType::regular )
{

AMREX_LAUNCH_HOST_DEVICE_LAMBDA ( vbx, thread_box,
{
eb_interp_cc2face(thread_box, ccfab,
AMREX_D_DECL(phi_x,phi_y,phi_z),
ncomp,
domain, d_bcs);
});

AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM
(xbx, txbx,
{
eb_interp_cc2face_x(txbx, ccfab, phi_x, ncomp, domain, d_bcs);
},
ybx, tybx,
{
eb_interp_cc2face_y(tybx, ccfab, phi_y, ncomp, domain, d_bcs);
},
zbx, tzbx,
{
eb_interp_cc2face_z(tzbx, ccfab, phi_z, ncomp, domain, d_bcs);
});
}
else
{
AMREX_D_TERM(Array4<Real const> const& apxfab = area[0]->const_array(mfi);,
Array4<Real const> const& apyfab = area[1]->const_array(mfi);,
Array4<Real const> const& apzfab = area[2]->const_array(mfi););
Array4<Real const> const& apzfab = area[2]->const_array(mfi));

AMREX_D_TERM(Array4<Real const> const& fcx = fcent[0]->const_array(mfi);,
Array4<Real const> const& fcy = fcent[1]->const_array(mfi);,
Array4<Real const> const& fcz = fcent[2]->const_array(mfi););
Array4<Real const> const& fcz = fcent[2]->const_array(mfi));

Array4<Real const> const& cvol = vfrac.const_array(mfi);
Array4<Real const> const& cct = ccent.const_array(mfi);

AMREX_LAUNCH_HOST_DEVICE_LAMBDA ( vbx, thread_box,
{
eb_interp_centroid2facecent(thread_box, ccfab,
AMREX_D_DECL(apxfab,apyfab,apzfab),
cvol,
cct,
AMREX_D_DECL(fcx,fcy,fcz),
AMREX_D_DECL(phi_x,phi_y,phi_z),
ncomp,
domain, d_bcs);
});

AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM
(xbx, txbx,
{
eb_interp_centroid2facecent_x(txbx, ccfab, apxfab, cvol, cct, fcx, phi_x,
ncomp, domain, d_bcs);
},
ybx, tybx,
{
eb_interp_centroid2facecent_y(tybx, ccfab, apyfab, cvol, cct, fcy, phi_y,
ncomp, domain, d_bcs);
},
zbx, tzbx,
{
eb_interp_centroid2facecent_z(tzbx, ccfab, apzfab, cvol, cct, fcz, phi_z,
ncomp, domain, d_bcs);
});
}
}
}
Expand Down
Loading