Skip to content

Commit

Permalink
Merge pull request seahorce-scidac#239 from hklion/land_sea_mask
Browse files Browse the repository at this point in the history 
Land sea mask
  • Loading branch information
@hklion
hklion authored Aug 3, 2024
2 parents 6c8e4e3 + 6900b4c commit 39cafda
Show file tree
Hide file tree
Showing 30 changed files with 424 additions and 195 deletions.
35 changes: 19 additions & 16 deletions Source/BoundaryConditions/BoundaryConditions_cons.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,10 +13,13 @@ using namespace amrex;
// time is the time at which the data should be filled
// bccomp is the index into both domain_bcs_type_bcr and bc_extdir_vals for icomp = 0 --
// so this follows the BCVars enum
// n_not_fill perimeter of cells in x and y where BCs are not applied for conditions other than ext_dir.
// Foextrap is done based on the values at dom_lo-n_not_fill and dom_hi+n_not_fill. Reflecting conditions
// are untested.
//
void REMORAPhysBCFunct::impose_cons_bcs (const Array4<Real>& dest_arr, const Box& bx, const Box& domain,
const GpuArray<Real,AMREX_SPACEDIM> /*dxInv*/,
int icomp, int ncomp, Real /*time*/, int bccomp)
const GpuArray<Real,AMREX_SPACEDIM> /*dxInv*/, const Array4<const Real>& mskr,
int icomp, int ncomp, Real /*time*/, int bccomp, int n_not_fill)
{
BL_PROFILE_VAR("impose_cons_bcs()",impose_cons_bcs);
const auto& dom_lo = amrex::lbound(domain);
Expand Down Expand Up @@ -64,12 +67,12 @@ void REMORAPhysBCFunct::impose_cons_bcs (const Array4<Real>& dest_arr, const Box
ParallelFor(
bx_xlo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
if (bc_ptr[n].lo(0) == REMORABCType::ext_dir) {
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][0];
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][0] * mskr(i,j,0);
}
},
bx_xhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
if (bc_ptr[n].hi(0) == REMORABCType::ext_dir) {
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][3];
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][3] * mskr(i,j,0);
}
}
);
Expand All @@ -82,12 +85,12 @@ void REMORAPhysBCFunct::impose_cons_bcs (const Array4<Real>& dest_arr, const Box
ParallelFor(
bx_ylo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
if (bc_ptr[n].lo(1) == REMORABCType::ext_dir) {
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][1];
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][1] * mskr(i,j,0);
}
},
bx_yhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
if (bc_ptr[n].hi(1) == REMORABCType::ext_dir) {
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][4];
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][4] * mskr(i,j,0);
}
}
);
Expand All @@ -99,12 +102,12 @@ void REMORAPhysBCFunct::impose_cons_bcs (const Array4<Real>& dest_arr, const Box
ParallelFor(
bx_zlo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
if (bc_ptr[n].lo(2) == REMORABCType::ext_dir) {
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][2];
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][2] * mskr(i,j,0);
}
},
bx_zhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
if (bc_ptr[n].hi(2) == REMORABCType::ext_dir) {
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][5];
dest_arr(i,j,k,icomp+n) = l_bc_extdir_vals_d[n][5] * mskr(i,j,0);
}
}
);
Expand All @@ -115,16 +118,16 @@ void REMORAPhysBCFunct::impose_cons_bcs (const Array4<Real>& dest_arr, const Box
if (!is_periodic_in_x)
{
// Populate ghost cells on lo-x and hi-x domain boundaries
Box bx_xlo(bx); bx_xlo.setBig (0,dom_lo.x-1);
Box bx_xlo(bx); bx_xlo.setBig (0,dom_lo.x-1-n_not_fill);
bx_xlo.setSmall(2,std::max(dom_lo.z,bx.smallEnd(2)));
bx_xlo.setBig (2,std::min(dom_hi.z,bx.bigEnd(2)));
Box bx_xhi(bx); bx_xhi.setSmall(0,dom_hi.x+1);
Box bx_xhi(bx); bx_xhi.setSmall(0,dom_hi.x+1+n_not_fill);
bx_xhi.setSmall(2,std::max(dom_lo.z,bx.smallEnd(2)));
bx_xhi.setBig (2,std::min(dom_hi.z,bx.bigEnd(2)));
ParallelFor(bx_xlo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int iflip = dom_lo.x - 1 - i;
if (bc_ptr[n].lo(0) == REMORABCType::foextrap) {
dest_arr(i,j,k,icomp+n) = dest_arr(dom_lo.x,j,k,icomp+n);
dest_arr(i,j,k,icomp+n) = dest_arr(dom_lo.x-n_not_fill,j,k,icomp+n);
} else if (bc_ptr[n].lo(0) == REMORABCType::reflect_even) {
dest_arr(i,j,k,icomp+n) = dest_arr(iflip,j,k,icomp+n);
} else if (bc_ptr[n].lo(0) == REMORABCType::reflect_odd) {
Expand All @@ -134,7 +137,7 @@ void REMORAPhysBCFunct::impose_cons_bcs (const Array4<Real>& dest_arr, const Box
bx_xhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int iflip = 2*dom_hi.x + 1 - i;
if (bc_ptr[n].hi(0) == REMORABCType::foextrap) {
dest_arr(i,j,k,icomp+n) = dest_arr(dom_hi.x,j,k,icomp+n);
dest_arr(i,j,k,icomp+n) = dest_arr(dom_hi.x+n_not_fill,j,k,icomp+n);
} else if (bc_ptr[n].hi(0) == REMORABCType::reflect_even) {
dest_arr(i,j,k,icomp+n) = dest_arr(iflip,j,k,icomp+n);
} else if (bc_ptr[n].hi(0) == REMORABCType::reflect_odd) {
Expand All @@ -147,17 +150,17 @@ void REMORAPhysBCFunct::impose_cons_bcs (const Array4<Real>& dest_arr, const Box
if (!is_periodic_in_y)
{
// Populate ghost cells on lo-y and hi-y domain boundaries
Box bx_ylo(bx); bx_ylo.setBig (1,dom_lo.y-1);
Box bx_ylo(bx); bx_ylo.setBig (1,dom_lo.y-1-n_not_fill);
bx_ylo.setSmall(2,std::max(dom_lo.z,bx.smallEnd(2)));
bx_ylo.setBig (2,std::min(dom_hi.z,bx.bigEnd(2)));
Box bx_yhi(bx); bx_yhi.setSmall(1,dom_hi.y+1);
Box bx_yhi(bx); bx_yhi.setSmall(1,dom_hi.y+1+n_not_fill);
bx_yhi.setSmall(2,std::max(dom_lo.z,bx.smallEnd(2)));
bx_yhi.setBig (2,std::min(dom_hi.z,bx.bigEnd(2)));
ParallelFor(
bx_ylo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int jflip = dom_lo.y - 1 - j;
if (bc_ptr[n].lo(1) == REMORABCType::foextrap) {
dest_arr(i,j,k,icomp+n) = dest_arr(i,dom_lo.y,k,icomp+n);
dest_arr(i,j,k,icomp+n) = dest_arr(i,dom_lo.y-n_not_fill,k,icomp+n);
} else if (bc_ptr[n].lo(1) == REMORABCType::reflect_even) {
dest_arr(i,j,k,icomp+n) = dest_arr(i,jflip,k,icomp+n);
} else if (bc_ptr[n].lo(1) == REMORABCType::reflect_odd) {
Expand All @@ -167,7 +170,7 @@ void REMORAPhysBCFunct::impose_cons_bcs (const Array4<Real>& dest_arr, const Box
bx_yhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int jflip = 2*dom_hi.y + 1 - j;
if (bc_ptr[n].hi(1) == REMORABCType::foextrap) {
dest_arr(i,j,k,icomp+n) = dest_arr(i,dom_hi.y,k,icomp+n);
dest_arr(i,j,k,icomp+n) = dest_arr(i,dom_hi.y+n_not_fill,k,icomp+n);
} else if (bc_ptr[n].hi(1) == REMORABCType::reflect_even) {
dest_arr(i,j,k,icomp+n) = dest_arr(i,jflip,k,icomp+n);
} else if (bc_ptr[n].hi(1) == REMORABCType::reflect_odd) {
Expand Down
19 changes: 10 additions & 9 deletions Source/BoundaryConditions/BoundaryConditions_netcdf.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ using namespace amrex;
*/

void
REMORA::fill_from_bdyfiles (MultiFab& mf_to_fill, const Real time, const int bdy_var_type, const int icomp_to_fill)
REMORA::fill_from_bdyfiles (MultiFab& mf_to_fill, const MultiFab& mf_mask, const Real time, const int bdy_var_type, const int icomp_to_fill)
{
int lev = 0;

Expand Down Expand Up @@ -99,36 +99,37 @@ REMORA::fill_from_bdyfiles (MultiFab& mf_to_fill, const Real time, const int bdy
Box xhi_yhi = xhi & yhi;

const Array4<Real>& dest_arr = mf_to_fill.array(mfi);
const Array4<const Real>& mask_arr = mf_mask.array(mfi);

if (!xlo.isEmpty()) {
ParallelFor(xlo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k,icomp+icomp_to_fill) = oma * bdatxlo_n (ubound(xlo).x,j,k,0)
+ alpha * bdatxlo_np1(ubound(xlo).x,j,k,0);
dest_arr(i,j,k,icomp+icomp_to_fill) = (oma * bdatxlo_n (ubound(xlo).x,j,k,0)
+ alpha * bdatxlo_np1(ubound(xlo).x,j,k,0)) * mask_arr(i,j,0);
});
}

if (!xhi.isEmpty()) {
ParallelFor(xhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k,icomp+icomp_to_fill) = oma * bdatxhi_n (lbound(xhi).x,j,k,0)
+ alpha * bdatxhi_np1(lbound(xhi).x,j,k,0);
dest_arr(i,j,k,icomp+icomp_to_fill) = (oma * bdatxhi_n (lbound(xhi).x,j,k,0)
+ alpha * bdatxhi_np1(lbound(xhi).x,j,k,0)) * mask_arr(i,j,0);
});
}

if (!ylo.isEmpty()) {
ParallelFor(ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k,icomp+icomp_to_fill) = oma * bdatylo_n (i,ubound(ylo).y,k,0)
+ alpha * bdatylo_np1(i,ubound(ylo).y,k,0);
dest_arr(i,j,k,icomp+icomp_to_fill) = (oma * bdatylo_n (i,ubound(ylo).y,k,0)
+ alpha * bdatylo_np1(i,ubound(ylo).y,k,0)) * mask_arr(i,j,0);
});
}

if (!yhi.isEmpty()) {
ParallelFor(yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k,icomp+icomp_to_fill) = oma * bdatyhi_n (i,lbound(yhi).y,k,0)
+ alpha * bdatyhi_np1(i,lbound(yhi).y,k,0);
dest_arr(i,j,k,icomp+icomp_to_fill) = (oma * bdatyhi_n (i,lbound(yhi).y,k,0)
+ alpha * bdatyhi_np1(i,lbound(yhi).y,k,0)) * mask_arr(i,j,0);
});
}

Expand Down
18 changes: 9 additions & 9 deletions Source/BoundaryConditions/BoundaryConditions_xvel.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@ using namespace amrex;
// so this follows the BCVars enum
//
void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box& bx, const Box& domain,
const GpuArray<Real,AMREX_SPACEDIM> /*dxInv*/,
const GpuArray<Real,AMREX_SPACEDIM> /*dxInv*/, const Array4<const Real>& msku,
Real /*time*/, int bccomp)
{
BL_PROFILE_VAR("impose_xvel_bcs()",impose_xvel_bcs);
Expand Down Expand Up @@ -60,7 +60,7 @@ void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box
bx_xlo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int iflip = dom_lo.x - i;
if (bc_ptr[n].lo(0) == REMORABCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][0];
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][0]*msku(i,j,0);
} else if (bc_ptr[n].lo(0) == REMORABCType::foextrap) {
dest_arr(i,j,k) = dest_arr(dom_lo.x,j,k);
} else if (bc_ptr[n].lo(0) == REMORABCType::reflect_even) {
Expand All @@ -72,14 +72,14 @@ void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box
// We only set the values on the domain faces themselves if EXT_DIR
bx_xlo_face, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
if (bc_ptr[n].lo(0) == REMORABCType::ext_dir)
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][0];
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][0]*msku(i,j,0);
}
);
ParallelFor(
bx_xhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int iflip = 2*(dom_hi.x + 1) - i;
if (bc_ptr[n].hi(0) == REMORABCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][3];
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][3]*msku(i,j,0);
} else if (bc_ptr[n].hi(0) == REMORABCType::foextrap) {
dest_arr(i,j,k) = dest_arr(dom_hi.x+1,j,k);
} else if (bc_ptr[n].hi(0) == REMORABCType::reflect_even) {
Expand All @@ -91,7 +91,7 @@ void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box
// We only set the values on the domain faces themselves if EXT_DIR
bx_xhi_face, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
if (bc_ptr[n].lo(3) == REMORABCType::ext_dir)
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][3];
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][3]*msku(i,j,0);
}
);
} // not periodic in x
Expand All @@ -105,7 +105,7 @@ void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box
bx_ylo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int jflip = dom_lo.y - 1 - j;
if (bc_ptr[n].lo(1) == REMORABCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][1];
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][1]*msku(i,j,0);
} else if (bc_ptr[n].lo(1) == REMORABCType::foextrap) {
dest_arr(i,j,k) = dest_arr(i,dom_lo.y,k);
} else if (bc_ptr[n].lo(1) == REMORABCType::reflect_even) {
Expand All @@ -117,7 +117,7 @@ void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box
bx_yhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int jflip = 2*dom_hi.y + 1 - j;
if (bc_ptr[n].hi(1) == REMORABCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][4];
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][4]*msku(i,j,0);
} else if (bc_ptr[n].hi(1) == REMORABCType::foextrap) {
dest_arr(i,j,k) = dest_arr(i,dom_hi.y,k);
} else if (bc_ptr[n].hi(1) == REMORABCType::reflect_even) {
Expand All @@ -137,7 +137,7 @@ void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box
bx_zlo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int kflip = dom_lo.z - 1 - k;
if (bc_ptr[n].lo(2) == REMORABCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][2];
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][2]*msku(i,j,0);
} else if (bc_ptr[n].lo(2) == REMORABCType::foextrap) {
dest_arr(i,j,k) = dest_arr(i,j,dom_lo.z);
} else if (bc_ptr[n].lo(2) == REMORABCType::reflect_even) {
Expand All @@ -149,7 +149,7 @@ void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box
bx_zhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int kflip = 2*dom_hi.z + 1 - k;
if (bc_ptr[n].hi(2) == REMORABCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][5];
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][5]*msku(i,j,0);
} else if (bc_ptr[n].hi(2) == REMORABCType::foextrap) {
dest_arr(i,j,k) = dest_arr(i,j,dom_hi.z);
} else if (bc_ptr[n].hi(2) == REMORABCType::reflect_even) {
Expand Down
Loading

0 comments on commit 39cafda

Please sign in to comment.