1 #ifndef AMREX_GEOMETRY_H_
2 #define AMREX_GEOMETRY_H_
3 #include <AMReX_Config.H>
21 class DistributionMapping;
31 Real
CellSize (
int dir)
const noexcept {
return dx[dir]; }
57 Real
dx[AMREX_SPACEDIM];
98 int const* is_per =
nullptr) noexcept;
126 static_cast<int>(
c_sys) };
130 static void Setup (
const RealBox* rb =
nullptr,
int coord = -1,
int const* isper =
nullptr) noexcept;
150 void define (const
Box& dom, const
RealBox* rb =
nullptr,
int coord = -1,
int const* is_per =
nullptr) noexcept;
167 void define (const
Box& dom, const
RealBox& rb,
int coord,
Array<
int,AMREX_SPACEDIM> const& is_per) noexcept;
239 #if AMREX_SPACEDIM == 1
241 auto coord = geomdata.
Coord();
251 Real rl = geomdata.
ProbLo()[0] +
static_cast<Real
>(point[0]) *
dx[0];
252 Real rr = rl +
dx[0];
254 constexpr Real
pi = Real(3.1415926535897932);
255 vol = (4.0_rt / 3.0_rt) *
pi *
dx[0] * (rl * rl + rl * rr + rr * rr);
258 #elif AMREX_SPACEDIM == 2
260 auto coord = geomdata.
Coord();
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];
273 constexpr Real
pi = Real(3.1415926535897932);
274 vol =
pi * (r_l + r_r) *
dx[0] *
dx[1];
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];
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];
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);
296 vol =
dx[0] *
dx[1] *
dx[2];
405 for (
int i = 0; i < AMREX_SPACEDIM; ++i) {
413 for (
int i = 0; i < AMREX_SPACEDIM; ++i) {
#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:549
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
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