Block-Structured AMR Software Framework
AMReX_Geometry.H
Go to the documentation of this file.
1 #ifndef AMREX_GEOMETRY_H_
2 #define AMREX_GEOMETRY_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_Array.H>
6 #include <AMReX_CoordSys.H>
8 #include <AMReX_RealBox.H>
9 #include <AMReX_Periodicity.H>
10 
11 #ifdef AMREX_USE_OMP
12 #include <omp.h>
13 #endif
14 
15 #include <iosfwd>
16 #include <map>
17 
18 namespace amrex {
19 
20 class MultiFab;
21 class DistributionMapping;
22 class BoxArray;
23 
25 {
28  const Real* CellSize () const noexcept { return dx; }
29  //
31  Real CellSize (int dir) const noexcept { return dx[dir]; }
34  const Real* ProbLo () const noexcept { return prob_domain.lo(); }
35  //
37  Real ProbLo (int dir) const noexcept { return prob_domain.lo(dir); }
40  const Real* ProbHi () const noexcept { return prob_domain.hi(); }
41  //
43  Real ProbHi (int dir) const noexcept { return prob_domain.hi(dir); }
46  const Box& Domain () const noexcept { return domain; }
49  int isPeriodic (const int i) const noexcept { return is_periodic[i]; }
52  int Coord () const noexcept { return coord; }
53 
54 public:
57  Real dx[AMREX_SPACEDIM];
58  int is_periodic[AMREX_SPACEDIM];
59  int coord;
60 };
61 
70 class Geometry
71  :
72  public CoordSys
73 {
74 public:
79  Geometry () noexcept;
80 
95  explicit Geometry (const Box& dom,
96  const RealBox* rb = nullptr,
97  int coord = -1,
98  int const* is_per = nullptr) noexcept;
113  Geometry (const Box& dom, const RealBox& rb, int coord,
114  Array<int,AMREX_SPACEDIM> const& is_per) noexcept;
115 
116  ~Geometry () = default;
117  Geometry (const Geometry& rhs) = default;
118  Geometry (Geometry&& rhs) noexcept = default;
119  Geometry& operator= (const Geometry& rhs) = default;
120  Geometry& operator= (Geometry&& rhs) noexcept = default;
121 
123  [[nodiscard]] GeometryData data() const noexcept {
124  return { prob_domain, domain, {AMREX_D_DECL(dx[0],dx[1],dx[2])},
126  static_cast<int>(c_sys) };
127  }
128 
130  static void Setup (const RealBox* rb = nullptr, int coord = -1, int const* isper = nullptr) noexcept;
131  static void ResetDefaultProbDomain (const RealBox& rb) noexcept;
132  static void ResetDefaultPeriodicity (const Array<int,AMREX_SPACEDIM>& is_per) noexcept;
133  static void ResetDefaultCoord (int coord) noexcept;
134 
150  void define (const Box& dom, const RealBox* rb = nullptr, int coord = -1, int const* is_per = nullptr) noexcept;
151 
167  void define (const Box& dom, const RealBox& rb, int coord, Array<int,AMREX_SPACEDIM> const& is_per) noexcept;
168 
170  [[nodiscard]] const RealBox& ProbDomain () const noexcept { return prob_domain; }
172  void ProbDomain (const RealBox& rb) noexcept
173  {
174  prob_domain = rb;
176  }
178  [[nodiscard]] const Real* ProbLo () const noexcept { return prob_domain.lo(); }
180  [[nodiscard]] const Real* ProbHi () const noexcept { return prob_domain.hi(); }
182  [[nodiscard]] Real ProbLo (int dir) const noexcept { return prob_domain.lo(dir); }
184  [[nodiscard]] Real ProbHi (int dir) const noexcept { return prob_domain.hi(dir); }
185 
186  [[nodiscard]] GpuArray<Real,AMREX_SPACEDIM> ProbLoArray () const noexcept {
187  return {{AMREX_D_DECL(prob_domain.lo(0),prob_domain.lo(1),prob_domain.lo(2))}};
188  }
189 
190  [[nodiscard]] GpuArray<Real,AMREX_SPACEDIM> ProbHiArray () const noexcept {
191  return {{AMREX_D_DECL(prob_domain.hi(0),prob_domain.hi(1),prob_domain.hi(2))}};
192  }
193 
195  return roundoff_lo;
196  }
197 
199  return roundoff_hi;
200  }
201 
203  [[nodiscard]] Real ProbSize () const noexcept
204  {
206  }
208  [[nodiscard]] Real ProbLength (int dir) const noexcept { return prob_domain.length(dir); }
210  [[nodiscard]] const Box& Domain () const noexcept { return domain; }
212  void Domain (const Box& bx) noexcept
213  {
214  AMREX_ASSERT(bx.cellCentered());
215  domain = bx;
217  }
219  void GetVolume (MultiFab& vol,
220  const BoxArray& grds,
221  const DistributionMapping& dm,
222  int grow) const;
224  void GetVolume (MultiFab& vol) const;
225 
226  void GetVolume (FArrayBox& vol,
227  const BoxArray& grds,
228  int idx,
229  int grow) const;
230 
233  static Real Volume (const IntVect& point, const GeometryData& geomdata)
234  {
235  const auto *dx = geomdata.CellSize();
236 
237  Real vol;
238 
239 #if AMREX_SPACEDIM == 1
240 
241  auto coord = geomdata.Coord();
242 
243  if (coord == CoordSys::cartesian) {
244  // Cartesian
245 
246  vol = dx[0];
247  }
248  else {
249  // Spherical
250 
251  Real rl = geomdata.ProbLo()[0] + static_cast<Real>(point[0]) * dx[0];
252  Real rr = rl + dx[0];
253 
254  constexpr Real pi = Real(3.1415926535897932);
255  vol = (4.0_rt / 3.0_rt) * pi * dx[0] * (rl * rl + rl * rr + rr * rr);
256  }
257 
258 #elif AMREX_SPACEDIM == 2
259 
260  auto coord = geomdata.Coord();
261 
262  if (coord == CoordSys::cartesian) {
263  // Cartesian
264 
265  vol = dx[0] * dx[1];
266  }
267  else if (coord == CoordSys::RZ) {
268  // Cylindrical
269 
270  Real r_l = geomdata.ProbLo()[0] + static_cast<Real>(point[0]) * dx[0];
271  Real r_r = geomdata.ProbLo()[0] + static_cast<Real>(point[0]+1) * dx[0];
272 
273  constexpr Real pi = Real(3.1415926535897932);
274  vol = pi * (r_l + r_r) * dx[0] * dx[1];
275  }
276  else {
277  // Spherical
278 
279  Real r_l = geomdata.ProbLo()[0] + static_cast<Real>(point[0]) * dx[0];
280  Real r_r = geomdata.ProbLo()[0] + static_cast<Real>(point[0]+1) * dx[0];
281 
282  Real theta_l = geomdata.ProbLo()[1] + static_cast<Real>(point[1]) * dx[1];
283  Real theta_r = geomdata.ProbLo()[1] + static_cast<Real>(point[1]+1) * dx[1];
284 
285  constexpr Real twoThirdsPi = static_cast<Real>(2.0 * 3.1415926535897932 / 3.0);
286  vol = twoThirdsPi * (std::cos(theta_l) - std::cos(theta_r)) * dx[0] *
287  (r_r*r_r + r_r*r_l + r_l*r_l);
288  }
289 
290 #else
291 
292  amrex::ignore_unused(point);
293 
294  // Cartesian
295 
296  vol = dx[0] * dx[1] * dx[2];
297 
298 #endif
299 
300  return vol;
301  }
302 
307  void GetDLogA (MultiFab& dloga,
308  const BoxArray& grds,
309  const DistributionMapping& dm,
310  int dir,
311  int grow) const;
316  void GetFaceArea (MultiFab& area,
317  const BoxArray& grds,
318  const DistributionMapping& dm,
319  int dir,
320  int grow) const;
322  void GetFaceArea (MultiFab& area,
323  int dir) const;
324 
325  void GetFaceArea (FArrayBox& area,
326  const BoxArray& grds,
327  int idx,
328  int dir,
329  int grow) const;
331  [[nodiscard]] bool isPeriodic (int dir) const noexcept { return is_periodic[dir]; }
333  [[nodiscard]] bool isAnyPeriodic () const noexcept
334  {
335  return AMREX_D_TERM(isPeriodic(0),||isPeriodic(1),||isPeriodic(2));
336  }
338  [[nodiscard]] bool isAllPeriodic () const noexcept
339  {
340  return AMREX_D_TERM(isPeriodic(0),&&isPeriodic(1),&&isPeriodic(2));
341  }
342  [[nodiscard]] Array<int,AMREX_SPACEDIM> isPeriodic () const noexcept {
343  return {{AMREX_D_DECL(static_cast<int>(is_periodic[0]),
344  static_cast<int>(is_periodic[1]),
345  static_cast<int>(is_periodic[2]))}};
346  }
347  [[nodiscard]] GpuArray<int,AMREX_SPACEDIM> isPeriodicArray () const noexcept {
348  return {{AMREX_D_DECL(static_cast<int>(is_periodic[0]),
349  static_cast<int>(is_periodic[1]),
350  static_cast<int>(is_periodic[2]))}};
351  }
353  [[nodiscard]] int period (int dir) const noexcept { BL_ASSERT(is_periodic[dir]); return domain.length(dir); }
354 
355  [[nodiscard]] Periodicity periodicity () const noexcept {
357  domain.length(1) * is_periodic[1],
358  domain.length(2) * is_periodic[2])));
359  }
360 
361  [[nodiscard]] Periodicity periodicity (const Box& b) const noexcept {
362  AMREX_ASSERT(b.cellCentered());
363  return Periodicity(IntVect(AMREX_D_DECL(b.length(0) * is_periodic[0],
364  b.length(1) * is_periodic[1],
365  b.length(2) * is_periodic[2])));
366  }
367 
376  void periodicShift (const Box& target,
377  const Box& src,
378  Vector<IntVect>& out) const noexcept;
379 
381  [[nodiscard]] Box growNonPeriodicDomain (IntVect const& ngrow) const noexcept;
383  [[nodiscard]] Box growNonPeriodicDomain (int ngrow) const noexcept;
385  [[nodiscard]] Box growPeriodicDomain (IntVect const& ngrow) const noexcept;
387  [[nodiscard]] Box growPeriodicDomain (int ngrow) const noexcept;
388 
395  is_periodic[1],
396  is_periodic[2])}};
398  is_periodic[1] = period[1];,
399  is_periodic[2] = period[2];);
400  return r;
401  }
402 
403  void coarsen (IntVect const& rr) {
404  domain.coarsen(rr);
405  for (int i = 0; i < AMREX_SPACEDIM; ++i) {
406  dx[i] = (ProbHi(i) - ProbLo(i)) / static_cast<Real>(domain.length(i));
407  inv_dx[i] = Real(1.)/dx[i];
408  }
409  }
410 
411  void refine (IntVect const& rr) {
412  domain.refine(rr);
413  for (int i = 0; i < AMREX_SPACEDIM; ++i) {
414  dx[i] = (ProbHi(i) - ProbLo(i)) / static_cast<Real>(domain.length(i));
415  inv_dx[i] = Real(1.)/dx[i];
416  }
417  }
418 
425  [[nodiscard]] bool outsideRoundoffDomain (AMREX_D_DECL(ParticleReal x, ParticleReal y, ParticleReal z)) const;
426 
433  [[nodiscard]] bool insideRoundoffDomain (AMREX_D_DECL(ParticleReal x, ParticleReal y, ParticleReal z)) const;
434 
439  void computeRoundoffDomain ();
440 
442  [[nodiscard]] GpuArray<ParticleReal, AMREX_SPACEDIM> const&
443  RoundOffLo () const { return roundoff_lo; }
444 
446  [[nodiscard]] GpuArray<ParticleReal, AMREX_SPACEDIM> const&
447  RoundOffHi () const { return roundoff_hi; }
448 
449 private:
450  void read_params ();
451 
452  // is_periodic and RealBox used to be static
453  bool is_periodic[AMREX_SPACEDIM] = {AMREX_D_DECL(false,false,false)};
455 
456  // Due to round-off errors, not all floating point numbers for which plo >= x < phi
457  // will map to a cell that is inside "domain". "roundoff_{lo,hi}" store
458  // positions that are very close to that in prob_domain, and for which all doubles and floats
459  // between them will map to a cell inside domain.
461 
462  //
464 
465  friend std::istream& operator>> (std::istream&, Geometry&);
466 };
467 
468 
470 std::ostream& operator<< (std::ostream&, const Geometry&);
472 std::istream& operator>> (std::istream&, Geometry&);
473 
474 inline
475 Geometry
476 coarsen (Geometry const& fine, IntVect const& rr) {
477  Geometry r{fine};
478  r.coarsen(rr);
479  return r;
480 }
481 
482 inline
483 Geometry
484 coarsen (Geometry const& fine, int rr) { return coarsen(fine, IntVect(rr)); }
485 
486 inline
487 Geometry
488 refine (Geometry const& crse, IntVect const& rr) {
489  Geometry r{crse};
490  r.refine(rr);
491  return r;
492 }
493 
494 inline
495 Geometry
496 refine (Geometry const& crse, int rr) { return refine(crse, IntVect(rr)); }
497 
498 inline
499 const Geometry&
501  return *(AMReX::top()->getDefaultGeometry());
502 }
503 
504 }
505 
506 #endif /*_GEOMETRY_H_*/
#define BL_ASSERT(EX)
Definition: AMReX_BLassert.H:39
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_INLINE
Definition: AMReX_Extension.H:127
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
Array4< Real > fine
Definition: AMReX_InterpFaceRegister.cpp:90
Array4< Real const > crse
Definition: AMReX_InterpFaceRegister.cpp:92
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
Geometry * getDefaultGeometry() noexcept
Definition: AMReX.H:279
static AMReX * top() noexcept
Definition: AMReX.H:268
A collection of Boxes stored in an Array.
Definition: AMReX_BoxArray.H:550
AMREX_GPU_HOST_DEVICE BoxND & refine(int ref_ratio) noexcept
Refine BoxND by given (positive) refinement ratio. NOTE: if type(dir) = CELL centered: lo <- lo*ratio...
Definition: AMReX_Box.H:684
AMREX_GPU_HOST_DEVICE IntVectND< dim > length() const noexcept
Return the length of the BoxND.
Definition: AMReX_Box.H:146
AMREX_GPU_HOST_DEVICE BoxND & coarsen(int ref_ratio) noexcept
Coarsen BoxND by given (positive) refinement ratio. NOTE: if type(dir) = CELL centered: lo <- lo/rati...
Definition: AMReX_Box.H:708
Coordinate System.
Definition: AMReX_CoordSys.H:24
Real inv_dx[AMREX_SPACEDIM]
Definition: AMReX_CoordSys.H:244
Real dx[AMREX_SPACEDIM]
Definition: AMReX_CoordSys.H:242
CoordType c_sys
Definition: AMReX_CoordSys.H:239
@ RZ
Definition: AMReX_CoordSys.H:27
@ cartesian
Definition: AMReX_CoordSys.H:27
Calculates the distribution of FABs to MPI processes.
Definition: AMReX_DistributionMapping.H:41
A Fortran Array of REALs.
Definition: AMReX_FArrayBox.H:229
Rectangular problem domain geometry.
Definition: AMReX_Geometry.H:73
GpuArray< ParticleReal, AMREX_SPACEDIM > roundoff_lo
Definition: AMReX_Geometry.H:460
Geometry(Geometry &&rhs) noexcept=default
const Real * ProbLo() const noexcept
Returns the lo end of the problem domain in each dimension.
Definition: AMReX_Geometry.H:178
Real ProbLo(int dir) const noexcept
Returns the lo end of the problem domain in specified direction.
Definition: AMReX_Geometry.H:182
~Geometry()=default
Array< int, AMREX_SPACEDIM > setPeriodicity(Array< int, AMREX_SPACEDIM > const &period) noexcept
Definition: AMReX_Geometry.H:393
GpuArray< ParticleReal, AMREX_SPACEDIM > const & RoundOffLo() const
Returns roundoff domain's lower end.
Definition: AMReX_Geometry.H:443
void GetVolume(MultiFab &vol, const BoxArray &grds, const DistributionMapping &dm, int grow) const
Define a multifab of areas and volumes with given grow factor.
Definition: AMReX_Geometry.cpp:207
void refine(IntVect const &rr)
Definition: AMReX_Geometry.H:411
void coarsen(IntVect const &rr)
Definition: AMReX_Geometry.H:403
AMREX_GPU_HOST_DEVICE static AMREX_INLINE Real Volume(const IntVect &point, const GeometryData &geomdata)
Return the volume of the specified cell.
Definition: AMReX_Geometry.H:233
GpuArray< Real, AMREX_SPACEDIM > ProbLoArray() const noexcept
Definition: AMReX_Geometry.H:186
friend std::istream & operator>>(std::istream &, Geometry &)
Nice ASCII input.
Definition: AMReX_Geometry.cpp:25
bool isAnyPeriodic() const noexcept
Is domain periodic in any direction?
Definition: AMReX_Geometry.H:333
Array< int, AMREX_SPACEDIM > isPeriodic() const noexcept
Definition: AMReX_Geometry.H:342
static void ResetDefaultProbDomain(const RealBox &rb) noexcept
Definition: AMReX_Geometry.cpp:183
Box domain
Definition: AMReX_Geometry.H:463
void GetFaceArea(MultiFab &area, const BoxArray &grds, const DistributionMapping &dm, int dir, int grow) const
Compute area of cell faces in given region and stuff stuff the results into the passed MultiFab.
Definition: AMReX_Geometry.cpp:298
Real ProbHi(int dir) const noexcept
Returns the hi end of the problem domain in specified direction.
Definition: AMReX_Geometry.H:184
Geometry(const Geometry &rhs)=default
int period(int dir) const noexcept
What's period in specified direction?
Definition: AMReX_Geometry.H:353
const RealBox & ProbDomain() const noexcept
Returns the problem domain.
Definition: AMReX_Geometry.H:170
GpuArray< int, AMREX_SPACEDIM > isPeriodicArray() const noexcept
Definition: AMReX_Geometry.H:347
void Domain(const Box &bx) noexcept
Sets our rectangular domain.
Definition: AMReX_Geometry.H:212
Periodicity periodicity() const noexcept
Definition: AMReX_Geometry.H:355
GpuArray< ParticleReal, AMREX_SPACEDIM > roundoff_hi
Definition: AMReX_Geometry.H:460
bool is_periodic[AMREX_SPACEDIM]
Definition: AMReX_Geometry.H:453
void GetDLogA(MultiFab &dloga, const BoxArray &grds, const DistributionMapping &dm, int dir, int grow) const
Compute d(log(A))/dr at cell centers in given region and stuff the results into the passed MultiFab.
Definition: AMReX_Geometry.cpp:263
Box growNonPeriodicDomain(IntVect const &ngrow) const noexcept
Return domain box with non-periodic directions grown by ngrow.
Definition: AMReX_Geometry.cpp:479
GpuArray< Real, AMREX_SPACEDIM > ProbHiArray() const noexcept
Definition: AMReX_Geometry.H:190
GpuArray< ParticleReal, AMREX_SPACEDIM > ProbLoArrayInParticleReal() const noexcept
Definition: AMReX_Geometry.H:194
const Box & Domain() const noexcept
Returns our rectangular domain.
Definition: AMReX_Geometry.H:210
void periodicShift(const Box &target, const Box &src, Vector< IntVect > &out) const noexcept
Compute Array of shifts which will translate src so that it will intersect target with non-zero inter...
Definition: AMReX_Geometry.cpp:369
const Real * ProbHi() const noexcept
Returns the hi end of the problem domain in each dimension.
Definition: AMReX_Geometry.H:180
static void Setup(const RealBox *rb=nullptr, int coord=-1, int const *isper=nullptr) noexcept
Read static values from ParmParse database.
Definition: AMReX_Geometry.cpp:108
GpuArray< ParticleReal, AMREX_SPACEDIM > ProbHiArrayInParticleReal() const noexcept
Definition: AMReX_Geometry.H:198
Real ProbLength(int dir) const noexcept
Returns length of problem domain in specified dimension.
Definition: AMReX_Geometry.H:208
void computeRoundoffDomain()
Compute the roundoff domain. Public because it contains an extended host / device lambda.
Definition: AMReX_Geometry.cpp:515
bool insideRoundoffDomain(AMREX_D_DECL(ParticleReal x, ParticleReal y, ParticleReal z)) const
Returns true if a point is inside the roundoff domain. All particles with positions inside the roundo...
Definition: AMReX_Geometry.cpp:700
GeometryData data() const noexcept
Returns non-static copy of geometry's stored data.
Definition: AMReX_Geometry.H:123
Geometry() noexcept
The default constructor.
Definition: AMReX_Geometry.cpp:45
void define(const Box &dom, const RealBox *rb=nullptr, int coord=-1, int const *is_per=nullptr) noexcept
Definition: AMReX_Geometry.cpp:70
Box growPeriodicDomain(IntVect const &ngrow) const noexcept
Return domain box with periodic directions grown by ngrow.
Definition: AMReX_Geometry.cpp:491
Geometry & operator=(const Geometry &rhs)=default
GpuArray< ParticleReal, AMREX_SPACEDIM > const & RoundOffHi() const
Returns roundoff domain's higher end.
Definition: AMReX_Geometry.H:447
Real ProbSize() const noexcept
Returns the overall size of the domain by multiplying the ProbLength's together.
Definition: AMReX_Geometry.H:203
void ProbDomain(const RealBox &rb) noexcept
Sets the problem domain.
Definition: AMReX_Geometry.H:172
bool isAllPeriodic() const noexcept
Is domain periodic in all directions?
Definition: AMReX_Geometry.H:338
bool isPeriodic(int dir) const noexcept
Is the domain periodic in the specified direction?
Definition: AMReX_Geometry.H:331
bool outsideRoundoffDomain(AMREX_D_DECL(ParticleReal x, ParticleReal y, ParticleReal z)) const
Returns true if a point is outside the roundoff domain. All particles with positions inside the round...
Definition: AMReX_Geometry.cpp:688
RealBox prob_domain
Definition: AMReX_Geometry.H:454
static void ResetDefaultPeriodicity(const Array< int, AMREX_SPACEDIM > &is_per) noexcept
Definition: AMReX_Geometry.cpp:192
static void ResetDefaultCoord(int coord) noexcept
Definition: AMReX_Geometry.cpp:199
Periodicity periodicity(const Box &b) const noexcept
Definition: AMReX_Geometry.H:361
A collection (stored as an array) of FArrayBox objects.
Definition: AMReX_MultiFab.H:38
This provides length of period for periodic domains. 0 means it is not periodic in that direction....
Definition: AMReX_Periodicity.H:17
A Box with real dimensions. A RealBox is OK iff volume >= 0.
Definition: AMReX_RealBox.H:21
AMREX_GPU_HOST_DEVICE const Real * hi() const &noexcept
Returns hide side.
Definition: AMReX_RealBox.H:51
AMREX_GPU_HOST_DEVICE const Real * lo() const &noexcept
Returns lo side.
Definition: AMReX_RealBox.H:46
AMREX_GPU_HOST_DEVICE Real length(int dir) const noexcept
Returns length in specified direction.
Definition: AMReX_RealBox.H:62
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition: AMReX_Vector.H:27
constexpr std::enable_if_t< std::is_floating_point_v< T >, T > pi()
Definition: AMReX_Math.H:62
Definition: AMReX_Amr.cpp:49
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
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
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 > refine(const BoxND< dim > &b, int ref_ratio) noexcept
Definition: AMReX_Box.H:1342
const Geometry & DefaultGeometry()
Definition: AMReX_Geometry.H:500
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
std::ostream & operator<<(std::ostream &os, AmrMesh const &amr_mesh)
Definition: AMReX_AmrMesh.cpp:1236
std::istream & operator>>(std::istream &is, BoxND< dim > &bx)
Read from istream.
Definition: AMReX_Box.H:1700
std::array< T, N > Array
Definition: AMReX_Array.H:24
Definition: AMReX_Geometry.H:25
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int Coord() const noexcept
Coordinates type.
Definition: AMReX_Geometry.H:52
Box domain
Definition: AMReX_Geometry.H:56
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real CellSize(int dir) const noexcept
Definition: AMReX_Geometry.H:31
int is_periodic[AMREX_SPACEDIM]
Definition: AMReX_Geometry.H:58
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const Real * ProbLo() const noexcept
Returns the lo end of the problem domain in each dimension.
Definition: AMReX_Geometry.H:34
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real ProbHi(int dir) const noexcept
Definition: AMReX_Geometry.H:43
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int isPeriodic(const int i) const noexcept
Returns whether the domain is periodic in the given direction.
Definition: AMReX_Geometry.H:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real ProbLo(int dir) const noexcept
Definition: AMReX_Geometry.H:37
int coord
Definition: AMReX_Geometry.H:59
Real dx[AMREX_SPACEDIM]
Definition: AMReX_Geometry.H:57
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const Box & Domain() const noexcept
Returns our rectangular domain.
Definition: AMReX_Geometry.H:46
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const Real * CellSize() const noexcept
Returns the cellsize for each coordinate direction.
Definition: AMReX_Geometry.H:28
RealBox prob_domain
Definition: AMReX_Geometry.H:55
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const Real * ProbHi() const noexcept
Returns the hi end of the problem domain in each dimension.
Definition: AMReX_Geometry.H:40