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.

Sign up for GitHub

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

Add option to interpolate data on faces linearly in the tangential di… #3794

Merged
merged 9 commits into from
Mar 12, 2024
Merged

Add option to interpolate data on faces linearly in the tangential di… #3794

Show file tree
Hide file tree
Changes from 6 commits
Commits
Show all changes
9 commits
Select commit Hold shift + click to select a range
543feb3
Add option to interpolate data on faces linearly in the tangential di…
asalmgren Mar 11, 2024
b1e5990
fix for compilation
asalmgren Mar 11, 2024
8d447be
remove if 0 accidentally left in
asalmgren Mar 11, 2024
a6466d2
fix oops
asalmgren Mar 11, 2024
6d7b1b8
add ignore_unused
asalmgren Mar 11, 2024
7f7aa67
allow FaceConservativeLinearInterpolater even if DIM == 1
asalmgren Mar 11, 2024
61f0e38
fixes based on PR feedback
asalmgren Mar 12, 2024
ab41da3
fix
asalmgren Mar 12, 2024
c41f8c9
revert to original versions of AMReX_InterpFaceReg*
asalmgren Mar 12, 2024
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
6 changes: 3 additions & 3 deletions Src/AmrCore/AMReX_InterpFaceReg_2D_C.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,12 +6,12 @@ namespace amrex {
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
void interp_face_reg (int i, int j, IntVect const& rr, Array4<Real> const& fine, int scomp,
Array4<Real const> const& crse, Array4<Real> const& slope, int ncomp,
Box const& domface, int idim)
Box const& per_grown_domain, int idim)
{
int ic = amrex::coarsen(i,rr[0]);
int jc = amrex::coarsen(j,rr[1]);
if (idim == 0) {
if (jc == domface.smallEnd(1) || jc == domface.bigEnd(1)) {
if (jc == per_grown_domain.smallEnd(1) || jc == per_grown_domain.bigEnd(1)) {
for (int n = 0; n < ncomp; ++n) {
fine(i,j,0,n+scomp) = crse(ic,jc,0,n);
}
Expand All @@ -35,7 +35,7 @@ void interp_face_reg (int i, int j, IntVect const& rr, Array4<Real> const& fine,
}
}
} else {
if (ic == domface.smallEnd(0) || ic == domface.bigEnd(0)) {
if (ic == per_grown_domain.smallEnd(0) || ic == per_grown_domain.bigEnd(0)) {
for (int n = 0; n < ncomp; ++n) {
fine(i,j,0,n+scomp) = crse(ic,jc,0,n);
}
Expand Down
15 changes: 7 additions & 8 deletions Src/AmrCore/AMReX_InterpFaceReg_3D_C.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ namespace amrex {
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
void interp_face_reg (int i, int j, int k, IntVect const& rr, Array4<Real> const& fine, int scomp,
Array4<Real const> const& crse, Array4<Real> const& slope, int ncomp,
Box const& domface, int idim)
Box const& per_grown_domain, int idim)
{
int ic = amrex::coarsen(i,rr[0]);
int jc = amrex::coarsen(j,rr[1]);
Expand All @@ -15,7 +15,7 @@ void interp_face_reg (int i, int j, int k, IntVect const& rr, Array4<Real> const
for (int n = 0; n < ncomp; ++n) {
fine(i,j,k,n+scomp) = crse(ic,jc,kc,n);
}
if (jc != domface.smallEnd(1) && jc != domface.bigEnd(1) && rr[1] > 1) {
if (jc != per_grown_domain.smallEnd(1) && jc != per_grown_domain.bigEnd(1) && rr[1] > 1) {
Real sfy = Real(1.0);
for (int n = 0; n < ncomp; ++n) {
Real dc = Real(0.5) * (crse(ic,jc+1,kc,n) - crse(ic,jc-1,kc,n));
Expand All @@ -34,8 +34,7 @@ void interp_face_reg (int i, int j, int k, IntVect const& rr, Array4<Real> const
fine(i,j,k,n+scomp) += yoff * slope(i,j,k,n) * sfy;
}
}

if (kc != domface.smallEnd(2) && kc != domface.bigEnd(2) && rr[2] > 1) {
if (kc != per_grown_domain.smallEnd(2) && kc != per_grown_domain.bigEnd(2) && rr[2] > 1) {
Real sfz = Real(1.0);
for (int n = 0; n < ncomp; ++n) {
Real dc = Real(0.5) * (crse(ic,jc,kc+1,n) - crse(ic,jc,kc-1,n));
Expand All @@ -58,7 +57,7 @@ void interp_face_reg (int i, int j, int k, IntVect const& rr, Array4<Real> const
for (int n = 0; n < ncomp; ++n) {
fine(i,j,k,n+scomp) = crse(ic,jc,kc,n);
}
if (ic != domface.smallEnd(0) && ic != domface.bigEnd(0) && rr[0] > 1) {
if (ic != per_grown_domain.smallEnd(0) && ic != per_grown_domain.bigEnd(0) && rr[0] > 1) {
Real sfx = Real(1.0);
for (int n = 0; n < ncomp; ++n) {
Real dc = Real(0.5) * (crse(ic+1,jc,kc,n) - crse(ic-1,jc,kc,n));
Expand All @@ -78,7 +77,7 @@ void interp_face_reg (int i, int j, int k, IntVect const& rr, Array4<Real> const
}
}

if (kc != domface.smallEnd(2) && kc != domface.bigEnd(2) && rr[2] > 1) {
if (kc != per_grown_domain.smallEnd(2) && kc != per_grown_domain.bigEnd(2) && rr[2] > 1) {
Real sfz = Real(1.0);
for (int n = 0; n < ncomp; ++n) {
Real dc = Real(0.5) * (crse(ic,jc,kc+1,n) - crse(ic,jc,kc-1,n));
Expand All @@ -101,7 +100,7 @@ void interp_face_reg (int i, int j, int k, IntVect const& rr, Array4<Real> const
for (int n = 0; n < ncomp; ++n) {
fine(i,j,k,n+scomp) = crse(ic,jc,kc,n);
}
if (ic != domface.smallEnd(0) && ic != domface.bigEnd(0) && rr[0] > 1) {
if (ic != per_grown_domain.smallEnd(0) && ic != per_grown_domain.bigEnd(0) && rr[0] > 1) {
Real sfx = Real(1.0);
for (int n = 0; n < ncomp; ++n) {
Real dc = Real(0.5) * (crse(ic+1,jc,kc,n) - crse(ic-1,jc,kc,n));
Expand All @@ -121,7 +120,7 @@ void interp_face_reg (int i, int j, int k, IntVect const& rr, Array4<Real> const
}
}

if (jc != domface.smallEnd(1) && jc != domface.bigEnd(1) && rr[1] > 1) {
if (jc != per_grown_domain.smallEnd(1) && jc != per_grown_domain.bigEnd(1) && rr[1] > 1) {
Real sfy = Real(1.0);
for (int n = 0; n < ncomp; ++n) {
Real dc = Real(0.5) * (crse(ic,jc+1,kc,n) - crse(ic,jc-1,kc,n));
Expand Down
12 changes: 10 additions & 2 deletions Src/AmrCore/AMReX_InterpFaceRegister.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,14 @@ void InterpFaceRegister::define (BoxArray const& fba, DistributionMapping const&

m_crse_geom = amrex::coarsen(m_fine_geom, m_ref_ratio);

// We don't need to worry about face-based domain because this is only used in the tangential interpolation
Box per_grown_domain = m_crse_geom.Domain();
for (int dim = 0; dim < AMREX_SPACEDIM; dim++) {
if (m_crse_geom.isPeriodic(dim)) {
per_grown_domain.grow(dim,1);
}
}

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
// We don't need to worry about face-based domain because this is only used in the tangential interpolation
Box per_grown_domain = m_crse_geom.Domain();
for (int dim = 0; dim < AMREX_SPACEDIM; dim++) {
if (m_crse_geom.isPeriodic(dim)) {
per_grown_domain.grow(dim,1);
}
}

Don't see it being used anywhere.

constexpr int crse_fine_face = 1;
constexpr int fine_fine_face = 0;
// First, set the face mask to 1 (i.e., coarse/fine boundary),
Expand Down Expand Up @@ -91,7 +99,7 @@ namespace {
Array4<Real> slope;
Array4<Real const> crse;
Array4<int const> mask;
Box domface;
Box per_grown_domain;
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
Box box() const noexcept { return Box(mask); }
};
Expand Down Expand Up @@ -151,7 +159,7 @@ InterpFaceRegister::interp (Array<MultiFab*, AMREX_SPACEDIM> const& fine, // NOL
{
if (tag.mask(i,j,k)) {
interp_face_reg(AMREX_D_DECL(i,j,k), rr, tag.fine, scomp, tag.crse,
tag.slope, ncomp, tag.domface, idim);
tag.slope, ncomp, tag.per_grown_domain, idim);
}
});
} else
Expand Down
152 changes: 150 additions & 2 deletions Src/AmrCore/AMReX_Interp_C.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,8 +15,8 @@ namespace amrex {

//
// Fill fine values with piecewise-constant interpolation of coarse data.
// Operate only on faces that overlap--ie, only fill the fine faces that make up each
// coarse face, leave the in-between faces alone.
// Operate only on faces that overlap -- i.e., only fill the fine faces that
// make up each coarse face, leave the in-between faces alone.
//
template<typename T>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void
Expand Down Expand Up @@ -84,6 +84,154 @@ face_linear_face_interp_z (int fi, int fj, int fk, int n, Array4<T> const& fine,
}
}

//
// Fill fine values with tangential interpolation of coarse data.
// Operate only on faces that overlap -- i.e., only fill the fine faces that
// make up each coarse face, leave the in-between faces alone.
//
template<typename T>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void
face_cons_linear_face_interp (int i, int j, int k, int n, Array4<T> const& fine,
Array4<T const> const& crse, Array4<int const> const& mask,
IntVect const& ratio, Box const& per_grown_domain, int dim) noexcept
{
int ci = amrex::coarsen(i, ratio[0]);

#if (AMREX_SPACEDIM == 1)
amrex::ignore_unused(per_grown_domain);
int cj = 0;
#else
int cj = amrex::coarsen(j, ratio[1]);
#endif

#if (AMREX_SPACEDIM == 3)
int ck = amrex::coarsen(k, ratio[2]);
#else
int ck = 0;
#endif

if (dim == 0 && ci*ratio[0] == i) {
// Check solve mask to ensure we don't overwrite valid fine data.
if (!mask || mask(ci, cj, ck, n)) {
fine(i, j, k, n) = crse(ci, cj, ck, n);
#if (AMREX_SPACEDIM >= 2)
if (cj != per_grown_domain.smallEnd(1) && cj != per_grown_domain.bigEnd(1) && ratio[1] > 1) {
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is the testing necessary? The CoarsenBox function grows the coarsened box by one cell. If the testing is necessary, it might be safer to test with > and < rather than !=.

Real sfy = Real(1.0);
Real dc = Real(0.5) * (crse(ci,cj+1,ck,n) - crse(ci,cj-1,ck,n));
Real df = Real(2.0) * (crse(ci,cj+1,ck,n) - crse(ci,cj ,ck,n));
Real db = Real(2.0) * (crse(ci,cj ,ck,n) - crse(ci,cj-1,ck,n));
Real sy = (df*db >= Real(0.0)) ?
amrex::min(std::abs(df),std::abs(db)) : Real(0.);
sy = std::copysign(Real(1.),dc)*amrex::min(sy,std::abs(dc));
if (dc != Real(0.0)) {
sfy = amrex::min(sfy, sy / dc);
}
Real slope = dc;
Real yoff = (static_cast<Real>(j - cj*ratio[1]) + Real(0.5)) / Real(ratio[1]) - Real(0.5);
fine(i,j,k,n) += yoff * slope * sfy;
} // jc
#if (AMREX_SPACEDIM == 3)
if (ck != per_grown_domain.smallEnd(2) && ck != per_grown_domain.bigEnd(2) && ratio[2] > 1) {
Real sfz = Real(1.0);
Real dc = Real(0.5) * (crse(ci,cj,ck+1,n) - crse(ci,cj,ck-1,n));
Real df = Real(2.0) * (crse(ci,cj,ck+1,n) - crse(ci,cj,ck ,n));
Real db = Real(2.0) * (crse(ci,cj,ck ,n) - crse(ci,cj,ck-1,n));
Real sz = (df*db >= Real(0.0)) ?
amrex::min(std::abs(df),std::abs(db)) : Real(0.);
sz = std::copysign(Real(1.),dc)*amrex::min(sz,std::abs(dc));
if (dc != Real(0.0)) {
sfz = amrex::min(sfz, sz / dc);
}
Real slope = dc;
Real zoff = (static_cast<Real>(k - ck*ratio[2]) + Real(0.5)) / Real(ratio[2]) - Real(0.5);
fine(i,j,k,n) += zoff * slope * sfz;
} // ck
#endif
#endif
} // mask
} // dim

#if (AMREX_SPACEDIM >= 2)
if (dim == 1 && cj*ratio[1] == j) {
// Check solve mask to ensure we don't overwrite valid fine data.
if (!mask || mask(ci, cj, ck, n)) {
fine(i, j, k, n) = crse(ci, cj, ck, n);
if (ci != per_grown_domain.smallEnd(0) && ci != per_grown_domain.bigEnd(0) && ratio[0] > 1) {
Real sfx = Real(1.0);
Real dc = Real(0.5) * (crse(ci+1,cj,ck,n) - crse(ci-1,cj,ck,n));
Real df = Real(2.0) * (crse(ci+1,cj,ck,n) - crse(ci ,cj,ck,n));
Real db = Real(2.0) * (crse(ci ,cj,ck,n) - crse(ci-1,cj,ck,n));
Real sx = (df*db >= Real(0.0)) ?
amrex::min(std::abs(df),std::abs(db)) : Real(0.);
sx = std::copysign(Real(1.),dc)*amrex::min(sx,std::abs(dc));
if (dc != Real(0.0)) {
sfx = amrex::min(sfx, sx / dc);
}
Real slope = dc;
Real xoff = (static_cast<Real>(i - ci*ratio[0]) + Real(0.5)) / Real(ratio[0]) - Real(0.5);
fine(i,j,k,n) += xoff * slope * sfx;
} // ci
#if (AMREX_SPACEDIM == 3)
if (ck != per_grown_domain.smallEnd(2) && ck != per_grown_domain.bigEnd(2) && ratio[2] > 1) {
Real sfz = Real(1.0);
Real dc = Real(0.5) * (crse(ci,cj,ck+1,n) - crse(ci,cj,ck-1,n));
Real df = Real(2.0) * (crse(ci,cj,ck+1,n) - crse(ci,cj,ck ,n));
Real db = Real(2.0) * (crse(ci,cj,ck ,n) - crse(ci,cj,ck-1,n));
Real sz = (df*db >= Real(0.0)) ?
amrex::min(std::abs(df),std::abs(db)) : Real(0.);
sz = std::copysign(Real(1.),dc)*amrex::min(sz,std::abs(dc));
if (dc != Real(0.0)) {
sfz = amrex::min(sfz, sz / dc);
}
Real slope = dc;
Real zoff = (static_cast<Real>(k - ck*ratio[2]) + Real(0.5)) / Real(ratio[2]) - Real(0.5);
fine(i,j,k,n) += zoff * slope * sfz;
} // ck
#endif // SPACEDIM >= 3
} // mask
} // dim == 1
#endif // SPACEDIM >= 2

#if (AMREX_SPACEDIM == 3)
if (dim == 2 && ck*ratio[2] == k) {
// Check solve mask to ensure we don't overwrite valid fine data.
if (!mask || mask(ci, cj, ck, n)) {
fine(i, j, k, n) = crse(ci, cj, ck, n);
if (ci != per_grown_domain.smallEnd(0) && ci != per_grown_domain.bigEnd(0) && ratio[0] > 1) {
Real sfx = Real(1.0);
Real dc = Real(0.5) * (crse(ci+1,cj,ck,n) - crse(ci-1,cj,ck,n));
Real df = Real(2.0) * (crse(ci+1,cj,ck,n) - crse(ci ,cj,ck,n));
Real db = Real(2.0) * (crse(ci ,cj,ck,n) - crse(ci-1,cj,ck,n));
Real sx = (df*db >= Real(0.0)) ?
amrex::min(std::abs(df),std::abs(db)) : Real(0.);
sx = std::copysign(Real(1.),dc)*amrex::min(sx,std::abs(dc));
if (dc != Real(0.0)) {
sfx = amrex::min(sfx, sx / dc);
}
Real slope = dc;
Real xoff = (static_cast<Real>(i - ci*ratio[0]) + Real(0.5)) / Real(ratio[0]) - Real(0.5);
fine(i,j,k,n) += xoff * slope * sfx;
} // ci
if (cj != per_grown_domain.smallEnd(1) && cj != per_grown_domain.bigEnd(1) && ratio[1] > 1) {
Real sfy = Real(1.0);
Real dc = Real(0.5) * (crse(ci,cj+1,ck,n) - crse(ci,cj-1,ck,n));
Real df = Real(2.0) * (crse(ci,cj+1,ck,n) - crse(ci,cj ,ck,n));
Real db = Real(2.0) * (crse(ci,cj ,ck,n) - crse(ci,cj-1,ck,n));
Real sy = (df*db >= Real(0.0)) ?
amrex::min(std::abs(df),std::abs(db)) : Real(0.);
sy = std::copysign(Real(1.),dc)*amrex::min(sy,std::abs(dc));
if (dc != Real(0.0)) {
sfy = amrex::min(sfy, sy / dc);
}
Real slope = dc;
Real yoff = (static_cast<Real>(j - cj*ratio[1]) + Real(0.5)) / Real(ratio[1]) - Real(0.5);
fine(i,j,k,n) += yoff * slope * sfy;
} // cj
} // mask
} // dim == 2
#endif
}

//
// Do linear in dir, pc transverse to dir, leave alone the fine values
// lining up with coarse edges--assume these have been set to hold the
Expand Down
Loading
Loading