REMORA
Regional Modeling of Oceans Refined Adaptively
Loading...
Searching...
No Matches
REMORA_FillPatch.cpp
Go to the documentation of this file.
1#include <REMORA.H>
5
6using namespace amrex;
7
8PhysBCFunctNoOp null_bc;
9
10/**
11 * Fill valid and ghost data in the MultiFab "mf"
12 * This version fills the MultiFab mf_to_fill in valid regions with
13 * the "state data" at the given time;
14 * values in mf when it is passed in are *not* used.
15 *
16 * @param[in ] lev level to FillPatch on
17 * @param[in ] time current time
18 * @param[inout] mf_to_fill MultiFab to fill
19 * @param[in ] mfs vector over levels of multifabs associated with mf_to_fill
20 * @param[in ] bccomp index into domain_bcs_type
21 * @param[in ] bdy_var_type when filling from NetCDF, which boundary data
22 * @param[in ] icomp component to fill
23 * @param[in ] fill_all whether to fill all components
24 * @param[in ] fill_set whether to use FillPatchers
25 * @param[in ] icomp_calc component to use in RHS boundary calculation
26 * @param[in ] dt_lev time step at this level
27 * @param[in ] mf_calc data to use in RHS boundary calculation
28 */
29void
30REMORA::FillPatch (int lev, Real time, MultiFab& mf_to_fill, Vector<MultiFab*> const& mfs,
31 const int bccomp,
32#ifdef REMORA_USE_NETCDF
33 const int bdy_var_type,
34#else
35 const int /*bdy_var_type*/,
36#endif
37 const int icomp,
38 const bool fill_all,
39 const bool fill_set,
40 const int n_not_fill,
41#ifdef REMORA_USE_NETCDF
42 const int icomp_calc,
43 const Real dt_lev,
44#else
45 const int /*icomp_calc*/,
46 const Real /*dt_lev*/,
47#endif
48 const MultiFab& mf_calc)
49{
50 BL_PROFILE_VAR("REMORA::FillPatch()",REMORA_FillPatch);
51 amrex::Interpolater* mapper = nullptr;
52
53 Box mf_box(mf_to_fill.boxArray()[0]);
54 bool is_2d = mf_box.length(2) == 1;
55
56 MultiFab* mask;
57
58 //
59 // ***************************************************************************
60 // The first thing we do is interpolate the momenta on the "valid" faces of
61 // the fine grids (where the interface is coarse/fine not fine/fine) -- this
62 // will not be over-written below because the FillPatch operators see these as
63 // valid faces.
64 // ***************************************************************************
65 if (lev>0 && fill_set) {
66 if (cf_set_width > 0 &&
67 mf_box.ixType() == IndexType(IntVect(0,0,0))) {
68 FPr_c[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
69 } else if (cf_set_width >= 0) {
70 if (!is_2d) {
71 if (mf_box.ixType() == IndexType(IntVect(1,0,0))) {
72 FPr_u[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
73 } else if (mf_box.ixType() == IndexType(IntVect(0,1,0))) {
74 FPr_v[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
75 } else if (mf_box.ixType() == IndexType(IntVect(0,0,1))) {
76 FPr_w[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
77 }
78 } else {
79 if (mf_box.ixType() == IndexType(IntVect(1,0,0))) {
80 FPr_ubar[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
81 } else if (mf_box.ixType() == IndexType(IntVect(0,1,0))) {
82 FPr_vbar[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
83 }
84 }
85 }
86 }
87
88 int ncomp;
89 if (fill_all) {
90 ncomp = mf_to_fill.nComp();
91 } else {
92 ncomp = 1;
93 }
94
95 if (mf_box.ixType() == IndexType(IntVect(0,0,0)))
96 {
97 mapper = &cell_cons_interp;
98 mask = vec_mskr[lev].get();
99 }
100 else if (mf_box.ixType() == IndexType(IntVect(1,0,0)))
101 {
102 mapper = &face_linear_interp;
103 mask = vec_msku[lev].get();
104 }
105 else if (mf_box.ixType() == IndexType(IntVect(0,1,0)))
106 {
107 mapper = &face_linear_interp;
108 mask = vec_mskv[lev].get();
109 }
110 else {
111 mapper = &face_linear_interp;
112 mask = vec_mskr[lev].get();
113 }
114
115 if (lev == 0)
116 {
117 Vector<MultiFab*> fmf = {mfs[lev], mfs[lev]};
118 Vector<Real> ftime = {t_old[lev], t_new[lev]};
119 amrex::FillPatchSingleLevel(mf_to_fill, time, fmf, ftime, icomp, icomp, ncomp,
120 geom[lev], null_bc, bccomp);
121 }
122 else
123 {
124 Vector<MultiFab*> fmf = {mfs[lev], mfs[lev]};
125 Vector<Real> ftime = {t_old[lev], t_new[lev]};
126 Vector<MultiFab*> cmf = {mfs[lev-1], mfs[lev-1]};
127 Vector<Real> ctime = {t_old[lev-1], t_new[lev-1]};
128
129 mfs[lev-1]->FillBoundary(geom[lev-1].periodicity());
130 amrex::FillPatchTwoLevels(mf_to_fill, mf_to_fill.nGrowVect(), IntVect(0,0,0),
131 time, cmf, ctime, fmf, ftime,
132 icomp, icomp, ncomp, geom[lev-1], geom[lev],
133 //null_bc, bccomp, null_bc, bccomp,
134 refRatio(lev-1),
135 mapper, domain_bcs_type, bccomp);
136 } // lev > 0
137
138 // This is currently unconditionally true, but do_bc flag should actually be passed in
139 const bool do_bc=true;
140 if (do_bc) {
141 // ***************************************************************************
142 // Physical bc's at domain boundary
143 // ***************************************************************************
144
145 // Enforce physical boundary conditions
146 (*physbcs[lev])(mf_to_fill,*mask,icomp,ncomp,mf_to_fill.nGrowVect(),time,bccomp,n_not_fill,mf_calc,
147 *vec_msku[lev].get(), *vec_mskv[lev].get());
148
149#ifdef REMORA_USE_NETCDF
150 // Fill the data which is stored in the boundary data read from netcdf files
151 if ( (solverChoice.boundary_from_netcdf) && (lev==0) &&
152 (bdy_var_type != BdyVars::null) )
153 {
154 fill_from_bdyfiles(mf_to_fill,*mask,time,bccomp,bdy_var_type, icomp,icomp_calc,mf_calc,dt_lev);
155 }
156#endif
157 // Fill corners of the domain with periodic data
158 if ( mf_box.ixType() == IndexType(IntVect(0,0,0)) ) {
159 mf_to_fill.EnforcePeriodicity(geom[lev].periodicity());
160 }
161
162 // Also enforce free-slip at top boundary (on xvel or yvel)
163 if ( (mf_box.ixType() == IndexType(IntVect(1,0,0))) ||
164 (mf_box.ixType() == IndexType(IntVect(0,1,0))) )
165 {
166 int khi = geom[lev].Domain().bigEnd(2);
167 for (MFIter mfi(mf_to_fill); mfi.isValid(); ++mfi)
168 {
169 Box gbx = mfi.growntilebox(); // Note this is face-centered since vel is
170 gbx.setSmall(2,khi+1);
171 if (gbx.ok()) {
172 Array4<Real> vel_arr = mf_to_fill.array(mfi);
173 ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
174 {
175 vel_arr(i,j,k) = vel_arr(i,j,khi);
176 });
177 }
178 }
179 }
180 }
181}
182
183/**
184 * Fill valid and ghost data in the MultiFab "mf"
185 * This version fills the MultiFab mf_to_fill in valid regions with
186 * the "state data" at the given time;
187 * values in mf when it is passed in are *not* used.
188 * Unlike FillPatch, FillPatchNoBC does not apply boundary conditions
189 *
190 * @param[in ] lev level to FillPatch on
191 * @param[in ] time current time
192 * @param[inout] mf_to_fill MultiFab to fill
193 * @param[in ] mfs vector over levels of multifabs associated with mf_to_fill
194 * @param[in ] bdy_var_type when filling from NetCDF, which boundary data
195 * @param[in ] icomp component to fill
196 * @param[in ] fill_all whether to fill all components
197 * @param[in ] fill_set whether to use FillPatchers
198 */
199void
200REMORA::FillPatchNoBC (int lev, Real time, MultiFab& mf_to_fill, Vector<MultiFab*> const& mfs,
201#ifdef REMORA_USE_NETCDF
202 const int /*bdy_var_type*/,
203#else
204 const int /*bdy_var_type*/,
205#endif
206 const int icomp,
207 const bool fill_all,
208 const bool fill_set)
209{
210 // HACK: Note that this is hacky; should be able to have a single call to FillPatch with a
211 // flag for bcs, but for some reason it was acting weird, so we're splitting this out into
212 // two functions with repeated code for the time being
213 BL_PROFILE_VAR("REMORA::FillPatchNoBC()",REMORA_FillPatch);
214 int bccomp;
215 amrex::Interpolater* mapper = nullptr;
216
217 Box mf_box(mf_to_fill.boxArray()[0]);
218 bool is_2d = mf_box.length(2) == 1;
219
220 //
221 // ***************************************************************************
222 // The first thing we do is interpolate the momenta on the "valid" faces of
223 // the fine grids (where the interface is coarse/fine not fine/fine) -- this
224 // will not be over-written below because the FillPatch operators see these as
225 // valid faces.
226 // ***************************************************************************
227 if (lev>0 && fill_set) {
228 if (cf_set_width > 0 &&
229 mf_box.ixType() == IndexType(IntVect(0,0,0))) {
230 FPr_c[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
231 } else if (cf_set_width >= 0) {
232 if (!is_2d) {
233 if (mf_box.ixType() == IndexType(IntVect(1,0,0))) {
234 FPr_u[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
235 } else if (mf_box.ixType() == IndexType(IntVect(0,1,0))) {
236 FPr_v[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
237 } else if (mf_box.ixType() == IndexType(IntVect(0,0,1))) {
238 FPr_w[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
239 }
240 } else {
241 if (mf_box.ixType() == IndexType(IntVect(1,0,0))) {
242 FPr_ubar[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
243 } else if (mf_box.ixType() == IndexType(IntVect(0,1,0))) {
244 FPr_vbar[lev-1].FillSet(mf_to_fill, time, null_bc, domain_bcs_type);
245 }
246 }
247 }
248 }
249
250 int ncomp;
251 if (fill_all) {
252 ncomp = mf_to_fill.nComp();
253 } else {
254 ncomp = 1;
255 }
256
257 if (mf_box.ixType() == IndexType(IntVect(0,0,0)))
258 {
259 bccomp = 0;
260 mapper = &cell_cons_interp;
261 }
262 else if (mf_box.ixType() == IndexType(IntVect(1,0,0)))
263 {
264 bccomp = BCVars::xvel_bc;
265 mapper = &face_linear_interp;
266 }
267 else if (mf_box.ixType() == IndexType(IntVect(0,1,0)))
268 {
269 bccomp = BCVars::yvel_bc;
270 mapper = &face_linear_interp;
271 }
272 else {
273 bccomp = BCVars::zvel_bc;
274 mapper = &face_linear_interp;
275 }
276
277 if (lev == 0)
278 {
279 Vector<MultiFab*> fmf = {mfs[lev], mfs[lev]};
280 Vector<Real> ftime = {t_old[lev], t_new[lev]};
281 amrex::FillPatchSingleLevel(mf_to_fill, time, fmf, ftime, icomp, icomp, ncomp,
282 geom[lev], null_bc, bccomp);
283 }
284 else
285 {
286 Vector<MultiFab*> fmf = {mfs[lev], mfs[lev]};
287 Vector<Real> ftime = {t_old[lev], t_new[lev]};
288 Vector<MultiFab*> cmf = {mfs[lev-1], mfs[lev-1]};
289 Vector<Real> ctime = {t_old[lev-1], t_new[lev-1]};
290
291 mfs[lev-1]->FillBoundary(geom[lev-1].periodicity());
292 amrex::FillPatchTwoLevels(mf_to_fill, mf_to_fill.nGrowVect(), IntVect(0,0,0),
293 time, cmf, ctime, fmf, ftime,
294 icomp, icomp, ncomp, geom[lev-1], geom[lev],
295 //null_bc, bccomp, null_bc, bccomp,
296 refRatio(lev-1),
297 mapper, domain_bcs_type, bccomp);
298 } // lev > 0
299}
300
301
302// utility to copy in data from old/new data into a struct that holds data for FillPatching
304REMORA::GetDataAtTime (int /*lev*/, Real /*time*/)
305{
306 BL_PROFILE_VAR("GetDataAtTime()",GetDataAtTime);
308
309// HACK HACK HACK
310#if 0
311 const Real teps = (t_new[lev] - t_old[lev]) * 1.e-3_rt;
312
313 if (time > t_new[lev] - teps && time < t_new[lev] + teps)
314 {
315 for (int i = 0; i < Vars::NumTypes; ++i) {
316 data.add_var(&vars_new[lev][i], data.non_owning);
317 }
318 data.set_time(t_new[lev]);
319 }
320 else if (time > t_old[lev] - teps && time < t_old[lev] + teps)
321 {
322 for (int i = 0; i < Vars::NumTypes; ++i) {
323 data.add_var(&vars_old[lev][i], data.non_owning);
324 }
325 data.set_time(t_old[lev]);
326 }
327 else if (time > t_old[lev] && time < t_new[lev])
328 {
329 // do first order interpolation in time between [t_old[lev], t_new[lev]]
330 // time interpolation includes the ghost cells
331 for (int i = 0; i < Vars::NumTypes; ++i) {
332 MultiFab* mf_tmp = new MultiFab(vars_new[lev][i].boxArray(),
333 vars_new[lev][i].DistributionMap(),
334 vars_new[lev][i].nComp(), vars_new[lev][i].nGrowVect());
335 mf_tmp->setVal(0.0_rt);
336
337 const Real dt_fraction = (time - t_old[lev]) / (t_new[lev] - t_old[lev]);
338 MultiFab::Saxpy(*mf_tmp, 1.0_rt - dt_fraction, vars_old[lev][i], 0, 0, mf_tmp->nComp(), mf_tmp->nGrowVect());
339 MultiFab::Saxpy(*mf_tmp, dt_fraction, vars_new[lev][i], 0, 0, mf_tmp->nComp(), mf_tmp->nGrowVect());
340
341 data.add_var(mf_tmp, data.owning);
342 }
343 data.set_time(time);
344 }
345 else
346 {
347 amrex::Error("Requested data at a time outside the interval [t_old, t_new]");
348 }
349
350 // We need to make sure to fill these before we compute the viscosity
351 for (int i = 0; i < Vars::NumTypes; ++i) {
352 data.get_var(Vars::xvel).FillBoundary(geom[lev].periodicity());
353 data.get_var(Vars::yvel).FillBoundary(geom[lev].periodicity());
354 data.get_var(Vars::zvel).FillBoundary(geom[lev].periodicity());
355 data.get_var(Vars::cons).FillBoundary(geom[lev].periodicity());
356 }
357#endif
358 return data;
359}
360
361
362/**
363 * Like FillCoarsePatch but uses piecewise constant interpolater
364 *
365 * @param[in ] lev level to FillPatch on
366 * @param[in ] time current time
367 * @param[inout] mf_to_fill MultiFab to fill
368 * @param[inout] mf_crse MultiFab of coarse data
369 * @param[in ] bccomp index into domain_bcs_type
370 * @param[in ] bdy_var_type when filling from NetCDF, which boundary data
371 * @param[in ] icomp component to fill
372 * @param[in ] fill_all whether to fill all components
373 * @param[in ] icomp_calc component to use in RHS boundary calculation
374 * @param[in ] dt_lev time step at this level
375 * @param[in ] mf_calc data to use in RHS boundary calculation
376 */
377void
378REMORA::FillCoarsePatchPC (int lev, Real time, MultiFab* mf_to_fill,
379 MultiFab* mf_crse,
380 const int bccomp,
381 const int bdy_var_type,
382 const int icomp,
383 const bool fill_all,
384 const int n_not_fill,
385 const int icomp_calc,
386 const Real dt_lev,
387 const MultiFab& mf_calc) {
388 amrex::Interpolater* mapper = &pc_interp;
389 FillCoarsePatchMap(lev, time, mf_to_fill, mf_crse, bccomp, bdy_var_type,icomp,
390 fill_all, n_not_fill, icomp_calc, dt_lev, mf_calc, mapper);
391}
392
393/**
394 * Fill an entire multifab by interpolating from the coarser level --
395 * this is used only when a new level of refinement is being created
396 * during a run (i.e not at initialization)
397 * This will never be used with static refinement.
398 *
399 * @param[in ] lev level to FillPatch on
400 * @param[in ] time current time
401 * @param[inout] mf_to_fill MultiFab to fill
402 * @param[inout] mf_crse MultiFab of coarse data
403 * @param[in ] bccomp index into domain_bcs_type
404 * @param[in ] bdy_var_type when filling from NetCDF, which boundary data
405 * @param[in ] icomp component to fill
406 * @param[in ] fill_all whether to fill all components
407 * @param[in ] icomp_calc component to use in RHS boundary calculation
408 * @param[in ] dt_lev time step at this level
409 * @param[in ] mf_calc data to use in RHS boundary calculation
410 */
411void
412REMORA::FillCoarsePatch (int lev, Real time, MultiFab* mf_to_fill,
413 MultiFab* mf_crse,
414 const int bccomp,
415 const int bdy_var_type,
416 const int icomp,
417 const bool fill_all,
418 const int n_not_fill,
419 const int icomp_calc,
420 const Real dt_lev,
421 const MultiFab& mf_calc) {
422 amrex::Interpolater* mapper = nullptr;
423 FillCoarsePatchMap(lev, time, mf_to_fill, mf_crse, bccomp, bdy_var_type,icomp,
424 fill_all, n_not_fill, icomp_calc, dt_lev, mf_calc, mapper);
425}
426
427/**
428 * Fill an entire multifab by interpolating from the coarser level --
429 * this is used only when a new level of refinement is being created
430 * during a run (i.e not at initialization)
431 * This will never be used with static refinement.
432 *
433 * @param[in ] lev level to FillPatch on
434 * @param[in ] time current time
435 * @param[inout] mf_to_fill MultiFab to fill
436 * @param[inout] mf_crse MultiFab of coarse data
437 * @param[in ] bccomp index into domain_bcs_type
438 * @param[in ] bdy_var_type when filling from NetCDF, which boundary data
439 * @param[in ] icomp component to fill
440 * @param[in ] fill_all whether to fill all components
441 * @param[in ] icomp_calc component to use in RHS boundary calculation
442 * @param[in ] dt_lev time step at this level
443 * @param[in ] mf_calc data to use in RHS boundary calculation
444 */
445void
446REMORA::FillCoarsePatchMap (int lev, Real time, MultiFab* mf_to_fill, MultiFab* mf_crse,
447 const int bccomp,
448#ifdef REMORA_USE_NETCDF
449 const int bdy_var_type,
450#else
451 const int /*bdy_var_type*/,
452#endif
453 const int icomp,
454 const bool fill_all,
455 const int n_not_fill,
456#ifdef REMORA_USE_NETCDF
457 const int icomp_calc,
458 const Real dt_lev,
459#else
460 const int /*icomp_calc*/,
461 const Real /*dt_lev*/,
462#endif
463 const MultiFab& mf_calc,
464 amrex::Interpolater* mapper)
465{
466 BL_PROFILE_VAR("FillCoarsePatch()",FillCoarsePatch);
467 AMREX_ASSERT(lev > 0);
468
469 MultiFab* mask;
470
471 int ncomp;
472 if (fill_all) {
473 ncomp = mf_to_fill->nComp();
474 } else {
475 ncomp = 1;
476 }
477
478 Box box_mf(mf_to_fill->boxArray()[0]);
479
480 if (box_mf.ixType() == IndexType(IntVect(0,0,0)))
481 {
482 mask = vec_mskr[lev].get();
483 }
484 else if (box_mf.ixType() == IndexType(IntVect(1,0,0)))
485 {
486 mask = vec_msku[lev].get();
487 }
488 else if (box_mf.ixType() == IndexType(IntVect(0,1,0)))
489 {
490 mask = vec_mskv[lev].get();
491 }
492 else if (box_mf.ixType() == IndexType(IntVect(0,0,1)))
493 {
494 mask = vec_mskr[lev].get();
495 } else {
496 amrex::Abort("Dont recognize this box type in REMORA_FillPatch");
497 }
498
499 if (mapper == nullptr) {
500 if (box_mf.ixType() == IndexType(IntVect(0,0,0)))
501 {
502 mapper = &cell_cons_interp;
503 }
504 else if (box_mf.ixType() == IndexType(IntVect(1,0,0)))
505 {
506 mapper = &face_linear_interp;
507 }
508 else if (box_mf.ixType() == IndexType(IntVect(0,1,0)))
509 {
510 mapper = &face_linear_interp;
511 }
512 else if (box_mf.ixType() == IndexType(IntVect(0,0,1)))
513 {
514 mapper = &face_linear_interp;
515 } else {
516 amrex::Abort("Dont recognize this box type in REMORA_FillPatch");
517 }
518 }
519
520#if 0
521 TimeInterpolatedData cdata = GetDataAtTime(lev-1, time);
522 TimeInterpolatedData fdata = GetDataAtTime(lev , time);
523 Vector<Real> ctime = {cdata.get_time()};
524
525 Vector<MultiFab*> cmf = {mf_crse};
526 Vector<Real> ctime = {time};
527
528 REMORAPhysBCFunct cphysbc(lev-1,geom[lev-1],
530 cdata,
532#ifdef REMORA_USE_NETCDF
533 ,ic_type,bdy_data_xlo,bdy_data_xhi,
535#endif
536 );
537 REMORAPhysBCFunct fphysbc(lev,geom[lev],
539 fdata,
541#ifdef REMORA_USE_NETCDF
542 ,ic_type,bdy_data_xlo,bdy_data_xhi,
544#endif
545 );
546#endif
547
548
549 mf_crse->FillBoundary(geom[lev-1].periodicity());
550 amrex::InterpFromCoarseLevel(*mf_to_fill, mf_to_fill->nGrowVect(), IntVect(0,0,0),
551 *mf_crse, 0, icomp, ncomp, geom[lev-1], geom[lev],
552 refRatio(lev-1),
553 mapper, domain_bcs_type, bccomp);
554
555 // ***************************************************************************
556 // Physical bc's at domain boundary
557 // ***************************************************************************
558
559 // Enforce physical boundary conditions
560 (*physbcs[lev])(*mf_to_fill,*mask,icomp,ncomp,mf_to_fill->nGrowVect(),
561 time,bccomp,n_not_fill,mf_calc,
562 *vec_msku[lev].get(), *vec_mskv[lev].get());
563
564#ifdef REMORA_USE_NETCDF
565 // Fill the data which is stored in the boundary data read from netcdf files
566 if ( (solverChoice.boundary_from_netcdf) && (lev==0) &&
567 (bdy_var_type != BdyVars::null) )
568 {
569 fill_from_bdyfiles(*mf_to_fill,*mask,time,bccomp,bdy_var_type, icomp,icomp_calc,mf_calc,dt_lev);
570 }
571#endif
572 // Fill corners of the domain with periodic data
573 if ( box_mf.ixType() == IndexType(IntVect(0,0,0)) ) {
574 mf_to_fill->EnforcePeriodicity(geom[lev].periodicity());
575 }
576
577 // Enforce free-slip at top boundary (on xvel or yvel)
578 if ( (box_mf.ixType() == IndexType(IntVect(1,0,0))) ||
579 (box_mf.ixType() == IndexType(IntVect(0,1,0))) )
580 {
581 int khi = geom[lev].Domain().bigEnd(2);
582 for (MFIter mfi(*mf_to_fill); mfi.isValid(); ++mfi)
583 {
584 Box gbx = mfi.growntilebox(); // Note this is face-centered since vel is
585 gbx.setSmall(2,khi+1);
586 if (gbx.ok()) {
587 Array4<Real> vel_arr = mf_to_fill->array(mfi);
588 ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
589 {
590 vel_arr(i,j,k) = vel_arr(i,j,khi);
591 });
592 }
593 }
594 }
595
596}
597
598/**
599 * @param[in ] lev level to fill at
600 * @param[in ] lev MultiFab of cell-centered velocities
601 */
602void
603REMORA::FillBdyCCVels (int lev, MultiFab& mf_cc_vel)
604{
605 // Impose bc's at domain boundaries
606// for (int lev = 0; lev <= finest_level; ++lev)
607// {
608 Box domain(Geom(lev).Domain());
609
610 int ihi = domain.bigEnd(0);
611 int jhi = domain.bigEnd(1);
612 int khi = domain.bigEnd(2);
613
614 // Impose periodicity first
615 mf_cc_vel.FillBoundary(geom[lev].periodicity());
616
617#ifdef _OPENMP
618#pragma omp parallel if (Gpu::notInLaunchRegion())
619#endif
620 for (MFIter mfi(mf_cc_vel, TilingIfNotGPU()); mfi.isValid(); ++mfi)
621 {
622 // Note that we don't fill corners here -- only the cells that share a face
623 // with interior cells -- this is all that is needed to calculate vorticity
624 const Box& bx = mfi.tilebox();
625 const Array4<Real>& vel_arr = mf_cc_vel.array(mfi);
626
627 if (!Geom(lev).isPeriodic(0)) {
628 // Low-x side
629 if (bx.smallEnd(0) <= domain.smallEnd(0)) {
630 Real mult = (phys_bc_type[BCVars::xvel_bc][0] == REMORA_BC::no_slip_wall) ? -1. : 1.;
631 ParallelFor(makeSlab(bx,0,0), [=] AMREX_GPU_DEVICE(int , int j, int k) noexcept
632 {
633 vel_arr(-1,j,k,1) = mult*vel_arr(0,j,k,1); // v
634 vel_arr(-1,j,k,2) = mult*vel_arr(0,j,k,2); // w
635 });
636 }
637
638 // High-x side
639 if (bx.bigEnd(0) >= domain.bigEnd(0)) {
640 Real mult = (phys_bc_type[BCVars::xvel_bc][3] == REMORA_BC::no_slip_wall) ? -1. : 1.;
641 ParallelFor(makeSlab(bx,0,0), [=] AMREX_GPU_DEVICE(int , int j, int k) noexcept
642 {
643 vel_arr(ihi+1,j,k,1) = mult*vel_arr(ihi,j,k,1); // v
644 vel_arr(ihi+1,j,k,2) = mult*vel_arr(ihi,j,k,2); // w
645 });
646 }
647 } // !periodic
648
649 if (!Geom(lev).isPeriodic(1)) {
650 // Low-y side
651 if (bx.smallEnd(1) <= domain.smallEnd(1)) {
652 Real mult = (phys_bc_type[BCVars::yvel_bc][1] == REMORA_BC::no_slip_wall) ? -1. : 1.;
653 ParallelFor(makeSlab(bx,1,0), [=] AMREX_GPU_DEVICE(int i, int , int k) noexcept
654 {
655 vel_arr(i,-1,k,0) = mult*vel_arr(i,0,k,0); // u
656 vel_arr(i,-1,k,2) = mult*vel_arr(i,0,k,2); // w
657 });
658 }
659
660 // High-y side
661 if (bx.bigEnd(1) >= domain.bigEnd(1)) {
662 Real mult = (phys_bc_type[BCVars::yvel_bc][4] == REMORA_BC::no_slip_wall) ? -1. : 1.;
663 ParallelFor(makeSlab(bx,1,0), [=] AMREX_GPU_DEVICE(int i, int , int k) noexcept
664 {
665 vel_arr(i,jhi+1,k,0) = mult*vel_arr(i,jhi,k,0); // u
666 vel_arr(i,jhi+1,k,2) = mult*-vel_arr(i,jhi,k,2); // w
667 });
668 }
669 } // !periodic
670
671 if (!Geom(lev).isPeriodic(2)) {
672 // Low-z side
673 if (bx.smallEnd(2) <= domain.smallEnd(2)) {
674 Real mult = (phys_bc_type[BCVars::zvel_bc][2] == REMORA_BC::no_slip_wall) ? -1. : 1.;
675 ParallelFor(makeSlab(bx,2,0), [=] AMREX_GPU_DEVICE(int i, int j, int) noexcept
676 {
677 vel_arr(i,j,-1,0) = mult*vel_arr(i,j,0,0); // u
678 vel_arr(i,j,-1,1) = mult*vel_arr(i,j,0,1); // v
679 });
680 }
681
682 // High-z side
683 if (bx.bigEnd(2) >= domain.bigEnd(2)) {
684 Real mult = (phys_bc_type[BCVars::zvel_bc][5] == REMORA_BC::no_slip_wall) ? -1. : 1.;
685 ParallelFor(makeSlab(bx,2,0), [=] AMREX_GPU_DEVICE(int i, int j, int) noexcept
686 {
687 vel_arr(i,j,khi+1,0) = mult*vel_arr(i,j,khi,0); // u
688 vel_arr(i,j,khi+1,1) = mult*vel_arr(i,j,khi,1); // v
689 });
690 }
691 } // !periodic
692 } // MFIter
693
694// } // lev
695}
PhysBCFunctNoOp null_bc
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=amrex::Real(0.0), const amrex::MultiFab &mf_calc=amrex::MultiFab())
fill an entire multifab by interpolating from the coarser level using the piecewise constant interpol...
amrex::Vector< amrex::BCRec > domain_bcs_type
vector (over BCVars) of BCRecs
Definition REMORA.H:1178
amrex::Vector< REMORAFillPatcher > FPr_v
Vector over levels of FillPatchers for v (3D)
Definition REMORA.H:1113
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=amrex::Real(0.0), const amrex::MultiFab &mf_calc=amrex::MultiFab())
fill an entire multifab by interpolating from the coarser level
amrex::GpuArray< amrex::GpuArray< REMORA_BC, AMREX_SPACEDIM *2 >, BCVars::NumTypes > phys_bc_type
Array holding the "physical" boundary condition types (e.g. "inflow")
Definition REMORA.H:1189
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskr
land/sea mask at cell centers (2D)
Definition REMORA.H:357
amrex::Vector< REMORAFillPatcher > FPr_u
Vector over levels of FillPatchers for u (3D)
Definition REMORA.H:1111
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_msku
land/sea mask at x-faces (2D)
Definition REMORA.H:359
void FillPatchNoBC(int lev, amrex::Real time, amrex::MultiFab &mf_to_be_filled, amrex::Vector< amrex::MultiFab * > const &mfs, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const bool fill_set=true)
Fill a new MultiFab by copying in phi from valid region and filling ghost cells without applying boun...
amrex::Vector< amrex::Vector< amrex::FArrayBox > > bdy_data_yhi
Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Northern boundary data from ...
Definition REMORA.H:1054
amrex::Real bdy_time_interval
Interval between boundary data times.
Definition REMORA.H:1059
amrex::Vector< REMORAFillPatcher > FPr_vbar
Vector over levels of FillPatchers for vbar (2D)
Definition REMORA.H:1121
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskv
land/sea mask at y-faces (2D)
Definition REMORA.H:361
amrex::Vector< std::unique_ptr< REMORAPhysBCFunct > > physbcs
Vector (over level) of functors to apply physical boundary conditions.
Definition REMORA.H:1166
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=true, const int n_not_fill=0, const int icomp_calc=0, const amrex::Real dt=amrex::Real(0.0), const amrex::MultiFab &mf_calc=amrex::MultiFab())
Fill a new MultiFab by copying in phi from valid region and filling ghost cells.
amrex::Vector< amrex::Vector< amrex::FArrayBox > > bdy_data_ylo
Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Southern boundary data from ...
Definition REMORA.H:1052
void fill_from_bdyfiles(amrex::MultiFab &mf_to_fill, const amrex::MultiFab &mf_mask, const amrex::Real time, const int bccomp, const int bdy_var_type, const int icomp_to_fill, const int icomp_calc=0, const amrex::MultiFab &mf_calc=amrex::MultiFab(), const amrex::Real=amrex::Real(0.0))
Fill boundary data from netcdf file.
void FillBdyCCVels(int lev, amrex::MultiFab &mf_cc_vel)
Fill the physical boundary conditions for cell-centered velocity (diagnostic only)
TimeInterpolatedData GetDataAtTime(int lev, amrex::Real time)
utility to copy in data from old and/or new state into another multifab
amrex::Vector< REMORAFillPatcher > FPr_c
Vector over levels of FillPatchers for scalars.
Definition REMORA.H:1109
amrex::Vector< REMORAFillPatcher > FPr_w
Vector over levels of FillPatchers for w.
Definition REMORA.H:1115
amrex::Vector< amrex::Real > t_new
new time at each level
Definition REMORA.H:1156
static SolverChoice solverChoice
Container for algorithmic choices.
Definition REMORA.H:1279
int cf_set_width
Width for fixing values at coarse-fine interface.
Definition REMORA.H:1106
amrex::Vector< amrex::Vector< amrex::FArrayBox > > bdy_data_xhi
Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Eastern boundary data from f...
Definition REMORA.H:1050
void FillCoarsePatchMap(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=amrex::Real(0.0), const amrex::MultiFab &mf_calc=amrex::MultiFab(), amrex::Interpolater *mapper=nullptr)
fill an entire multifab by interpolating from the coarser level, explicitly specifying interpolator t...
amrex::Array< amrex::Array< amrex::Real, AMREX_SPACEDIM *2 >, AMREX_SPACEDIM+NCONS+8 > m_bc_extdir_vals
Array holding the Dirichlet values at walls which need them.
Definition REMORA.H:1186
amrex::Vector< REMORAFillPatcher > FPr_ubar
Vector over levels of FillPatchers for ubar (2D)
Definition REMORA.H:1119
amrex::Vector< amrex::Vector< amrex::FArrayBox > > bdy_data_xlo
Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Western boundary data from f...
Definition REMORA.H:1048
amrex::Vector< amrex::Real > t_old
old time at each level
Definition REMORA.H:1158
amrex::Gpu::DeviceVector< amrex::BCRec > domain_bcs_type_d
GPU vector (over BCVars) of BCRecs.
Definition REMORA.H:1180
void set_time(amrex::Real time)
void add_var(amrex::MultiFab *var_data, int own_data)
amrex::MultiFab & get_var(int var_idx)