Block-Structured AMR Software Framework
AMReX_EB2_Level.H
Go to the documentation of this file.
1 #ifndef AMREX_EB2_LEVEL_H_
2 #define AMREX_EB2_LEVEL_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_ParmParse.H>
6 #include <AMReX_Geometry.H>
7 #include <AMReX_MultiFab.H>
8 #include <AMReX_iMultiFab.H>
9 #include <AMReX_LayoutData.H>
10 #include <AMReX_VisMF.H>
11 #include <AMReX_Array.H>
12 #include <AMReX_EBCellFlag.H>
13 #include <AMReX_MultiCutFab.H>
14 #include <AMReX_EB2_MultiGFab.H>
15 #include <AMReX_EB2_C.H>
17 
18 #ifdef AMREX_USE_OMP
19 #include <omp.h>
20 #endif
21 
22 #include <unordered_map>
23 #include <limits>
24 #include <cmath>
25 #include <type_traits>
26 
27 namespace amrex::EB2 {
28 
29 class IndexSpace;
30 template <typename G> class GShopLevel;
31 
32 class Level
33 {
34 public:
35 
36  bool isAllRegular () const noexcept { return m_allregular; }
37  bool isOK () const noexcept { return m_ok; }
38  void fillEBCellFlag (FabArray<EBCellFlagFab>& cellflag, const Geometry& geom) const;
39  void fillVolFrac (MultiFab& vfrac, const Geometry& geom) const;
40  void fillCentroid (MultiCutFab& centroid, const Geometry& geom) const;
41  void fillCentroid ( MultiFab& centroid, const Geometry& geom) const;
42  void fillBndryArea (MultiCutFab& bndryarea, const Geometry& geom) const;
43  void fillBndryArea ( MultiFab& bndryarea, const Geometry& geom) const;
44  void fillBndryCent (MultiCutFab& bndrycent, const Geometry& geom) const;
45  void fillBndryCent ( MultiFab& bndrycent, const Geometry& geom) const;
46  void fillBndryNorm (MultiCutFab& bndrynorm, const Geometry& geom) const;
47  void fillBndryNorm ( MultiFab& bndrynorm, const Geometry& geom) const;
48  void fillAreaFrac (Array<MultiCutFab*,AMREX_SPACEDIM> const& areafrac, const Geometry& geom) const;
49  void fillAreaFrac (Array< MultiFab*,AMREX_SPACEDIM> const& areafrac, const Geometry& geom) const;
50  void fillFaceCent (Array<MultiCutFab*,AMREX_SPACEDIM> const& facecent, const Geometry& geom) const;
51  void fillFaceCent (Array< MultiFab*,AMREX_SPACEDIM> const& facecent, const Geometry& geom) const;
52  void fillEdgeCent (Array<MultiCutFab*,AMREX_SPACEDIM> const& edgecent, const Geometry& geom) const;
53  void fillEdgeCent (Array< MultiFab*,AMREX_SPACEDIM> const& edgecent, const Geometry& geom) const;
54  void fillLevelSet (MultiFab& levelset, const Geometry& geom) const;
55 
56  const BoxArray& boxArray () const noexcept { return m_grids; }
57  const DistributionMapping& DistributionMap () const noexcept { return m_dmap; }
58 
59  Level (IndexSpace const* is, const Geometry& geom) : m_geom(geom), m_parent(is) {}
60  void prepareForCoarsening (const Level& rhs, int max_grid_size, IntVect const& ngrow);
61 
62  const Geometry& Geom () const noexcept { return m_geom; }
63  IndexSpace const* getEBIndexSpace () const noexcept { return m_parent; }
64 
65  IntVect const& nGrowVect () const noexcept { return m_ngrow; }
66 
67  void write_to_chkpt_file (const std::string& fname, bool extend_domain_face, int max_grid_size) const;
68 
69  bool hasEBInfo () const noexcept { return m_has_eb_info; }
70  void fillCutCellMask (iMultiFab& cutcellmask, const Geometry& geom) const;
71 
72 // public: // for cuda
73  int coarsenFromFine (Level& fineLevel, bool fill_boundary);
74  void buildCellFlag ();
75  void buildCutCellMask (Level const& fine_level);
76 
77 protected:
78 
96  bool m_allregular = false;
97  bool m_ok = false;
98  bool m_has_eb_info = true;
100 
101 private:
102  template <typename G> friend class GShopLevel;
103 
104  // Need this function to work around a gcc bug.
105  void setRegularLevel () {
106  m_allregular = true;
107  m_ok = true;
108  m_has_eb_info = false;
109  }
110 };
111 
112 template <typename G>
114  : public Level
115 {
116 public:
117  GShopLevel (IndexSpace const* is, G const& gshop, const Geometry& geom, int max_grid_size,
118  int ngrow, bool extend_domain_face, int num_crse_opt);
119  GShopLevel (IndexSpace const* is, int ilev, int max_grid_size, int ngrow,
120  const Geometry& geom, GShopLevel<G>& fineLevel);
121  GShopLevel (IndexSpace const* is, const Geometry& geom);
122  void define_fine (G const& gshop, const Geometry& geom,
123  int max_grid_size, int ngrow, bool extend_domain_face, int num_crse_opt);
124 
125  static GShopLevel<G>
126  makeAllRegular(IndexSpace const* is, const Geometry& geom)
127  {
128  GShopLevel<G> r(is, geom);
129  r.setRegularLevel();
130  return r;
131  }
132 };
133 
134 template <typename G>
136  : Level(is, geom)
137 {}
138 
139 template <typename G>
140 GShopLevel<G>::GShopLevel (IndexSpace const* is, G const& gshop, const Geometry& geom,
141  int max_grid_size, int ngrow, bool extend_domain_face, int num_crse_opt)
142  : Level(is, geom)
143 {
144  if (std::is_same<typename G::FunctionType, AllRegularIF>::value) {
145  m_allregular = true;
146  m_ok = true;
147  return;
148  }
149 
150  define_fine(gshop, geom, max_grid_size, ngrow, extend_domain_face, num_crse_opt);
151 }
152 
153 template <typename G>
154 void
155 GShopLevel<G>::define_fine (G const& gshop, const Geometry& geom,
156  int max_grid_size, int ngrow, bool extend_domain_face, int num_crse_opt)
157 {
158  if (amrex::Verbose() > 0 && extend_domain_face == false) {
159  amrex::Print() << "AMReX WARNING: extend_domain_face=false is not recommended!\n";
160  }
161 
162  BL_PROFILE("EB2::GShopLevel()-fine");
163 
164 #ifdef AMREX_USE_FLOAT
165  Real small_volfrac = 1.e-5_rt;
166 #else
167  Real small_volfrac = 1.e-14;
168 #endif
169  bool cover_multiple_cuts = false;
170  int maxiter = 32;
171  {
172  ParmParse pp("eb2");
173  pp.queryAdd("small_volfrac", small_volfrac);
174  pp.queryAdd("cover_multiple_cuts", cover_multiple_cuts);
175  pp.queryAdd("maxiter", maxiter);
176  }
177  maxiter = std::min(100000, maxiter);
178 
179  // make sure ngrow is multiple of 16
180  m_ngrow = IntVect{static_cast<int>(std::ceil(ngrow/16.)) * 16};
181 
182  Box const& domain = geom.Domain();
183  Box domain_grown = domain;
184  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
185  if (geom.isPeriodic(idim)) {
186  m_ngrow[idim] = 0;
187  } else {
188  m_ngrow[idim] = std::min(m_ngrow[idim], domain_grown.length(idim));
189  }
190  }
191  domain_grown.grow(m_ngrow);
192  Box bounding_box = (extend_domain_face) ? domain : domain_grown;
193  bounding_box.surroundingNodes();
194 
195  BoxList cut_boxes;
196  BoxList covered_boxes;
197 
198  const int nprocs = ParallelDescriptor::NProcs();
199  const int iproc = ParallelDescriptor::MyProc();
200 
201  num_crse_opt = std::max(0,std::min(8,num_crse_opt));
202  for (int clev = num_crse_opt; clev >= 0; --clev) {
203  IntVect crse_ratio(1 << clev);
204  if (domain.coarsenable(crse_ratio)) {
205  Box const& crse_bounding_box = amrex::coarsen(bounding_box, crse_ratio);
206  Geometry const& crse_geom = amrex::coarsen(geom, crse_ratio);
207  BoxList test_boxes;
208  if (cut_boxes.isEmpty()) {
209  covered_boxes.clear();
210  test_boxes = BoxList(crse_geom.Domain());
211  test_boxes.maxSize(max_grid_size);
212  } else {
213  test_boxes.swap(cut_boxes);
214  test_boxes.coarsen(crse_ratio);
215  test_boxes.maxSize(max_grid_size);
216  }
217 
218  const Long nboxes = test_boxes.size();
219  const auto& boxes = test_boxes.data();
220  for (Long i = iproc; i < nboxes; i += nprocs) {
221  const Box& vbx = boxes[i];
222  const Box& gbx = amrex::surroundingNodes(amrex::grow(vbx,1));
223  auto box_type = gshop.getBoxType(gbx&crse_bounding_box,crse_geom,RunOn::Gpu);
224  if (box_type == gshop.allcovered) {
225  covered_boxes.push_back(amrex::refine(vbx, crse_ratio));
226  } else if (box_type == gshop.mixedcells) {
227  cut_boxes.push_back(amrex::refine(vbx, crse_ratio));
228  }
229  }
230 
231  amrex::AllGatherBoxes(cut_boxes.data());
232  }
233  }
234 
235  amrex::AllGatherBoxes(covered_boxes.data());
236 
237  if (m_ngrow != 0) {
238  auto grow_at_domain_boundary = [&] (BoxList& bl)
239  {
240  const IntVect& domlo = domain.smallEnd();
241  const IntVect& domhi = domain.bigEnd();
242  for (auto& b : bl) {
243  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
244  if (m_ngrow[idim] != 0) {
245  if (b.smallEnd(idim) == domlo[idim]) {
246  b.growLo(idim,m_ngrow[idim]);
247  }
248  if (b.bigEnd(idim) == domhi[idim]) {
249  b.growHi(idim,m_ngrow[idim]);
250  }
251  }
252  }
253  }
254  };
255  grow_at_domain_boundary(covered_boxes);
256  grow_at_domain_boundary(cut_boxes);
257  }
258 
259  if ( cut_boxes.isEmpty() &&
260  !covered_boxes.isEmpty())
261  {
262  amrex::Abort("AMReX_EB2_Level.H: Domain is completely covered");
263  }
264 
265  if (!covered_boxes.isEmpty()) {
266  if (num_crse_opt > 2) { // don't want the box too big
267  covered_boxes.maxSize(max_grid_size*4);
268  }
269  m_covered_grids = BoxArray(std::move(covered_boxes));
270  }
271 
272  if (cut_boxes.isEmpty()) {
273  m_grids = BoxArray();
274  m_dmap = DistributionMapping();
275  m_allregular = true;
276  m_ok = true;
277  return;
278  }
279 
280  m_grids = BoxArray(std::move(cut_boxes));
281  m_dmap = DistributionMapping(m_grids);
282 
283  m_mgf.define(m_grids, m_dmap);
284  const int ng = GFab::ng;
285  MFInfo mf_info;
286  mf_info.SetTag("EB2::Level");
287  m_cellflag.define(m_grids, m_dmap, 1, ng, mf_info);
288  m_volfrac.define(m_grids, m_dmap, 1, ng, mf_info);
289  m_centroid.define(m_grids, m_dmap, AMREX_SPACEDIM, ng, mf_info);
290  m_bndryarea.define(m_grids, m_dmap, 1, ng, mf_info);
291  m_bndrycent.define(m_grids, m_dmap, AMREX_SPACEDIM, ng, mf_info);
292  m_bndrynorm.define(m_grids, m_dmap, AMREX_SPACEDIM, ng, mf_info);
293  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
294  m_areafrac[idim].define(amrex::convert(m_grids, IntVect::TheDimensionVector(idim)),
295  m_dmap, 1, ng, mf_info);
296  m_facecent[idim].define(amrex::convert(m_grids, IntVect::TheDimensionVector(idim)),
297  m_dmap, AMREX_SPACEDIM-1, ng, mf_info);
298  IntVect edge_type{1}; edge_type[idim] = 0;
299  m_edgecent[idim].define(amrex::convert(m_grids, edge_type), m_dmap, 1, ng, mf_info);
300  }
301 
302  const auto dx = geom.CellSizeArray();
303  const auto problo = geom.ProbLoArray();
304 
305  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
306  if (!extend_domain_face || geom.isPeriodic(idim)) {
307  bounding_box.grow(idim,GFab::ng);
308  }
309  }
310 
311  RunOn gshop_run_on = (Gpu::inLaunchRegion() && gshop.isGPUable())
313 
314  bool hybrid = Gpu::inLaunchRegion() && (gshop_run_on == RunOn::Cpu);
315  amrex::ignore_unused(hybrid);
316 
317  int iter = 0;
318  for (; iter < maxiter; ++iter)
319  {
320  int nsmallcells = 0;
321  int nmulticuts = 0;
322 #ifdef AMREX_USE_OMP
323 #pragma omp parallel if (Gpu::notInLaunchRegion()) reduction(+:nsmallcells,nmulticuts)
324 #endif
325  {
326 #if (AMREX_SPACEDIM == 3)
327  Array<BaseFab<Real>, AMREX_SPACEDIM> M2;
328  EBCellFlagFab cellflagtmp;
329 #endif
330  for (MFIter mfi(m_mgf); mfi.isValid(); ++mfi)
331  {
332  auto& gfab = m_mgf[mfi];
333  const Box& vbx = gfab.validbox();
334 
335  auto& levelset = gfab.getLevelSet();
336  if (iter == 0) {
337  gshop.fillFab(levelset, geom, gshop_run_on, bounding_box);
338 #ifdef AMREX_USE_GPU
339  if (hybrid) {
340  levelset.prefetchToDevice();
341  }
342 #endif
343  }
344 
345  auto& cellflag = m_cellflag[mfi];
346 
347  gfab.buildTypes(cellflag);
348 
349  Array4<Real const> const& clst = levelset.const_array();
350  Array4<Real > const& lst = levelset.array();
351  Array4<EBCellFlag> const& cfg = m_cellflag.array(mfi);
352  Array4<Real> const& vfr = m_volfrac.array(mfi);
353  Array4<Real> const& ctr = m_centroid.array(mfi);
354  Array4<Real> const& bar = m_bndryarea.array(mfi);
355  Array4<Real> const& bct = m_bndrycent.array(mfi);
356  Array4<Real> const& bnm = m_bndrynorm.array(mfi);
357  AMREX_D_TERM(Array4<Real> const& apx = m_areafrac[0].array(mfi);,
358  Array4<Real> const& apy = m_areafrac[1].array(mfi);,
359  Array4<Real> const& apz = m_areafrac[2].array(mfi););
360  AMREX_D_TERM(Array4<Real> const& fcx = m_facecent[0].array(mfi);,
361  Array4<Real> const& fcy = m_facecent[1].array(mfi);,
362  Array4<Real> const& fcz = m_facecent[2].array(mfi););
363 
364  auto& facetype = gfab.getFaceType();
365  AMREX_D_TERM(Array4<Type_t> const& ftx = facetype[0].array();,
366  Array4<Type_t> const& fty = facetype[1].array();,
367  Array4<Type_t> const& ftz = facetype[2].array(););
368 
369  int nmc = 0;
370  int nsm = 0;
371 
372 #if (AMREX_SPACEDIM == 3)
373  auto& edgetype = gfab.getEdgeType();
374  Array4<Type_t const> const& xdg = edgetype[0].const_array();
375  Array4<Type_t const> const& ydg = edgetype[1].const_array();
376  Array4<Type_t const> const& zdg = edgetype[2].const_array();
377 
378  Array4<Real> const& xip = m_edgecent[0].array(mfi);
379  Array4<Real> const& yip = m_edgecent[1].array(mfi);
380  Array4<Real> const& zip = m_edgecent[2].array(mfi);
381 
382  if (iter == 0) {
383 #ifdef AMREX_USE_GPU
384  if (hybrid) {
386  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
387  edgetype[idim].prefetchToHost();
388  m_edgecent[idim][mfi].prefetchToHost();
389  }
390  }
391 #endif
392 
393  gshop.getIntercept({xip,yip,zip}, {xdg,ydg,zdg}, clst,
394  geom, gshop_run_on, bounding_box);
395 
396 #ifdef AMREX_USE_GPU
397  if (hybrid) {
398  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
399  edgetype[idim].prefetchToDevice();
400  m_edgecent[idim][mfi].prefetchToDevice();
401  }
402  }
403 #endif
404  } else {
405  gshop.updateIntercept({xip,yip,zip}, {xdg,ydg,zdg}, clst, geom);
406  }
407 
408  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
409  const Box& b = facetype[idim].box();
410  M2[idim].resize(b,3);
411  }
412  Array4<Real> const& xm2 = M2[0].array();
413  Array4<Real> const& ym2 = M2[1].array();
414  Array4<Real> const& zm2 = M2[2].array();
415 
416  nmc = build_faces(vbx, cfg, ftx, fty, ftz, xdg, ydg, zdg, lst,
417  xip, yip, zip, apx, apy, apz, fcx, fcy, fcz,
418  xm2, ym2, zm2, dx, problo, cover_multiple_cuts);
419 
420  cellflagtmp.resize(m_cellflag[mfi].box());
421  Elixir cellflagtmp_eli = cellflagtmp.elixir();
422  Array4<EBCellFlag> const& cfgtmp = cellflagtmp.array();
423 
424  build_cells(vbx, cfg, ftx, fty, ftz, apx, apy, apz,
425  fcx, fcy, fcz, xm2, ym2, zm2, dx, vfr, ctr,
426  bar, bct, bnm, cfgtmp, lst,
427  small_volfrac, geom, extend_domain_face, cover_multiple_cuts,
428  nsm, nmc);
429 
430  // Because it is used in a synchronous reduction kernel in
431  // build_cells, we do not need to worry about M2's lifetime.
432  // But we still need to use Elixir to extend the life of
433  // cellflagtmp.
434 
435 #elif (AMREX_SPACEDIM == 2)
436  Array4<Real> const& xip = m_edgecent[0].array(mfi);
437  Array4<Real> const& yip = m_edgecent[1].array(mfi);
438 
439  if (iter == 0) {
440 #ifdef AMREX_USE_GPU
441  if (hybrid) {
443  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
444  facetype[idim].prefetchToHost();
445  m_edgecent[idim][mfi].prefetchToHost();
446  }
447  }
448 #endif
449 
450  // yes, factype[1] and then [0]
451  gshop.getIntercept({xip,yip},
452  {facetype[1].const_array(), facetype[0].const_array()},
453  clst, geom, gshop_run_on, bounding_box);
454 
455 #ifdef AMREX_USE_GPU
456  if (hybrid) {
457  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
458  facetype[idim].prefetchToDevice();
459  m_edgecent[idim][mfi].prefetchToDevice();
460  }
461  }
462 #endif
463  } else {
464  gshop.updateIntercept({xip,yip},
465  {facetype[1].const_array(), facetype[0].const_array()},
466  clst, geom);
467  }
468 
469  nmc = build_faces(vbx, cfg, ftx, fty, lst, xip, yip, apx, apy, fcx, fcy,
470  dx, problo, cover_multiple_cuts, nsm);
471 
472  build_cells(vbx, cfg, ftx, fty, apx, apy, dx, vfr, ctr,
473  bar, bct, bnm, lst, small_volfrac, geom, extend_domain_face,
474  nsm, nmc);
475 #endif
476  nsmallcells += nsm;
477  nmulticuts += nmc;
478  }
479  }
480 
481  ParallelAllReduce::Sum<int>({nsmallcells,nmulticuts}, ParallelContext::CommunicatorSub());
482  if (nsmallcells == 0 && nmulticuts == 0) {
483  break;
484  } else {
485  auto ls = m_mgf.getLevelSet();
486  // This is an alias MultiFab, therefore FillBoundary on it is fine.
487  ls.FillBoundary(geom.periodicity());
488  if (amrex::Verbose() > 0) {
489  if (nsmallcells) {
490  amrex::Print() << "AMReX EB: Iter. " << iter+1 << " fixed " << nsmallcells
491  << " small cells" << '\n';
492  }
493  if (nmulticuts) {
494  amrex::Print() << "AMReX EB: Iter. " << iter+1 << " fixed " << nmulticuts
495  << " multicuts" << '\n';
496  }
497  }
498  }
499  }
500 
501  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < maxiter, "EB: failed to fix small cells");
502 
503 #ifdef AMREX_USE_OMP
504 #pragma omp parallel if (Gpu::notInLaunchRegion())
505 #endif
506  for (MFIter mfi(m_mgf); mfi.isValid(); ++mfi)
507  {
508  auto& gfab = m_mgf[mfi];
509  auto const& levelset = gfab.getLevelSet();
510  Array4<Real const> const& clst = levelset.const_array();
511 
512  AMREX_D_TERM(Array4<Real> const& xip = m_edgecent[0].array(mfi);,
513  Array4<Real> const& yip = m_edgecent[1].array(mfi);,
514  Array4<Real> const& zip = m_edgecent[2].array(mfi);)
515 #if (AMREX_SPACEDIM == 3)
516  auto const& edgetype = gfab.getEdgeType();
517  Array4<Type_t const> const& xdg = edgetype[0].const_array();
518  Array4<Type_t const> const& ydg = edgetype[1].const_array();
519  Array4<Type_t const> const& zdg = edgetype[2].const_array();
520 
521  intercept_to_edge_centroid(xip, yip, zip, xdg, ydg, zdg, clst, dx, problo);
522 
523 #elif (AMREX_SPACEDIM == 2)
524  auto& facetype = gfab.getFaceType();
525  Array4<Type_t const> const& ftx = facetype[0].const_array();
526  Array4<Type_t const> const& fty = facetype[1].const_array();
527  // fty then ftx
528  intercept_to_edge_centroid(xip, yip, fty, ftx, clst, dx, problo);
529 #endif
530  }
531 
532  m_levelset = m_mgf.getLevelSet();
533 
534  m_ok = true;
535 }
536 
537 
538 template <typename G>
539 GShopLevel<G>::GShopLevel (IndexSpace const* is, int /*ilev*/, int max_grid_size, int /*ngrow*/,
540  const Geometry& geom, GShopLevel<G>& fineLevel)
541  : Level(is, geom)
542 {
543  if (fineLevel.isAllRegular()) {
544  m_allregular = true;
545  m_ok = true;
546  return;
547  }
548 
549  BL_PROFILE("EB2::GShopLevel()-coarse");
550 
551  const BoxArray& fine_grids = fineLevel.m_grids;
552  const BoxArray& fine_covered_grids = fineLevel.m_covered_grids;
553 
554  const int coarse_ratio = 2;
555  const int min_width = 8;
556  bool coarsenable = fine_grids.coarsenable(coarse_ratio, min_width)
557  && (fine_covered_grids.empty() || fine_covered_grids.coarsenable(coarse_ratio));
558 
559  m_ngrow = amrex::coarsen(fineLevel.m_ngrow,2);
560  if (amrex::scale(m_ngrow,2) != fineLevel.m_ngrow) {
562  }
563 
564  if (coarsenable)
565  {
566  int ierr = coarsenFromFine(fineLevel, true);
567  m_ok = (ierr == 0);
568  }
569  else
570  {
571  Level fine_level_2(is, fineLevel.m_geom);
572  fine_level_2.prepareForCoarsening(fineLevel, max_grid_size, amrex::scale(m_ngrow,2));
573  int ierr = coarsenFromFine(fine_level_2, false);
574  m_ok = (ierr == 0);
575  }
576 }
577 
578 }
579 
580 #endif
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define AMREX_ALWAYS_ASSERT_WITH_MESSAGE(EX, MSG)
Definition: AMReX_BLassert.H:49
amrex::ParmParse pp
Input file parser instance for the given namespace.
Definition: AMReX_HypreIJIface.cpp:15
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
void resize(const Box &b, int N=1, Arena *ar=nullptr)
This function resizes a BaseFab so it covers the Box b with N components.
Definition: AMReX_BaseFab.H:2098
AMREX_FORCE_INLINE Array4< T const > array() const noexcept
Definition: AMReX_BaseFab.H:379
Elixir elixir() noexcept
Definition: AMReX_BaseFab.H:2140
A collection of Boxes stored in an Array.
Definition: AMReX_BoxArray.H:549
bool coarsenable(int refinement_ratio, int min_width=1) const
Coarsen each Box in the BoxArray to the specified ratio.
bool empty() const noexcept
Return whether the BoxArray is empty.
Definition: AMReX_BoxArray.H:602
A class for managing a List of Boxes that share a common IndexType. This class implements operations ...
Definition: AMReX_BoxList.H:52
void clear()
Remove all Boxes from this BoxList.
Vector< Box > & data() noexcept
Returns a reference to the Vector<Box>.
Definition: AMReX_BoxList.H:215
BoxList & coarsen(int ratio)
Coarsen each Box in the BoxList by the ratio.
bool isEmpty() const noexcept
Is this BoxList empty?
Definition: AMReX_BoxList.H:135
BoxList & maxSize(int chunk)
Forces each Box in the BoxList to have sides of length <= chunk.
void swap(BoxList &rhs) noexcept
Definition: AMReX_BoxList.H:219
void push_back(const Box &bn)
Append a Box to this BoxList.
Definition: AMReX_BoxList.H:93
Long size() const noexcept
The number of Boxes in this BoxList.
Definition: AMReX_BoxList.H:113
AMREX_GPU_HOST_DEVICE BoxND & surroundingNodes() noexcept
Convert to NODE type in all directions.
Definition: AMReX_Box.H:946
AMREX_GPU_HOST_DEVICE const IntVectND< dim > & smallEnd() const &noexcept
Get the smallend of the BoxND.
Definition: AMReX_Box.H:105
AMREX_GPU_HOST_DEVICE BoxND & grow(int i) noexcept
Definition: AMReX_Box.H:627
AMREX_GPU_HOST_DEVICE const IntVectND< dim > & bigEnd() const &noexcept
Get the bigend.
Definition: AMReX_Box.H:116
AMREX_GPU_HOST_DEVICE IntVectND< dim > length() const noexcept
Return the length of the BoxND.
Definition: AMReX_Box.H:146
AMREX_GPU_HOST_DEVICE bool coarsenable(const IntVectND< dim > &refrat, const IntVectND< dim > &min_width) const noexcept
Definition: AMReX_Box.H:745
GpuArray< Real, AMREX_SPACEDIM > CellSizeArray() const noexcept
Definition: AMReX_CoordSys.H:76
Calculates the distribution of FABs to MPI processes.
Definition: AMReX_DistributionMapping.H:41
static constexpr int ng
Definition: AMReX_EB2_MultiGFab.H:19
Definition: AMReX_EB2_Level.H:115
GShopLevel(IndexSpace const *is, const Geometry &geom)
Definition: AMReX_EB2_Level.H:135
void define_fine(G const &gshop, const Geometry &geom, int max_grid_size, int ngrow, bool extend_domain_face, int num_crse_opt)
Definition: AMReX_EB2_Level.H:155
static GShopLevel< G > makeAllRegular(IndexSpace const *is, const Geometry &geom)
Definition: AMReX_EB2_Level.H:126
GShopLevel(IndexSpace const *is, int ilev, int max_grid_size, int ngrow, const Geometry &geom, GShopLevel< G > &fineLevel)
Definition: AMReX_EB2_Level.H:539
GShopLevel(IndexSpace const *is, G const &gshop, const Geometry &geom, int max_grid_size, int ngrow, bool extend_domain_face, int num_crse_opt)
Definition: AMReX_EB2_Level.H:140
Definition: AMReX_EB2.H:26
Definition: AMReX_EB2_Level.H:33
MultiFab m_bndryarea
Definition: AMReX_EB2_Level.H:89
bool isOK() const noexcept
Definition: AMReX_EB2_Level.H:37
MultiFab m_volfrac
Definition: AMReX_EB2_Level.H:87
MultiGFab m_mgf
Definition: AMReX_EB2_Level.H:84
bool m_allregular
Definition: AMReX_EB2_Level.H:96
iMultiFab m_cutcellmask
Definition: AMReX_EB2_Level.H:95
IndexSpace const * getEBIndexSpace() const noexcept
Definition: AMReX_EB2_Level.H:63
void prepareForCoarsening(const Level &rhs, int max_grid_size, IntVect const &ngrow)
Definition: AMReX_EB2_Level.cpp:10
void fillCutCellMask(iMultiFab &cutcellmask, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:915
Array< MultiFab, AMREX_SPACEDIM > m_areafrac
Definition: AMReX_EB2_Level.H:92
MultiFab m_centroid
Definition: AMReX_EB2_Level.H:88
IntVect const & nGrowVect() const noexcept
Definition: AMReX_EB2_Level.H:65
FabArray< EBCellFlagFab > m_cellflag
Definition: AMReX_EB2_Level.H:86
Array< MultiFab, AMREX_SPACEDIM > m_facecent
Definition: AMReX_EB2_Level.H:93
BoxArray m_covered_grids
Definition: AMReX_EB2_Level.H:82
void fillLevelSet(MultiFab &levelset, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:880
Geometry m_geom
Definition: AMReX_EB2_Level.H:79
void fillFaceCent(Array< MultiCutFab *, AMREX_SPACEDIM > const &facecent, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:775
IntVect m_ngrow
Definition: AMReX_EB2_Level.H:80
IndexSpace const * m_parent
Definition: AMReX_EB2_Level.H:99
BoxArray m_grids
Definition: AMReX_EB2_Level.H:81
Array< MultiFab, AMREX_SPACEDIM > m_edgecent
Definition: AMReX_EB2_Level.H:94
bool hasEBInfo() const noexcept
Definition: AMReX_EB2_Level.H:69
void write_to_chkpt_file(const std::string &fname, bool extend_domain_face, int max_grid_size) const
Definition: AMReX_EB2_Level.cpp:924
void fillEBCellFlag(FabArray< EBCellFlagFab > &cellflag, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:426
int coarsenFromFine(Level &fineLevel, bool fill_boundary)
Definition: AMReX_EB2_Level.cpp:84
MultiFab m_bndrynorm
Definition: AMReX_EB2_Level.H:91
bool m_has_eb_info
Definition: AMReX_EB2_Level.H:98
MultiFab m_bndrycent
Definition: AMReX_EB2_Level.H:90
MultiFab m_levelset
Definition: AMReX_EB2_Level.H:85
void setRegularLevel()
Definition: AMReX_EB2_Level.H:105
bool isAllRegular() const noexcept
Definition: AMReX_EB2_Level.H:36
void fillEdgeCent(Array< MultiCutFab *, AMREX_SPACEDIM > const &edgecent, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:813
void fillBndryCent(MultiCutFab &bndrycent, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:581
void fillBndryNorm(MultiCutFab &bndrynorm, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:605
const Geometry & Geom() const noexcept
Definition: AMReX_EB2_Level.H:62
const BoxArray & boxArray() const noexcept
Definition: AMReX_EB2_Level.H:56
void fillVolFrac(MultiFab &vfrac, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:477
Level(IndexSpace const *is, const Geometry &geom)
Definition: AMReX_EB2_Level.H:59
const DistributionMapping & DistributionMap() const noexcept
Definition: AMReX_EB2_Level.H:57
friend class GShopLevel
Definition: AMReX_EB2_Level.H:102
bool m_ok
Definition: AMReX_EB2_Level.H:97
void fillCentroid(MultiCutFab &centroid, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:534
DistributionMapping m_dmap
Definition: AMReX_EB2_Level.H:83
void fillBndryArea(MultiCutFab &bndryarea, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:558
void buildCellFlag()
Definition: AMReX_EB2_Level.cpp:400
void buildCutCellMask(Level const &fine_level)
Definition: AMReX_EB2_Level.cpp:934
void fillAreaFrac(Array< MultiCutFab *, AMREX_SPACEDIM > const &areafrac, const Geometry &geom) const
Definition: AMReX_EB2_Level.cpp:629
Definition: AMReX_EB2_MultiGFab.H:65
Definition: AMReX_EBCellFlag.H:287
An Array of FortranArrayBox(FAB)-like Objects.
Definition: AMReX_FabArray.H:344
Rectangular problem domain geometry.
Definition: AMReX_Geometry.H:73
GpuArray< Real, AMREX_SPACEDIM > ProbLoArray() const noexcept
Definition: AMReX_Geometry.H:186
Periodicity periodicity() const noexcept
Definition: AMReX_Geometry.H:355
const Box & Domain() const noexcept
Returns our rectangular domain.
Definition: AMReX_Geometry.H:210
bool isPeriodic(int dir) const noexcept
Is the domain periodic in the specified direction?
Definition: AMReX_Geometry.H:331
Definition: AMReX_GpuElixir.H:13
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE IntVectND< dim > TheZeroVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition: AMReX_IntVect.H:672
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE IntVectND< dim > TheDimensionVector(int d) noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition: AMReX_IntVect.H:691
Definition: AMReX_MFIter.H:57
bool isValid() const noexcept
Is the iterator valid i.e. is it associated with a FAB?
Definition: AMReX_MFIter.H:141
Definition: AMReX_MultiCutFab.H:81
A collection (stored as an array) of FArrayBox objects.
Definition: AMReX_MultiFab.H:38
Parse Parameters From Command Line and Input Files.
Definition: AMReX_ParmParse.H:320
int queryAdd(const char *name, T &ref)
If name is found, the value in the ParmParse database will be stored in the ref argument....
Definition: AMReX_ParmParse.H:993
This class provides the user with a few print options.
Definition: AMReX_Print.H:35
Definition: AMReX_iMultiFab.H:32
Definition: AMReX_FabArrayBase.H:32
AMREX_EXPORT int max_grid_size
Definition: AMReX_EB2.cpp:23
int build_faces(Box const &bx, Array4< EBCellFlag > const &cell, Array4< Type_t > const &fx, Array4< Type_t > const &fy, Array4< Real > const &levset, Array4< Real const > const &interx, Array4< Real const > const &intery, Array4< Real > const &apx, Array4< Real > const &apy, Array4< Real > const &fcx, Array4< Real > const &fcy, GpuArray< Real, AMREX_SPACEDIM > const &dx, GpuArray< Real, AMREX_SPACEDIM > const &problo, bool cover_multiple_cuts, int &nsmallfaces) noexcept
Definition: AMReX_EB2_2D_C.cpp:198
AMREX_EXPORT bool extend_domain_face
Definition: AMReX_EB2.cpp:24
void build_cells(Box const &bx, Array4< EBCellFlag > const &cell, Array4< Type_t > const &fx, Array4< Type_t > const &fy, Array4< Real > const &apx, Array4< Real > const &apy, GpuArray< Real, AMREX_SPACEDIM > const &dx, Array4< Real > const &vfrac, Array4< Real > const &vcent, Array4< Real > const &barea, Array4< Real > const &bcent, Array4< Real > const &bnorm, Array4< Real > const &levset, Real small_volfrac, Geometry const &geom, bool extend_domain_face, int &nsmallcells, int const nmulticuts) noexcept
Definition: AMReX_EB2_2D_C.cpp:352
void intercept_to_edge_centroid(AMREX_D_DECL(Array4< Real > const &excent, Array4< Real > const &eycent, Array4< Real > const &ezcent), AMREX_D_DECL(Array4< Type_t const > const &fx, Array4< Type_t const > const &fy, Array4< Type_t const > const &fz), Array4< Real const > const &levset, GpuArray< Real, AMREX_SPACEDIM > const &dx, GpuArray< Real, AMREX_SPACEDIM > const &problo) noexcept
Definition: AMReX_EB2_ND_C.cpp:5
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
bool inLaunchRegion() noexcept
Definition: AMReX_GpuControl.H:86
MPI_Comm CommunicatorSub() noexcept
sub-communicator for current frame
Definition: AMReX_ParallelContext.H:70
int MyProc() noexcept
return the rank number local to the current Parallel Context
Definition: AMReX_ParallelDescriptor.H:125
int NProcs() noexcept
return the number of MPI ranks local to the current Parallel Context
Definition: AMReX_ParallelDescriptor.H:243
@ min
Definition: AMReX_ParallelReduce.H:18
@ max
Definition: AMReX_ParallelReduce.H:17
RunOn
Definition: AMReX_GpuControl.H:69
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > grow(const BoxND< dim > &b, int i) noexcept
Grow BoxND in all directions by given amount.
Definition: AMReX_Box.H:1211
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > convert(const BoxND< dim > &b, const IntVectND< dim > &typ) noexcept
Returns a BoxND with different type.
Definition: AMReX_Box.H:1435
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition: AMReX.H:111
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > surroundingNodes(const BoxND< dim > &b, int dir) noexcept
Returns a BoxND with NODE based coordinates in direction dir that encloses BoxND b....
Definition: AMReX_Box.H:1399
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE IntVectND< dim > scale(const IntVectND< dim > &p, int s) noexcept
Returns a IntVectND obtained by multiplying each of the components of this IntVectND by s.
Definition: AMReX_IntVect.H:1006
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > refine(const BoxND< dim > &b, int ref_ratio) noexcept
Definition: AMReX_Box.H:1342
int Verbose() noexcept
Definition: AMReX.cpp:164
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > coarsen(const BoxND< dim > &b, int ref_ratio) noexcept
Coarsen BoxND by given (positive) refinement ratio. NOTE: if type(dir) = CELL centered: lo <- lo/rati...
Definition: AMReX_Box.H:1304
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
void AllGatherBoxes(Vector< Box > &bxs, int n_extra_reserve)
Definition: AMReX_Box.cpp:120
std::array< T, N > Array
Definition: AMReX_Array.H:24
Definition: AMReX_Array4.H:61
FabArray memory allocation information.
Definition: AMReX_FabArray.H:66
MFInfo & SetTag() noexcept
Definition: AMReX_FabArray.H:79