Block-Structured AMR Software Framework
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Loading...
Searching...
No Matches
AMReX_ParGDB.H
Go to the documentation of this file.
1#ifndef AMREX_ParGDB_H_
2#define AMREX_ParGDB_H_
3#include <AMReX_Config.H>
4
5#include <AMReX_Array.H>
6#include <AMReX_Vector.H>
7#include <AMReX_Geometry.H>
8#include <AMReX_MultiFab.H>
9
10namespace amrex {
11
13{
14public:
15
16 ParGDBBase () noexcept = default;
17 virtual ~ParGDBBase () = default;
18 ParGDBBase (ParGDBBase const&) noexcept = default;
19 ParGDBBase (ParGDBBase &&) noexcept = default;
20 ParGDBBase& operator= (ParGDBBase const&) noexcept = default;
21 ParGDBBase& operator= (ParGDBBase &&) noexcept = default;
22
23 [[nodiscard]] virtual const Geometry& ParticleGeom (int level) const = 0;
24 [[nodiscard]] virtual const Geometry& Geom (int level) const = 0;
25
26 [[nodiscard]] virtual const Vector<Geometry>& ParticleGeom () const = 0;
27 [[nodiscard]] virtual const Vector<Geometry>& Geom () const = 0;
28
29 [[nodiscard]] virtual const DistributionMapping& ParticleDistributionMap (int level) const = 0;
30 [[nodiscard]] virtual const DistributionMapping& DistributionMap (int level) const = 0;
31
32 [[nodiscard]] virtual const Vector<DistributionMapping>& ParticleDistributionMap () const = 0;
33 [[nodiscard]] virtual const Vector<DistributionMapping>& DistributionMap () const = 0;
34
35 [[nodiscard]] virtual const BoxArray& ParticleBoxArray (int level) const = 0;
36 [[nodiscard]] virtual const BoxArray& boxArray (int level) const = 0;
37
38 [[nodiscard]] virtual const Vector<BoxArray>& ParticleBoxArray () const = 0;
39 [[nodiscard]] virtual const Vector<BoxArray>& boxArray () const = 0;
40
41 virtual void SetParticleBoxArray (int level, const BoxArray& new_ba) = 0;
42 virtual void SetParticleDistributionMap (int level, const DistributionMapping& new_dm) = 0;
43 virtual void SetParticleGeometry (int level, const Geometry& new_geom) = 0;
44
45 virtual void ClearParticleBoxArray (int level) = 0;
46 virtual void ClearParticleDistributionMap (int level) = 0;
47 virtual void ClearParticleGeometry (int level) = 0;
48
49 [[nodiscard]] virtual bool LevelDefined (int level) const = 0;
50 [[nodiscard]] virtual int finestLevel () const = 0;
51 [[nodiscard]] virtual int maxLevel () const = 0;
52
53 [[nodiscard]] virtual IntVect refRatio (int level) const = 0;
54 [[nodiscard]] virtual int MaxRefRatio (int level) const = 0;
55
56 [[nodiscard]] virtual Vector<IntVect> refRatio () const = 0;
57
58 template <class MF>
59 [[nodiscard]] bool OnSameGrids (int level, const MF& mf) const;
60};
61
62//
64class ParGDB
65 : public ParGDBBase
66{
67public:
68
69 ParGDB () = default;
70
71 ParGDB (const Geometry & geom,
72 const DistributionMapping & dmap,
73 const BoxArray & ba);
74
75 ParGDB (const Vector<Geometry> & geom,
76 const Vector<DistributionMapping> & dmap,
77 const Vector<BoxArray> & ba,
78 const Vector<int> & rr);
79
80 ParGDB (const Vector<Geometry> & geom,
81 const Vector<DistributionMapping> & dmap,
82 const Vector<BoxArray> & ba,
83 const Vector<IntVect> & rr);
84
85 [[nodiscard]] const Geometry& ParticleGeom (int level) const override;
86 [[nodiscard]] const Geometry& Geom (int level) const override;
87
88 [[nodiscard]] const Vector<Geometry>& ParticleGeom () const override;
89 [[nodiscard]] const Vector<Geometry>& Geom () const override;
90
92 (int level) const override;
93 [[nodiscard]] const DistributionMapping& DistributionMap
94 (int level) const override;
95
96 [[nodiscard]] const Vector<DistributionMapping>& ParticleDistributionMap () const override;
97 [[nodiscard]] const Vector<DistributionMapping>& DistributionMap () const override;
98
99 [[nodiscard]] const BoxArray& ParticleBoxArray (int level) const override;
100 [[nodiscard]] const BoxArray& boxArray (int level) const override;
101
102 [[nodiscard]] const Vector<BoxArray>& ParticleBoxArray () const override;
103 [[nodiscard]] const Vector<BoxArray>& boxArray () const override;
104
105 void SetParticleBoxArray (int level, const BoxArray& new_ba) override;
106 void SetParticleDistributionMap (int level, const DistributionMapping& new_dm) override;
107 void SetParticleGeometry (int level, const Geometry& new_geom) override;
108
109 void ClearParticleBoxArray (int level) override;
110 void ClearParticleDistributionMap (int level) override;
111 void ClearParticleGeometry (int level) override;
112
113 [[nodiscard]] bool LevelDefined (int level) const override;
114 [[nodiscard]] int finestLevel () const override;
115 [[nodiscard]] int maxLevel () const override;
116
117 [[nodiscard]] IntVect refRatio (int level) const override;
118 [[nodiscard]] int MaxRefRatio (int level) const override;
119
120 [[nodiscard]] Vector<IntVect> refRatio () const override;
121
122protected:
123
129};
130
131template <class MF>
132bool
133ParGDBBase::OnSameGrids (int level, const MF& mf) const
134{
135 return (mf.DistributionMap() == ParticleDistributionMap(level) &&
136 mf.boxArray().CellEqual(ParticleBoxArray(level)));
137}
138
139inline
141 const DistributionMapping & dmap,
142 const BoxArray & ba)
143 :
144 m_geom(1,geom),
145 m_dmap(1,dmap),
146 m_ba(1,ba),
147 m_nlevels(1)
148{ }
149
150inline
152 const Vector<DistributionMapping> & dmap,
153 const Vector<BoxArray> & ba,
154 const Vector<IntVect> & rr)
155 :
156 m_geom(geom),
157 m_dmap(dmap),
158 m_ba(ba),
159 m_rr(rr),
160 m_nlevels(static_cast<int>(ba.size()))
161{ }
162
163inline
165 const Vector<DistributionMapping> & dmap,
166 const Vector<BoxArray> & ba,
167 const Vector<int> & rr)
168 :
169 m_geom(geom),
170 m_dmap(dmap),
171 m_ba(ba),
172 m_nlevels(static_cast<int>(ba.size()))
173{
174 for (int level : rr)
175 {
176 m_rr.push_back(level*IntVect::TheUnitVector());
177 }
178}
179
180inline
181const Geometry&
182ParGDB::Geom (int level) const
183{
184 return m_geom[level];
185}
186
187inline
188const Geometry&
189ParGDB::ParticleGeom (int level) const
190{
191 return m_geom[level];
192}
193
194inline
195const Vector<Geometry>&
197{
198 return m_geom;
199}
200
201inline
202const Vector<Geometry>&
204{
205 return m_geom;
206}
207
208inline
211{
212 return m_dmap[level];
213}
214
215inline
217ParGDB::DistributionMap (int level) const
218{
219 return m_dmap[level];
220}
221
222inline
225{
226 return m_dmap;
227}
228
229inline
232{
233 return m_dmap;
234}
235
236inline
237const BoxArray&
239{
240 return m_ba[level];
241}
242
243inline
244const BoxArray&
245ParGDB::boxArray (int level) const
246{
247 return m_ba[level];
248}
249
250inline
251const Vector<BoxArray>&
253{
254 return m_ba;
255}
256
257inline
258const Vector<BoxArray>&
260{
261 return m_ba;
262}
263
264inline
265void
266ParGDB::SetParticleBoxArray (int level, const BoxArray& new_ba)
267{
268 AMREX_ASSERT(level < m_nlevels);
269 m_ba[level] = new_ba;
270}
271
272inline
273void
275{
276 AMREX_ASSERT(level < m_nlevels);
277 m_dmap[level] = new_dm;
278}
279
280inline
281void
282ParGDB::SetParticleGeometry (int level, const Geometry& new_geom)
283{
284 AMREX_ASSERT(level < m_nlevels);
285 m_geom[level] = new_geom;
286}
287
288inline
289void
291{
292 AMREX_ASSERT(level < m_nlevels);
293 m_ba[level] = BoxArray();
294}
295
296inline
297void
303
304inline
305void
307{
308 AMREX_ASSERT(level < m_nlevels);
309 m_geom[level] = Geometry();
310}
311
312inline
313bool
314ParGDB::LevelDefined (int level) const
315{
316 return (level < m_nlevels);
317}
318
319inline
320int
322{
323 return m_nlevels-1;
324}
325
326inline
327int
329{
330 return m_nlevels-1;
331}
332
333inline
335ParGDB::refRatio (int level) const
336{
337 return m_rr[level]*IntVect::TheUnitVector();
338}
339
340inline
341int
342ParGDB::MaxRefRatio (int /*level*/) const
343{
344 int max_ref_ratio = 0;
345 for (int lev = 0; lev < m_nlevels-1; lev++) {
346 max_ref_ratio = std::max(max_ref_ratio, m_rr[lev].max());
347 }
348 return max_ref_ratio;
349}
350
351inline
354{
355 return m_rr;
356}
357
358}
359
360#endif
#define AMREX_ASSERT(EX)
Definition AMReX_BLassert.H:38
A collection of Boxes stored in an Array.
Definition AMReX_BoxArray.H:550
Calculates the distribution of FABs to MPI processes.
Definition AMReX_DistributionMapping.H:41
Rectangular problem domain geometry.
Definition AMReX_Geometry.H:73
AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE constexpr IntVectND< dim > TheUnitVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition AMReX_IntVect.H:680
Definition AMReX_ParGDB.H:13
virtual const Vector< Geometry > & ParticleGeom() const =0
virtual void SetParticleDistributionMap(int level, const DistributionMapping &new_dm)=0
virtual void SetParticleGeometry(int level, const Geometry &new_geom)=0
virtual Vector< IntVect > refRatio() const =0
virtual int maxLevel() const =0
virtual const Vector< DistributionMapping > & DistributionMap() const =0
ParGDBBase() noexcept=default
virtual const Vector< DistributionMapping > & ParticleDistributionMap() const =0
virtual bool LevelDefined(int level) const =0
virtual void ClearParticleDistributionMap(int level)=0
bool OnSameGrids(int level, const MF &mf) const
Definition AMReX_ParGDB.H:133
virtual const Vector< BoxArray > & boxArray() const =0
virtual const Vector< Geometry > & Geom() const =0
virtual void ClearParticleGeometry(int level)=0
virtual int finestLevel() const =0
virtual void SetParticleBoxArray(int level, const BoxArray &new_ba)=0
virtual void ClearParticleBoxArray(int level)=0
virtual const Vector< BoxArray > & ParticleBoxArray() const =0
virtual int MaxRefRatio(int level) const =0
we use this for non-Amr particle code
Definition AMReX_ParGDB.H:66
void ClearParticleBoxArray(int level) override
Definition AMReX_ParGDB.H:290
void SetParticleBoxArray(int level, const BoxArray &new_ba) override
Definition AMReX_ParGDB.H:266
void SetParticleDistributionMap(int level, const DistributionMapping &new_dm) override
Definition AMReX_ParGDB.H:274
const Vector< DistributionMapping > & DistributionMap() const override
Definition AMReX_ParGDB.H:231
void ClearParticleGeometry(int level) override
Definition AMReX_ParGDB.H:306
const Vector< DistributionMapping > & ParticleDistributionMap() const override
Definition AMReX_ParGDB.H:224
void SetParticleGeometry(int level, const Geometry &new_geom) override
Definition AMReX_ParGDB.H:282
int m_nlevels
Definition AMReX_ParGDB.H:128
int MaxRefRatio(int level) const override
Definition AMReX_ParGDB.H:342
const Vector< BoxArray > & ParticleBoxArray() const override
Definition AMReX_ParGDB.H:252
Vector< IntVect > refRatio() const override
Definition AMReX_ParGDB.H:353
const Vector< Geometry > & ParticleGeom() const override
Definition AMReX_ParGDB.H:203
const Vector< BoxArray > & boxArray() const override
Definition AMReX_ParGDB.H:259
int maxLevel() const override
Definition AMReX_ParGDB.H:328
const Vector< Geometry > & Geom() const override
Definition AMReX_ParGDB.H:196
bool LevelDefined(int level) const override
Definition AMReX_ParGDB.H:314
Vector< Geometry > m_geom
Definition AMReX_ParGDB.H:124
void ClearParticleDistributionMap(int level) override
Definition AMReX_ParGDB.H:298
Vector< IntVect > m_rr
Definition AMReX_ParGDB.H:127
int finestLevel() const override
Definition AMReX_ParGDB.H:321
Vector< DistributionMapping > m_dmap
Definition AMReX_ParGDB.H:125
Vector< BoxArray > m_ba
Definition AMReX_ParGDB.H:126
ParGDB()=default
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition AMReX_Vector.H:27
Definition AMReX_Amr.cpp:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr const T & max(const T &a, const T &b) noexcept
Definition AMReX_Algorithm.H:35
const int[]
Definition AMReX_BLProfiler.cpp:1664