Skip to content

Commit

Permalink
FaceLinear Interpolater for face data (AMReX-Codes#1118)
Browse files Browse the repository at this point in the history 
It's linear in the nodal direction and constant in others.  This preserves the divergence, but it's only first-order.
  • Loading branch information
@WeiqunZhang @dwillcox
WeiqunZhang authored and dwillcox committed Oct 3, 2020
1 parent 6bfdd15 commit 558b644
Show file tree
Hide file tree
Showing 5 changed files with 207 additions and 1 deletion.
10 changes: 10 additions & 0 deletions Src/AmrCore/AMReX_Interp_1D_C.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -296,6 +296,16 @@ nodebilin_interp (Box const& bx, Array4<T> const& fine, const int fcomp, const i
}
}

template<typename T>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void
face_linear_interp_x (int i, int /*j*/, int /*k*/, int n, Array4<T> const& fine,
Array4<T const> const& crse, IntVect const& ratio) noexcept
{
int ii = amrex::coarsen(i,ratio[0]);
Real const w = static_cast<Real>(i-ii*ratio[0]) * (1.0_rt/ratio[0]);
fine(i,0,0,n) = (1.0_rt-w) * crse(ii,0,0,n) + w * crse(ii+1,0,0,n);
}

}

#endif
22 changes: 22 additions & 0 deletions Src/AmrCore/AMReX_Interp_2D_C.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -428,6 +428,28 @@ nodebilin_interp (Box const& bx, Array4<T> const& fine, const int fcomp, const i
}
}

template<typename T>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void
face_linear_interp_x (int i, int j, int /*k*/, int n, Array4<T> const& fine,
Array4<T const> const& crse, IntVect const& ratio) noexcept
{
int ii = amrex::coarsen(i,ratio[0]);
int jj = amrex::coarsen(j,ratio[1]);
Real const w = static_cast<Real>(i-ii*ratio[0]) * (1.0_rt/ratio[0]);
fine(i,j,0,n) = (1.0_rt-w) * crse(ii,jj,0,n) + w * crse(ii+1,jj,0,n);
}

template<typename T>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void
face_linear_interp_y (int i, int j, int /*k*/, int n, Array4<T> const& fine,
Array4<T const> const& crse, IntVect const& ratio) noexcept
{
int ii = amrex::coarsen(i,ratio[0]);
int jj = amrex::coarsen(j,ratio[1]);
Real const w = static_cast<Real>(j-jj*ratio[1]) * (1.0_rt/ratio[1]);
fine(i,j,0,n) = (1.0_rt-w) * crse(ii,jj,0,n) + w * crse(ii,jj+1,0,n);
}

}

#endif
36 changes: 36 additions & 0 deletions Src/AmrCore/AMReX_Interp_3D_C.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -566,6 +566,42 @@ nodebilin_interp (Box const& bx, Array4<T> const& fine, const int fcomp, const i
}
}

template<typename T>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void
face_linear_interp_x (int i, int j, int k, int n, Array4<T> const& fine,
Array4<T const> const& crse, IntVect const& ratio) noexcept
{
int ii = amrex::coarsen(i,ratio[0]);
int jj = amrex::coarsen(j,ratio[1]);
int kk = amrex::coarsen(k,ratio[2]);
Real const w = static_cast<Real>(i-ii*ratio[0]) * (1.0_rt/ratio[0]);
fine(i,j,k,n) = (1.0_rt-w) * crse(ii,jj,kk,n) + w * crse(ii+1,jj,kk,n);
}

template<typename T>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void
face_linear_interp_y (int i, int j, int k, int n, Array4<T> const& fine,
Array4<T const> const& crse, IntVect const& ratio) noexcept
{
int ii = amrex::coarsen(i,ratio[0]);
int jj = amrex::coarsen(j,ratio[1]);
int kk = amrex::coarsen(k,ratio[2]);
Real const w = static_cast<Real>(j-jj*ratio[1]) * (1.0_rt/ratio[1]);
fine(i,j,k,n) = (1.0_rt-w) * crse(ii,jj,kk,n) + w * crse(ii,jj+1,kk,n);
}

template<typename T>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void
face_linear_interp_z (int i, int j, int k, int n, Array4<T> const& fine,
Array4<T const> const& crse, IntVect const& ratio) noexcept
{
int ii = amrex::coarsen(i,ratio[0]);
int jj = amrex::coarsen(j,ratio[1]);
int kk = amrex::coarsen(k,ratio[2]);
Real const w = static_cast<Real>(k-kk*ratio[2]) * (1.0_rt/ratio[2]);
fine(i,j,k,n) = (1.0_rt-w) * crse(ii,jj,kk,n) + w * crse(ii,jj,kk+1,n);
}

}

#endif
66 changes: 66 additions & 0 deletions Src/AmrCore/AMReX_Interpolater.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -646,10 +646,76 @@ public:
};
#endif

/**
* \brief Bilinear interpolation on face data.
*
* Bilinear interpolation on data.
*/

class FaceLinear
:
public Interpolater
{
public:

/**
* \brief The destructor.
*/
virtual ~FaceLinear () override;

/**
* \brief Returns coarsened box given fine box and refinement ratio.
*
* \param fine
* \param ratio
*/
virtual Box CoarseBox (const Box& fine,
int ratio) override;

/**
* \brief Returns coarsened box given fine box and refinement ratio.
*
* \param fine
* \param ratio
*/
virtual Box CoarseBox (const Box& fine,
const IntVect& ratio) override;

/**
* \brief Coarse to fine interpolation in space.
*
* \param crse
* \param crse_comp
* \param fine
* \param fine_comp
* \param ncomp
* \param fine_region
* \param ratio
* \param crse_geom
* \param fine_geom
* \param bcr
* \param actual_comp
* \param actual_state
*/
virtual void interp (const FArrayBox& crse,
int crse_comp,
FArrayBox& fine,
int fine_comp,
int ncomp,
const Box& fine_region,
const IntVect& ratio,
const Geometry& crse_geom,
const Geometry& fine_geom,
Vector<BCRec> const& bcr,
int actual_comp,
int actual_state,
RunOn gpu_or_cpu) override;
};

//! CONSTRUCT A GLOBAL OBJECT OF EACH VERSION.
extern PCInterp pc_interp;
extern NodeBilinear node_bilinear_interp;
extern FaceLinear face_linear_interp;
extern CellConservativeLinear lincc_interp;
extern CellConservativeLinear cell_cons_interp;

Expand Down
74 changes: 73 additions & 1 deletion Src/AmrCore/AMReX_Interpolater.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,8 @@
namespace amrex {

//
// PCInterp, NodeBilinear, and CellConservativeLinear are supported for all dimensions on cpu and gpu.
// PCInterp, NodeBilinear, FaceLinear, and CellConservativeLinear are supported for all dimensions
// on cpu and gpu.
//
// CellConsertiveProtected only works in 2D and 3D on cpu.
//
Expand All @@ -29,6 +30,7 @@ namespace amrex {
//
PCInterp pc_interp;
NodeBilinear node_bilinear_interp;
FaceLinear face_linear_interp;
CellConservativeLinear lincc_interp;
CellConservativeLinear cell_cons_interp(0);

Expand Down Expand Up @@ -141,6 +143,76 @@ NodeBilinear::interp (const FArrayBox& crse,
});
}

Box
FaceLinear::CoarseBox (const Box& fine, int ratio)
{
return CoarseBox(fine, IntVect(ratio));
}

Box
FaceLinear::CoarseBox (const Box& fine, const IntVect& ratio)
{
Box b = amrex::coarsen(fine,ratio);
for (int i = 0; i < AMREX_SPACEDIM; i++) {
if (b.type(i) == IndexType::NODE && b.length(i) < 2) {
// Don't want degenerate boxes in nodal direction.
b.growHi(i,1);
}
}
return b;
}

void
FaceLinear::interp (const FArrayBox& crse,
int crse_comp,
FArrayBox& fine,
int fine_comp,
int ncomp,
const Box& fine_region,
const IntVect& ratio,
const Geometry& /*crse_geom */,
const Geometry& /*fine_geom */,
Vector<BCRec> const& /*bcr*/,
int /*actual_comp*/,
int /*actual_state*/,
RunOn runon)
{
BL_PROFILE("FaceLinear::interp()");

AMREX_ASSERT(AMREX_D_TERM(fine_region.type(0),+fine_region.type(1),+fine_region.type(2)) == 1);

Array4<Real> const& fine_arr = fine.array(fine_comp);
Array4<Real const> const& crse_arr = crse.const_array(crse_comp);

if (fine_region.type(0) == IndexType::NODE)
{
AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(runon,fine_region,ncomp,i,j,k,n,
{
face_linear_interp_x(i,j,k,n,fine_arr,crse_arr,ratio);
});
}
#if (AMREX_SPACEDIM >= 2)
else if (fine_region.type(1) == IndexType::NODE)
{
AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(runon,fine_region,ncomp,i,j,k,n,
{
face_linear_interp_y(i,j,k,n,fine_arr,crse_arr,ratio);
});
}
#if (AMREX_SPACEDIM == 3)
else
{
AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(runon,fine_region,ncomp,i,j,k,n,
{
face_linear_interp_z(i,j,k,n,fine_arr,crse_arr,ratio);
});
}
#endif
#endif
}

FaceLinear::~FaceLinear () {}

#ifndef BL_NO_FORT
CellBilinear::~CellBilinear () {}

Expand Down

0 comments on commit 558b644

Please sign in to comment.