8#include <AMReX_buildInfo.H>
25 const DistributionMapping& dm)
27 BoxList bl2d = ba.boxList();
28 for (
auto& b : bl2d) {
31 BoxArray ba2d(std::move(bl2d));
33 amrex::Print() <<
"Making level " << lev <<
" from coarse" << std::endl;
34 amrex::Print() <<
"GRIDS AT LEVEL " << lev <<
" ARE " << ba << std::endl;
39 xvel_new[lev] =
new MultiFab(convert(ba, IntVect(1,0,0)), dm, 1,
xvel_new[lev-1]->nGrowVect());
40 xvel_old[lev] =
new MultiFab(convert(ba, IntVect(1,0,0)), dm, 1,
xvel_new[lev-1]->nGrowVect());
42 yvel_new[lev] =
new MultiFab(convert(ba, IntVect(0,1,0)), dm, 1,
yvel_new[lev-1]->nGrowVect());
43 yvel_old[lev] =
new MultiFab(convert(ba, IntVect(0,1,0)), dm, 1,
yvel_new[lev-1]->nGrowVect());
45 zvel_new[lev] =
new MultiFab(convert(ba, IntVect(0,0,1)), dm, 1,
zvel_new[lev-1]->nGrowVect());
46 zvel_old[lev] =
new MultiFab(convert(ba, IntVect(0,0,1)), dm, 1,
zvel_new[lev-1]->nGrowVect());
52 vec_ubar[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,3,IntVect(
NGROW,
NGROW,0)));
53 vec_vbar[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,3,IntVect(
NGROW,
NGROW,0)));
55 vec_ru[lev].reset(
new MultiFab(convert(ba,IntVect(1,0,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
56 vec_rv[lev].reset(
new MultiFab(convert(ba,IntVect(0,1,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
58 vec_ru2d[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
59 vec_rv2d[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
101 for (
int icomp=0; icomp<3; icomp++) {
107 for (
int icomp=0; icomp<2; icomp++) {
134#ifdef REMORA_USE_NETCDF
148#ifdef REMORA_USE_PARTICLES
165 BoxArray ba_old(
cons_new[lev]->boxArray());
166 DistributionMapping dm_old(
cons_new[lev]->DistributionMap());
168 BoxList bl2d = ba.boxList();
169 for (
auto& b : bl2d) {
172 BoxArray ba2d(std::move(bl2d));
174 amrex::Print() <<
"Remaking level " << lev << std::endl;
175 amrex::Print() <<
"GRIDS AT LEVEL " << lev <<
" ARE " << ba << std::endl;
191 MultiFab tmp_cons_new(ba, dm,
ncons, ngrow_state);
192 MultiFab tmp_cons_old(ba, dm,
ncons, ngrow_state);
194 MultiFab tmp_xvel_new(convert(ba, IntVect(1,0,0)), dm, 1, ngrow_vels);
195 MultiFab tmp_xvel_old(convert(ba, IntVect(1,0,0)), dm, 1, ngrow_vels);
197 MultiFab tmp_yvel_new(convert(ba, IntVect(0,1,0)), dm, 1, ngrow_vels);
198 MultiFab tmp_yvel_old(convert(ba, IntVect(0,1,0)), dm, 1, ngrow_vels);
200 MultiFab tmp_zvel_new(convert(ba, IntVect(0,0,1)), dm, 1, IntVect(ngrow_vels,ngrow_vels,0));
201 MultiFab tmp_zvel_old(convert(ba, IntVect(0,0,1)), dm, 1, IntVect(ngrow_vels,ngrow_vels,0));
203 MultiFab tmp_Zt_avg1_new(ba2d, dm, 1, IntVect(ngrow_zeta,ngrow_zeta,0));
204 MultiFab tmp_h(ba2d, dm, 2, IntVect(ngrow_h,ngrow_h,0));
206 MultiFab tmp_ubar_new(convert(ba2d, IntVect(1,0,0)), dm, 3, IntVect(ngrow_velbar,ngrow_velbar,0));
208 MultiFab tmp_vbar_new(convert(ba2d, IntVect(0,1,0)), dm, 3, IntVect(ngrow_velbar,ngrow_velbar,0));
210 MultiFab tmp_ru_new(convert(ba, IntVect(1,0,0)),dm,2,IntVect(
NGROW,
NGROW,0));
211 MultiFab tmp_rv_new(convert(ba, IntVect(0,1,0)),dm,2,IntVect(
NGROW,
NGROW,0));
213 MultiFab tmp_ru2d_new(convert(ba2d, IntVect(1,0,0)),dm,2,IntVect(
NGROW,
NGROW,0));
214 MultiFab tmp_rv2d_new(convert(ba2d, IntVect(0,1,0)),dm,2,IntVect(
NGROW,
NGROW,0));
218 tmp_cons_new.setVal(
zero);
219 tmp_xvel_new.setVal(
zero);
220 tmp_yvel_new.setVal(
zero);
221 tmp_zvel_new.setVal(
zero);
223 tmp_cons_old.setVal(
zero);
224 tmp_xvel_old.setVal(
zero);
225 tmp_yvel_old.setVal(
zero);
226 tmp_zvel_old.setVal(
zero);
228 tmp_ru_new.setVal(
zero);
229 tmp_rv_new.setVal(
zero);
231 tmp_ru2d_new.setVal(
zero);
232 tmp_rv2d_new.setVal(
zero);
234 tmp_ubar_new.setVal(
zero);
235 tmp_vbar_new.setVal(
zero);
240 FillPatch(lev,
time, tmp_xvel_new,
xvel_new,
xvel_bc(),
BdyVars::u,0,
true,
false,0,0,
zero,tmp_xvel_new);
241 FillPatch(lev,
time, tmp_yvel_new,
yvel_new,
yvel_bc(),
BdyVars::v,0,
true,
false,0,0,
zero,tmp_yvel_new);
245 for (
int icomp=0; icomp<3; icomp++) {
249 for (
int icomp=0; icomp<2; icomp++) {
257 MultiFab::Copy(tmp_cons_old,tmp_cons_new,0,0,
ncons,tmp_cons_new.nGrowVect());
258 MultiFab::Copy(tmp_xvel_old,tmp_xvel_new,0,0, 1,tmp_xvel_new.nGrowVect());
259 MultiFab::Copy(tmp_yvel_old,tmp_yvel_new,0,0, 1,tmp_yvel_new.nGrowVect());
260 MultiFab::Copy(tmp_zvel_old,tmp_zvel_new,0,0, 1,tmp_zvel_new.nGrowVect());
262 std::swap(tmp_cons_new, *
cons_new[lev]);
263 std::swap(tmp_cons_old, *
cons_old[lev]);
264 std::swap(tmp_xvel_new, *
xvel_new[lev]);
265 std::swap(tmp_xvel_old, *
xvel_old[lev]);
266 std::swap(tmp_yvel_new, *
yvel_new[lev]);
267 std::swap(tmp_yvel_old, *
yvel_old[lev]);
268 std::swap(tmp_zvel_new, *
zvel_new[lev]);
269 std::swap(tmp_zvel_old, *
zvel_old[lev]);
271 std::swap(tmp_ubar_new, *
vec_ubar[lev]);
272 std::swap(tmp_vbar_new, *
vec_vbar[lev]);
273 std::swap(tmp_ru_new, *
vec_ru[lev]);
274 std::swap(tmp_rv_new, *
vec_rv[lev]);
275 std::swap(tmp_ru2d_new, *
vec_ru2d[lev]);
276 std::swap(tmp_rv2d_new, *
vec_rv2d[lev]);
282 std::swap(tmp_h, *
vec_h[lev]);
306#ifdef REMORA_USE_NETCDF
314 bool ba_changed = (ba != ba_old);
315 bool dm_changed = (dm != dm_old);
316 if (ba_changed || dm_changed) {
321#ifdef REMORA_USE_PARTICLES
322 particleData.Redistribute();
339 const DistributionMapping& dm)
342 SetBoxArray(lev, ba);
343 SetDistributionMap(lev, dm);
345 BoxList bl2d = ba.boxList();
346 for (
auto& b : bl2d) {
349 BoxArray ba2d(std::move(bl2d));
351 amrex::Print() <<
"Making level " << lev <<
" from scratch" << std::endl;
352 amrex::Print() <<
"GRIDS AT LEVEL " << lev <<
" ARE " << ba << std::endl;
367 xvel_new[lev] =
new MultiFab(convert(ba, IntVect(1,0,0)), dm, 1, ngrow_vels);
368 xvel_old[lev] =
new MultiFab(convert(ba, IntVect(1,0,0)), dm, 1, ngrow_vels);
370 yvel_new[lev] =
new MultiFab(convert(ba, IntVect(0,1,0)), dm, 1, ngrow_vels);
371 yvel_old[lev] =
new MultiFab(convert(ba, IntVect(0,1,0)), dm, 1, ngrow_vels);
373 zvel_new[lev] =
new MultiFab(convert(ba, IntVect(0,0,1)), dm, 1, IntVect(ngrow_vels,ngrow_vels,0));
374 zvel_old[lev] =
new MultiFab(convert(ba, IntVect(0,0,1)), dm, 1, IntVect(ngrow_vels,ngrow_vels,0));
380 vec_ubar[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,3,IntVect(
NGROW,
NGROW,0)));
381 vec_vbar[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,3,IntVect(
NGROW,
NGROW,0)));
383 vec_ru[lev].reset(
new MultiFab(convert(ba,IntVect(1,0,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
384 vec_rv[lev].reset(
new MultiFab(convert(ba,IntVect(0,1,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
386 vec_ru2d[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
387 vec_rv2d[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
394#ifdef REMORA_USE_PARTICLES
399 particleData.Redistribute();
532 BoxList bl2d = ba.boxList();
533 for (
auto& b : bl2d) {
537 BoxArray ba2d(std::move(bl2d));
539 vec_msku[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,1,IntVect(
NGROW+1,
NGROW+1,0)));
540 vec_mskv[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,1,IntVect(
NGROW+1,
NGROW+1,0)));
541 vec_mskp[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,1,0)),dm,1,IntVect(
NGROW+1,
NGROW+1,0)));
566 BoxList bl2d = ba.boxList();
567 for (
auto& b : bl2d) {
570 BoxArray ba2d(std::move(bl2d));
572 BoxList bl1d = ba.boxList();
573 for (
auto& b : bl1d) {
577 BoxArray ba1d(std::move(bl1d));
580 ba_nd.surroundingNodes();
582 ba_w.surroundingNodes(2);
585 vec_z_w[lev].reset (
new MultiFab(convert(ba,IntVect(0,0,1)),dm,1,IntVect(
NGROW+1,
NGROW+1,0)));
589 vec_Huon[lev].reset (
new MultiFab(convert(ba,IntVect(1,0,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
590 vec_Hvom[lev].reset (
new MultiFab(convert(ba,IntVect(0,1,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
592 vec_Akv[lev].reset (
new MultiFab(convert(ba,IntVect(0,0,1)),dm,1,IntVect(
NGROW,
NGROW,0)));
601 vec_rufrc[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
602 vec_rvfrc[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,2,IntVect(
NGROW,
NGROW,0)));
604 vec_sustr[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
605 vec_svstr[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
619 vec_bustr[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
620 vec_bvstr[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
632 vec_rubar[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,4,IntVect(
NGROW,
NGROW,0)));
633 vec_rvbar[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,4,IntVect(
NGROW,
NGROW,0)));
646 vec_xu[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
647 vec_yu[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
648 vec_xv[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
649 vec_yv[lev].reset(
new MultiFab(convert(ba2d,IntVect(0,1,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
650 vec_xp[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,1,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
651 vec_yp[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,1,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
654 vec_dndx[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
655 vec_dmde[lev].reset(
new MultiFab(convert(ba2d,IntVect(1,0,0)),dm,1,IntVect(
NGROW,
NGROW,0)));
663 vec_bvf[lev].reset(
new MultiFab(convert(ba,IntVect(0,0,1)),dm,1,IntVect(
NGROW,
NGROW,0)));
665 vec_tke[lev].reset(
new MultiFab(convert(ba,IntVect(0,0,1)),dm,3,IntVect(
NGROW,
NGROW,0)));
666 vec_gls[lev].reset(
new MultiFab(convert(ba,IntVect(0,0,1)),dm,3,IntVect(
NGROW,
NGROW,0)));
668 vec_Akk[lev].reset(
new MultiFab(convert(ba,IntVect(0,0,1)),dm,1,IntVect(
NGROW,
NGROW,0)));
669 vec_Akp[lev].reset(
new MultiFab(convert(ba,IntVect(0,0,1)),dm,1,IntVect(
NGROW,
NGROW,0)));
798 const auto dxi = Geom(lev).InvCellSize();
799 vec_pm[lev]->setVal(
dxi[0]);
vec_pm[lev]->FillBoundary(geom[lev].periodicity());
800 vec_pn[lev]->setVal(
dxi[1]);
vec_pn[lev]->FillBoundary(geom[lev].periodicity());
803 vec_pm[lev]->FillBoundary(geom[lev].periodicity());
804 vec_pn[lev]->FillBoundary(geom[lev].periodicity());
807#ifdef REMORA_USE_NETCDF
812 Real dummy_time =
zero;
816 int rrx = ref_ratio[lev-1][0];
817 int rry = ref_ratio[lev-1][1];
819 for ( MFIter mfi(*
cons_new[lev], TilingIfNotGPU()); mfi.isValid(); ++mfi )
821 Array4<Real>
const& pm =
vec_pm[lev]->array(mfi);
822 Array4<Real>
const& pn =
vec_pn[lev]->array(mfi);
823 Box ubx = mfi.growntilebox(IntVect(
NGROW+1,
NGROW+2,0));
824 Box vbx = mfi.growntilebox(IntVect(
NGROW+2,
NGROW+1,0));
825 ParallelFor(makeSlab(ubx,2,0), [=] AMREX_GPU_DEVICE (
int i,
int j,
int ) {
826 pm(i,j,0) = pm(i,j,0) * (
rrx);
828 ParallelFor(makeSlab(vbx,2,0), [=] AMREX_GPU_DEVICE (
int i,
int j,
int ) {
829 pn(i,j,0) = pn(i,j,0) * (
rry);
849 for ( MFIter mfi(*
vec_xr[lev], TilingIfNotGPU()); mfi.isValid(); ++mfi )
851 Array4<const Real>
const& pm =
vec_pm[lev]->const_array(mfi);
852 Array4<const Real>
const& pn =
vec_pn[lev]->const_array(mfi);
853 Array4<Real>
const& xr =
vec_xr[lev]->array(mfi);
854 Array4<Real>
const& yr =
vec_yr[lev]->array(mfi);
855 Array4<Real>
const& xu =
vec_xu[lev]->array(mfi);
856 Array4<Real>
const& yu =
vec_yu[lev]->array(mfi);
857 Array4<Real>
const& xv =
vec_xv[lev]->array(mfi);
858 Array4<Real>
const& yv =
vec_yv[lev]->array(mfi);
859 Array4<Real>
const& xp =
vec_xp[lev]->array(mfi);
860 Array4<Real>
const& yp =
vec_yp[lev]->array(mfi);
862 Box bx = mfi.growntilebox(IntVect(
NGROW,
NGROW,0));
863 ParallelFor(bx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int)
865 xr(i,j,0) = (i + Real(0.5)) / pm(i,j,0);
866 yr(i,j,0) = (j + Real(0.5)) / pn(i,j,0);
869 ParallelFor(grow(convert(bx,IntVect(1,0,0)),IntVect(-1,0,0)), [=] AMREX_GPU_DEVICE (
int i,
int j,
int)
871 xu(i,j,0) = i / pm(i,j,0);
872 yu(i,j,0) = (j + Real(0.5)) / pn(i,j,0);
875 ParallelFor(grow(convert(bx,IntVect(0,1,0)),IntVect(0,-1,0)), [=] AMREX_GPU_DEVICE (
int i,
int j,
int)
877 xv(i,j,0) = (i + Real(0.5)) / pm(i,j,0);
878 yv(i,j,0) = j / pn(i,j,0);
881 ParallelFor(grow(convert(bx,IntVect(1,1,0)),IntVect(-1,-1,0)), [=] AMREX_GPU_DEVICE (
int i,
int j,
int)
883 xp(i,j,0) = i / pm(i,j,0);
884 yp(i,j,0) = j / pn(i,j,0);
894 for ( MFIter mfi(*
vec_dndx[lev], TilingIfNotGPU()); mfi.isValid(); ++mfi )
896 Array4<const Real>
const& pm =
vec_pm[lev]->const_array(mfi);
897 Array4<const Real>
const& pn =
vec_pn[lev]->const_array(mfi);
898 Array4<Real>
const& dndx =
vec_dndx[lev]->array(mfi);
899 Array4<Real>
const& dmde =
vec_dmde[lev]->array(mfi);
901 Box bx = mfi.growntilebox(IntVect(
NGROW,
NGROW,0));
902 ParallelFor(bx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int)
904 dndx(i,j,0) = Real(0.5) * (
one / pn(i+1,j ,0) -
one / pn(i-1,j ,0));
905 dmde(i,j,0) = Real(0.5) * (
one / pm(i ,j+1,0) -
one / pn(i ,j-1,0));
916 BL_PROFILE(
"REMORA::set_zeta_to_Ztavg()");
917 std::unique_ptr<MultiFab>& mf_zeta =
vec_zeta[lev];
918 std::unique_ptr<MultiFab>& mf_Zt_avg1 =
vec_Zt_avg1[lev];
920 for ( MFIter mfi(*
vec_zeta[lev], TilingIfNotGPU()); mfi.isValid(); ++mfi )
922 Array4<Real>
const& Zt_avg1 = (mf_Zt_avg1)->array(mfi);
923 Array4<const Real>
const& evap =
vec_evap[lev]->const_array(mfi);
924 Array4<const Real>
const& rain =
vec_rain[lev]->const_array(mfi);
925 Array4<const Real>
const& EminusP =
vec_EminusP[lev]->const_array(mfi);
929 Box bx2 = mfi.growntilebox(IntVect(
NGROW,
NGROW,0));
931 Real cff =
dt[lev] /
rhow;
932 Real dt_lev =
dt[lev];
934 ParallelFor(bx2, [=] AMREX_GPU_DEVICE (
int i,
int j,
int )
936 if (use_EminusP_from_input) {
938 Zt_avg1(i,j,0) = Zt_avg1(i,j,0) - EminusP(i,j,0) * dt_lev;
940 Zt_avg1(i,j,0) = Zt_avg1(i,j,0) - (evap(i,j,0) - rain(i,j,0)) * cff;
945 Gpu::streamSynchronize();
947 vec_Zt_avg1[lev]->FillBoundary(geom[lev].periodicity());
949 for ( MFIter mfi(*
vec_zeta[lev], TilingIfNotGPU()); mfi.isValid(); ++mfi )
951 Box bx3 = mfi.tilebox(); bx3.grow(IntVect(
NGROW+1,
NGROW+1,0));
952 Array4<Real>
const& zeta = mf_zeta->array(mfi);
953 Array4<Real>
const& Zt_avg1 = (mf_Zt_avg1)->array(mfi);
955 ParallelFor(bx3, 3, [=] AMREX_GPU_DEVICE (
int i,
int j,
int ,
int n)
957 zeta(i,j,0,n) = Zt_avg1(i,j,0);
968 for ( MFIter mfi(*
vec_mskr[lev], TilingIfNotGPU()); mfi.isValid(); ++mfi )
970 Array4<const Real>
const& mskr =
vec_mskr[lev]->const_array(mfi);
971 Array4< Real>
const& mskp =
vec_mskp[lev]->array(mfi);
973 Box bx = mfi.tilebox(); bx.grow(IntVect(1,1,0)); bx.makeSlab(2,0);
978 ParallelFor(bx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int)
980 if ((mskr(i-1,j,0) > Real(0.5)) and (mskr(i,j,0) > Real(0.5)) and (mskr(i-1,j-1,0) > Real(0.5)) and (mskr(i,j-1,0) > Real(0.5))) {
982 }
else if ((mskr(i-1,j,0) < Real(0.5)) and (mskr(i,j,0) > Real(0.5)) and (mskr(i-1,j-1,0) > Real(0.5)) and (mskr(i,j-1,0) > Real(0.5))) {
984 }
else if ((mskr(i-1,j,0) > Real(0.5)) and (mskr(i,j,0) < Real(0.5)) and (mskr(i-1,j-1,0) > Real(0.5)) and (mskr(i,j-1,0) > Real(0.5))) {
986 }
else if ((mskr(i-1,j,0) > Real(0.5)) and (mskr(i,j,0) > Real(0.5)) and (mskr(i-1,j-1,0) < Real(0.5)) and (mskr(i,j-1,0) > Real(0.5))) {
988 }
else if ((mskr(i-1,j,0) > Real(0.5)) and (mskr(i,j,0) > Real(0.5)) and (mskr(i-1,j-1,0) > Real(0.5)) and (mskr(i,j-1,0) < Real(0.5))) {
990 }
else if ((mskr(i-1,j,0) > Real(0.5)) and (mskr(i,j,0) < Real(0.5)) and (mskr(i-1,j-1,0) > Real(0.5)) and (mskr(i,j-1,0) < Real(0.5))) {
992 }
else if ((mskr(i-1,j,0) < Real(0.5)) and (mskr(i,j,0) > Real(0.5)) and (mskr(i-1,j-1,0) < Real(0.5)) and (mskr(i,j-1,0) > Real(0.5))) {
994 }
else if ((mskr(i-1,j,0) > Real(0.5)) and (mskr(i,j,0) > Real(0.5)) and (mskr(i-1,j-1,0) < Real(0.5)) and (mskr(i,j-1,0) < Real(0.5))) {
996 }
else if ((mskr(i-1,j,0) < Real(0.5)) and (mskr(i,j,0) < Real(0.5)) and (mskr(i-1,j-1,0) > Real(0.5)) and (mskr(i,j-1,0) > Real(0.5))) {
1012 for ( MFIter mfi(*
vec_mskr[lev], TilingIfNotGPU()); mfi.isValid(); ++mfi )
1014 Array4<const Real>
const& mskr =
vec_mskr[lev]->const_array(mfi);
1015 Array4< Real>
const& msku =
vec_msku[lev]->array(mfi);
1016 Array4< Real>
const& mskv =
vec_mskv[lev]->array(mfi);
1017 Array4< Real>
const& mskp =
vec_mskp[lev]->array(mfi);
1019 Box bx = mfi.tilebox(); bx.grow(IntVect(1,1,0)); bx.makeSlab(2,0);
1021 ParallelFor(bx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int)
1023 msku(i,j,0) = mskr(i-1,j ,0) * mskr(i,j,0);
1024 mskv(i,j,0) = mskr(i ,j-1,0) * mskr(i,j,0);
1025 mskp(i,j,0) = mskr(i-1,j-1,0) * mskr(i,j,0) * mskr(i-1,j,0) * mskr(i,j-1,0);
1036 for ( MFIter mfi(*
vec_mskr3d[lev], TilingIfNotGPU()); mfi.isValid(); ++mfi )
1038 Array4<const Real>
const& mskr =
vec_mskr[lev]->const_array(mfi);
1039 Array4< Real>
const& mskr3d =
vec_mskr3d[lev]->array(mfi);
1041 Box bx = mfi.tilebox(); bx.grow(IntVect(1,1,0));
1043 ParallelFor(bx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k)
1045 mskr3d(i,j,k) = mskr(i,j,0);
constexpr amrex::Real two
constexpr amrex::Real bogus_large_value
constexpr amrex::Real one
constexpr amrex::Real zero
constexpr amrex::Real rhow
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_EminusP
evaporation minus precipitation [kg/m^2/s], defined at rho-points
void set_grid_vars_averaged_down(int lev)
Set pm/pn by averaging down from higher-resolution grid.
int ncons
Number of conserved scalars in the state (temperature + salt + passive scalars)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_zeta_full_domain
high resolution initial free surface height (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rv2d
v velocity RHS (2D, includes horizontal and vertical advection)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_evap
evaporation rate [kg/m^2/s]
int zvel_bc() const noexcept
int xvel_bc() const noexcept
amrex::Vector< amrex::BCRec > domain_bcs_type
vector (over BCVars) of BCRecs
void calculate_nodal_masks(int lev)
Calculate u-, v-, and psi-point masks based on rho-point masks after analytic initialization.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_fcor
coriolis factor (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_btflux
Bottom tracer flux; input arrays.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rubar
barotropic x velocity for the RHS (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xvel_full_domain
multilevel data container for high res initial x velocities (u in ROMS)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_h
multilevel data container for current step's z velocities (largely unused; W stored separately)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_pm
horizontal scaling factor: 1 / dx (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ZoBot
Bottom roughness length [m], defined at rho points.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_DU_avg2
correct time average of barotropic x velocity flux for coupling (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_lrflx
longwave radiation
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yv
y_grid on v-points (2D)
amrex::Vector< amrex::MultiFab * > cons_new
multilevel data container for current step's scalar data: temperature, salinity, passive tracer
virtual void MakeNewLevelFromCoarse(int lev, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm) override
Make a new level using provided BoxArray and DistributionMapping and fill with interpolated coarse le...
void stretch_transform(int lev)
Calculate vertical stretched coordinates.
void set_curvilinear_terms_from_grid_scale(int lev)
Set curvilinear derivative terms on level lev based on pm and pn.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_vwind
Wind in the v direction, defined at rho-points.
std::unique_ptr< ProblemBase > prob
Pointer to container of analytical functions for problem definition.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskr
land/sea mask at cell centers (2D)
void Construct_REMORAFillPatchers(int lev)
Construct FillPatchers.
void init_grid_vars_from_netcdf(int lev)
Grid variable initialization from NetCDF file.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rain
precipitation rate [kg/m^2/s]
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_tke
Turbulent kinetic energy.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_stflx
Surface tracer flux; working arrays.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_gls
Turbulent generic length scale.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_sustr
Surface stress in the u direction.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yp
y_grid on psi-points (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xr
x_grid on rho points (2D)
int yvel_bc() const noexcept
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ru2d
u velocity RHS (2D, includes horizontal and vertical advection)
virtual void ClearLevel(int lev) override
Delete level data Overrides the pure virtual function in AmrCore.
amrex::Vector< amrex::MultiFab * > zvel_new
multilevel data container for current step's z velocities (largely unused; W stored separately)
AMREX_FORCE_INLINE int ComputeGhostCells(const int &spatial_order)
Helper function to determine number of ghost cells.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_sstore
additional scratch space for calculations on temp, salt, etc
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xv
x_grid on v-points (2D)
void init_only(int lev, amrex::Real time)
Init (NOT restart or regrid)
void init_set_vmix(int lev)
Initialize vertical mixing coefficients from file or analytic.
void set_grid_scale(int lev)
Set pm and pn arrays and x/y coords on level lev.
void set_coriolis(int lev)
Initialize Coriolis factor from file or analytic.
int foextrap_bc() const noexcept
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Lscale
Vertical mixing turbulent length scale.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Hz
Width of cells in the vertical (z-) direction (3D, Hz in ROMS)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akt
Vertical diffusion coefficient (3D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_msku
land/sea mask at x-faces (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rvfrc
v velocity RHS, integrated, including advection and bottom/surface stresses (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_pm_full_domain
horizontal scaling factor: 1 / dx (2D) on whole domain
amrex::Vector< amrex::MultiFab * > xvel_old
multilevel data container for last step's x velocities (u in ROMS)
amrex::Vector< amrex::MultiFab * > yvel_new
multilevel data container for current step's y velocities (v in ROMS)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_uwind
Wind in the u direction, defined at rho-points.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rufrc
u velocity RHS, integrated, including advection and bottom/surface stresses (2D)
int zeta_bc() const noexcept
void Define_REMORAFillPatchers(int lev)
Define FillPatchers.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_shflx
sensible heat flux
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_visc2_p
Harmonic viscosity defined on the psi points (corners of horizontal grid cells)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_dmde
d(1/m)/d(eta)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rvbar
barotropic y velocity for the RHS (2D)
amrex::Vector< amrex::MultiFab * > zvel_old
multilevel data container for last step's z velocities (largely unused; W stored separately)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_z_r
z coordinates at rho points (cell centers)
amrex::Vector< amrex::MultiFab * > xvel_new
multilevel data container for current step's x velocities (u in ROMS)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_lhflx
latent heat flux
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskp
land/sea mask at cell corners (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_bvf
Brunt-Vaisala frequency (3D)
virtual void RemakeLevel(int lev, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm) override
Remake an existing level using provided BoxArray and DistributionMapping and fill with existing fine ...
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_lonp
longitude on psi-points (2D, degrees east); only filled when the grid NetCDF file carries lon_psi
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskv
land/sea mask at y-faces (2D)
void init_masks(int lev, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm)
Allocate MultiFabs for masks.
amrex::Vector< std::unique_ptr< REMORAPhysBCFunct > > physbcs
Vector (over level) of functors to apply physical boundary conditions.
void fill_3d_masks(int lev)
Copy maskr to all z levels.
void FillCoarsePatch(int lev, amrex::Real time, amrex::MultiFab *mf_fine, amrex::MultiFab *mf_crse, const int bccomp, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const int n_not_fill=0, const int icomp_calc=0, const amrex::Real dt=zero, const amrex::MultiFab &mf_calc=amrex::MultiFab())
fill an entire multifab by interpolating from the coarser level
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rhoS
density perturbation
void init_riv_pos_from_netcdf(int lev)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_visc2_r
Harmonic viscosity defined on the rho points (centers)
void update_mskp(int lev)
Set psi-point mask to be consistent with rho-point mask.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yvel_full_domain
multilevel data container for high res initial y velocities (v in ROMS)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_svstr
Surface stress in the v direction.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Huon
u-volume flux (3D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_h_full_domain
Bathymetry data on the whole domain at each potential level.
amrex::Vector< amrex::MultiFab * > yvel_old
multilevel data container for last step's y velocities (v in ROMS)
void init_stuff(int lev, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm)
Allocate MultiFabs for state and evolution variables.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rhoA
vertically-averaged density
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_DV_avg1
time average of barotropic y velocity flux
int hires_init_level
Which level the high resolution initialization data is at.
amrex::Vector< amrex::Real > t_new
new time at each level
static SolverChoice solverChoice
Container for algorithmic choices.
void set_grid_coords_from_grid_scale(int lev)
Set x/y coords on level lev based on pm and pn.
void resize_stuff(int lev)
Resize variable containers to accommodate data on levels 0 to max_lev.
amrex::Vector< std::unique_ptr< amrex::iMultiFab > > vec_river_position
iMultiFab for river positions; contents are indices of rivers
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akk
Turbulent kinetic energy vertical diffusion coefficient.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rdrag2
Quadratic drag coefficient [unitless], defined at rho points.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ru
u velocity RHS (3D, includes horizontal and vertical advection)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xp
x_grid on psi-points (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_longwave_down
Downward longwave radiation.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_zeta
free surface height (2D)
int ubar_bc() const noexcept
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_vbar
barotropic y velocity (2D)
void FillCoarsePatchPC(int lev, amrex::Real time, amrex::MultiFab *mf_fine, amrex::MultiFab *mf_crse, const int bccomp, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const int n_not_fill=0, const int icomp_calc=0, const amrex::Real dt=zero, const amrex::MultiFab &mf_calc=amrex::MultiFab())
fill an entire multifab by interpolating from the coarser level using the piecewise constant interpol...
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_DU_avg1
time average of barotropic x velocity flux (2D)
void set_zeta_to_Ztavg(int lev)
Set zeta components to be equal to time-averaged Zt_avg1.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_alpha
Thermal expansion coefficient (3D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ubar
barotropic x velocity (2D)
amrex::Vector< amrex::MultiFab * > cons_old
multilevel data container for last step's scalar data: temperature, salinity, passive tracer
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_bustr
Bottom stress in the u direction.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_cons_full_domain
multilevel data container for high res initial data: temperature, salinity, passive tracer
void FillPatch(int lev, amrex::Real time, amrex::MultiFab &mf_to_be_filled, amrex::Vector< amrex::MultiFab * > const &mfs, const int bccomp, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const bool fill_set=false, const int n_not_fill=0, const int icomp_calc=0, const amrex::Real dt=zero, const amrex::MultiFab &mf_calc=amrex::MultiFab())
Fill a new MultiFab by copying in phi from valid region and filling ghost cells.
int vbar_bc() const noexcept
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_DV_avg2
correct time average of barotropic y velocity flux for coupling (2D)
void set_hmixcoef(int lev)
Initialize horizontal mixing coefficients.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yu
y_grid on u-points (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_dndx
d(1/n)/d(xi)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_bvstr
Bottom stress in the v direction.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rzeta
free surface height for the RHS (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_z_phys_nd
z coordinates at psi points (cell nodes)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xu
x_grid on u-points (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_pn_full_domain
horizontal scaling factor: 1 / dy (2D) on whole domain
amrex::Vector< amrex::Vector< std::unique_ptr< amrex::MultiFab > > > vec_nudg_coeff
Climatology nudging coefficients.
void set_weights(int lev)
Set weights for averaging 3D variables to 2D.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_pn
horizontal scaling factor: 1 / dy (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akv
Vertical viscosity coefficient (3D)
virtual void MakeNewLevelFromScratch(int lev, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm) override
Make a new level from scratch using provided BoxArray and DistributionMapping. Only used during initi...
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_stflux
Surface tracer flux; input arrays.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rdrag
Linear drag coefficient [m/s], defined at rho points.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_cloud
cloud cover fraction [0-1], defined at rho-points
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Zt_avg1
Average of the free surface, zeta (2D)
std::string restart_chkfile
If set, restart from this checkpoint file.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rv
v velocity RHS (3D, includes horizontal and vertical advection)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_btflx
Bottom tracer flux; working arrays.
int cf_width
Nudging width at coarse-fine interface.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_beta
Saline contraction coefficient (3D)
amrex::Vector< amrex::Real > t_old
old time at each level
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskr3d
land/sea mask at cell centers, copied to all z levels (3D)
amrex::Vector< amrex::GpuArray< amrex::Real, AMREX_SPACEDIM *2 > > m_bc_extdir_vals
Array holding the Dirichlet values at walls which need them.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_srflx
Shortwave radiation flux [W/m²], defined at rho-points.
void set_bathymetry_averaged_down(int lev)
Copy over bathymetry data that has been averaged down from high resolution input netcdf file.
amrex::Vector< amrex::Real > dt
time step at each level
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Pair
Air pressure [mb], defined at rho-points.
amrex::Gpu::DeviceVector< amrex::BCRec > domain_bcs_type_d
GPU vector (over BCVars) of BCRecs.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_qair
Specific humidity [kg/kg], defined at rho-points.
int hires_grid_level
Which level the high resolution bathymetry is at.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_z_w
z coordinates at w points (faces between z-cells)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Tair
Air temperature [°C], defined at rho-points.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yr
y_grid on rho points (2D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Hvom
v-volume flux (3D)
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akp
Turbulent length scale vertical diffusion coefficient.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_diff2
Harmonic diffusivity for temperature / salinity.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_latp
latitude on psi-points (2D, degrees north); only filled when the grid NetCDF file carries lat_psi
static constexpr int cons_bc
static constexpr int Temp_bc_comp
@ EminusP
evaporation minus precipitation [m/s]
bool use_curvilinear_grid
BottomStressType bottom_stress_type
VertMixingType vert_mixing_type
std::array< BulkForcingType, BulkFlux::NumTypes > bulk_flux_type
GridScaleType grid_scale_type