Block-Structured AMR Software Framework
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Loading...
Searching...
No Matches
AMReX_ParticleContainerBase.H
Go to the documentation of this file.
1#ifndef AMREX_PARTICLECONTAINERBASE_H_
2#define AMREX_PARTICLECONTAINERBASE_H_
3#include <AMReX_Config.H>
4
5#include <AMReX_Extension.H>
6#include <AMReX_INT.H>
7#include <AMReX_IntVect.H>
8#include <AMReX_ParGDB.H>
9#include <AMReX_Geometry.H>
11#include <AMReX_BoxArray.H>
12#include <AMReX_Vector.H>
13#include <AMReX_ParticleUtil.H>
14#include <AMReX_MultiFab.H>
16#include <AMReX_DenseBins.H>
17
18#include <string>
19
20namespace amrex {
21
23{
24public:
25
27
29 :
30 m_gdb(gdb)
31 {}
32
34 const DistributionMapping & dmap,
35 const BoxArray & ba)
36 :
37 m_gdb_object(std::make_unique<ParGDB>(geom,dmap,ba)),
38 m_gdb(static_cast<ParGDBBase*>(m_gdb_object.get()))
39 {
40 }
41
43 const Vector<DistributionMapping> & dmap,
44 const Vector<BoxArray> & ba,
45 const Vector<int> & rr)
46 :
47 m_gdb_object(std::make_unique<ParGDB>(geom,dmap,ba,rr)),
48 m_gdb(static_cast<ParGDBBase*>(m_gdb_object.get()))
49 {
50 }
51
53 const Vector<DistributionMapping> & dmap,
54 const Vector<BoxArray> & ba,
55 const Vector<IntVect> & rr)
56 :
58 (std::make_unique<ParGDB>(geom,dmap,ba, [&]() -> Vector<int> {
59 Vector<int> ref_ratio;
60 for (auto const& iv : rr)
61 {
62#if AMREX_SPACEDIM > 1
63 AMREX_ASSERT(iv[0] == iv[1]);
64#endif
65#if AMREX_SPACEDIM > 2
66 AMREX_ASSERT(iv[0] == iv[2]);
67#endif
68 ref_ratio.push_back(iv[0]);
69 }
70 return ref_ratio; }() )),
71 m_gdb(static_cast<ParGDBBase*>(m_gdb_object.get()))
72 {
73 }
74
75 virtual ~ParticleContainerBase () = default;
76
79
82
83 void Define (ParGDBBase* gdb) { m_gdb = gdb;}
84
85 void Define (const Geometry & geom,
86 const DistributionMapping & dmap,
87 const BoxArray & ba);
88
89 void Define (const Vector<Geometry> & geom,
90 const Vector<DistributionMapping> & dmap,
91 const Vector<BoxArray> & ba,
92 const Vector<int> & rr);
93
94 void Define (const Vector<Geometry> & geom,
95 const Vector<DistributionMapping> & dmap,
96 const Vector<BoxArray> & ba,
97 const Vector<IntVect> & rr);
98
99 bool isDefined () const { return m_gdb != nullptr; }
100
101 virtual void reserveData ();
102 virtual void resizeData ();
103 void RedefineDummyMF (int lev);
104
105 MFIter MakeMFIter (int lev, const MFItInfo& info) const {
106 AMREX_ASSERT(m_dummy_mf[lev] != nullptr);
107 return MFIter(*m_dummy_mf[lev], info);
108 }
109
110 MFIter MakeMFIter (int lev) const {
111 AMREX_ASSERT(m_dummy_mf[lev] != nullptr);
113 }
114
115 MFIter MakeMFIter (int lev, bool tile) const {
116 AMREX_ASSERT(m_dummy_mf[lev] != nullptr);
117 return MFIter(*m_dummy_mf[lev], tile ? tile_size : IntVect::TheZeroVector());
118 }
119
128 void SetParGDB (const Geometry & geom,
129 const DistributionMapping & dmap,
130 const BoxArray & ba);
131
141 void SetParGDB (const Vector<Geometry> & geom,
142 const Vector<DistributionMapping> & dmap,
143 const Vector<BoxArray> & ba,
144 const Vector<IntVect> & rr);
145
157 void SetParGDB (const Vector<Geometry> & geom,
158 const Vector<DistributionMapping> & dmap,
159 const Vector<BoxArray> & ba,
160 const Vector<int> & rr);
161
169 void SetParticleBoxArray (int lev, BoxArray new_ba);
170
178 void SetParticleDistributionMap (int lev, DistributionMapping new_dmap);
179
187 void SetParticleGeometry (int lev, Geometry new_geom);
188
193 const BoxArray& ParticleBoxArray (int lev) const
194 { return m_gdb->ParticleBoxArray(lev); }
195
201 { return m_gdb->ParticleDistributionMap(lev); }
202
207 const Geometry& Geom (int lev) const { return m_gdb->ParticleGeom(lev); }
208
213 const Geometry& ParticleGeom (int lev) const { return m_gdb->ParticleGeom(lev); }
214
216 int finestLevel () const { return m_gdb->finestLevel(); }
217
219 int maxLevel () const { return m_gdb->maxLevel(); }
220
222 int numLevels() const { return finestLevel() + 1; }
223
225 const ParGDBBase* GetParGDB () const { return m_gdb; }
226
228 ParGDBBase* GetParGDB () { return m_gdb; }
229
230 int Verbose () const { return m_verbose; }
231
233
234 [[nodiscard]] int stableRedistribute () const {return m_stable_redistribute; }
235
236 void setStableRedistribute (int stable) { m_stable_redistribute = stable; }
237
238 const ParticleBufferMap& BufferMap () const {return m_buffer_map;}
239
240 Vector<int> NeighborProcs(int ngrow) const
241 {
242 return computeNeighborProcs(this->GetParGDB(), ngrow);
243 }
244
245 template <class MF>
246 bool OnSameGrids (int level, const MF& mf) const { return m_gdb->OnSameGrids(level, mf); }
247
248 static const std::string& CheckpointVersion ();
249 static const std::string& PlotfileVersion ();
250 static const std::string& DataPrefix ();
251 static int MaxReaders ();
252 static Long MaxParticlesPerRead ();
253 static const std::string& AggregationType ();
254 static int AggregationBuffer ();
255
261
262protected:
263
264 void BuildRedistributeMask (int lev, int nghost=1) const;
265 void defineBufferMap () const;
266
267 int m_verbose{0};
269 std::unique_ptr<ParGDB> m_gdb_object = std::make_unique<ParGDB>();
270 ParGDBBase* m_gdb{nullptr};
272
273 mutable std::unique_ptr<iMultiFab> redistribute_mask_ptr;
274 mutable int redistribute_mask_nghost = std::numeric_limits<int>::min();
277
278};
279
280} // namespace amrex
281
282#endif // AMREX_PARTICLECONTAINERBASE_H_
#define AMREX_ASSERT(EX)
Definition AMReX_BLassert.H:38
#define AMREX_EXPORT
Definition AMReX_Extension.H:191
Definition AMReX_ParticleLocator.H:250
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 > TheZeroVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition AMReX_IntVect.H:670
Definition AMReX_MFIter.H:57
Definition AMReX_ParGDB.H:13
virtual int maxLevel() const =0
virtual const DistributionMapping & ParticleDistributionMap(int level) const =0
virtual const BoxArray & ParticleBoxArray(int level) const =0
bool OnSameGrids(int level, const MF &mf) const
Definition AMReX_ParGDB.H:133
virtual int finestLevel() const =0
virtual const Geometry & ParticleGeom(int level) const =0
we use this for non-Amr particle code
Definition AMReX_ParGDB.H:66
Definition AMReX_ParticleBufferMap.H:53
Definition AMReX_ParticleContainerBase.H:23
int stableRedistribute() const
Definition AMReX_ParticleContainerBase.H:234
ParticleContainerBase & operator=(const ParticleContainerBase &)=delete
amrex::Vector< int > neighbor_procs
Definition AMReX_ParticleContainerBase.H:275
bool isDefined() const
Definition AMReX_ParticleContainerBase.H:99
void BuildRedistributeMask(int lev, int nghost=1) const
Definition AMReX_ParticleContainerBase.cpp:278
std::unique_ptr< iMultiFab > redistribute_mask_ptr
Definition AMReX_ParticleContainerBase.H:273
ParticleContainerBase(const Vector< Geometry > &geom, const Vector< DistributionMapping > &dmap, const Vector< BoxArray > &ba, const Vector< int > &rr)
Definition AMReX_ParticleContainerBase.H:42
const DistributionMapping & ParticleDistributionMap(int lev) const
Get the DistributionMapping for a given level.
Definition AMReX_ParticleContainerBase.H:200
Vector< int > NeighborProcs(int ngrow) const
Definition AMReX_ParticleContainerBase.H:240
ParGDBBase * m_gdb
Definition AMReX_ParticleContainerBase.H:270
void SetVerbose(int verbose)
Definition AMReX_ParticleContainerBase.H:232
const BoxArray & ParticleBoxArray(int lev) const
Get the BoxArray for a given level.
Definition AMReX_ParticleContainerBase.H:193
virtual void reserveData()
Definition AMReX_ParticleContainerBase.cpp:41
ParticleBufferMap m_buffer_map
Definition AMReX_ParticleContainerBase.H:276
MFIter MakeMFIter(int lev, bool tile) const
Definition AMReX_ParticleContainerBase.H:115
void setStableRedistribute(int stable)
Definition AMReX_ParticleContainerBase.H:236
MFIter MakeMFIter(int lev) const
Definition AMReX_ParticleContainerBase.H:110
static const std::string & PlotfileVersion()
Definition AMReX_ParticleContainerBase.cpp:166
static const std::string & AggregationType()
Definition AMReX_ParticleContainerBase.cpp:238
const Geometry & Geom(int lev) const
Get the Geometry for a given level.
Definition AMReX_ParticleContainerBase.H:207
ParGDBBase * GetParGDB()
Get the ParGDB object used to define this container.
Definition AMReX_ParticleContainerBase.H:228
static const std::string & CheckpointVersion()
Definition AMReX_ParticleContainerBase.cpp:150
void Define(ParGDBBase *gdb)
Definition AMReX_ParticleContainerBase.H:83
static AMREX_EXPORT bool memEfficientSort
Definition AMReX_ParticleContainerBase.H:258
ParticleContainerBase(const ParticleContainerBase &)=delete
ParticleContainerBase(const Geometry &geom, const DistributionMapping &dmap, const BoxArray &ba)
Definition AMReX_ParticleContainerBase.H:33
int numLevels() const
the number of defined levels in the ParticleContainer
Definition AMReX_ParticleContainerBase.H:222
static const std::string & DataPrefix()
Definition AMReX_ParticleContainerBase.cpp:181
ParticleContainerBase(ParGDBBase *gdb)
Definition AMReX_ParticleContainerBase.H:28
std::unique_ptr< ParGDB > m_gdb_object
Definition AMReX_ParticleContainerBase.H:269
static Long MaxParticlesPerRead()
Definition AMReX_ParticleContainerBase.cpp:213
int m_verbose
Definition AMReX_ParticleContainerBase.H:267
void SetParticleBoxArray(int lev, BoxArray new_ba)
Set the particle BoxArray. If the container was previously set to to track the AMR hierarchy of an Am...
Definition AMReX_ParticleContainerBase.cpp:112
const ParGDBBase * GetParGDB() const
Get the ParGDB object used to define this container (const version)
Definition AMReX_ParticleContainerBase.H:225
bool OnSameGrids(int level, const MF &mf) const
Definition AMReX_ParticleContainerBase.H:246
virtual ~ParticleContainerBase()=default
static AMREX_EXPORT bool do_tiling
Definition AMReX_ParticleContainerBase.H:256
AmrParticleLocator< DenseBins< Box > > m_particle_locator
Definition AMReX_ParticleContainerBase.H:260
int finestLevel() const
the finest level actually defined for the ParticleContainer
Definition AMReX_ParticleContainerBase.H:216
int redistribute_mask_nghost
Definition AMReX_ParticleContainerBase.H:274
void defineBufferMap() const
Definition AMReX_ParticleContainerBase.cpp:73
static int AggregationBuffer()
Definition AMReX_ParticleContainerBase.cpp:258
void RedefineDummyMF(int lev)
Definition AMReX_ParticleContainerBase.cpp:55
int maxLevel() const
the finest allowed level in the ParticleContainer, whether it is defined or not.
Definition AMReX_ParticleContainerBase.H:219
static AMREX_EXPORT bool use_comms_arena
Definition AMReX_ParticleContainerBase.H:259
Vector< std::unique_ptr< MultiFab > > m_dummy_mf
Definition AMReX_ParticleContainerBase.H:271
MFIter MakeMFIter(int lev, const MFItInfo &info) const
Definition AMReX_ParticleContainerBase.H:105
const ParticleBufferMap & BufferMap() const
Definition AMReX_ParticleContainerBase.H:238
void SetParticleGeometry(int lev, Geometry new_geom)
Set the particle Geometry. If the container was previously set to to track the AMR hierarchy of an Am...
Definition AMReX_ParticleContainerBase.cpp:138
void SetParticleDistributionMap(int lev, DistributionMapping new_dmap)
Set the particle DistributionMapping. If the container was previously set to to track the AMR hierarc...
Definition AMReX_ParticleContainerBase.cpp:125
const Geometry & ParticleGeom(int lev) const
Get the particle Geometry for a given level.
Definition AMReX_ParticleContainerBase.H:213
static int MaxReaders()
Definition AMReX_ParticleContainerBase.cpp:191
void SetParGDB(const Geometry &geom, const DistributionMapping &dmap, const BoxArray &ba)
Set the particle Geometry, DistributionMapping, and BoxArray. If the container was previously set to ...
Definition AMReX_ParticleContainerBase.cpp:83
int Verbose() const
Definition AMReX_ParticleContainerBase.H:230
int m_stable_redistribute
Definition AMReX_ParticleContainerBase.H:268
static AMREX_EXPORT IntVect tile_size
Definition AMReX_ParticleContainerBase.H:257
ParticleContainerBase(ParticleContainerBase &&)=default
ParticleContainerBase(const Vector< Geometry > &geom, const Vector< DistributionMapping > &dmap, const Vector< BoxArray > &ba, const Vector< IntVect > &rr)
Definition AMReX_ParticleContainerBase.H:52
virtual void resizeData()
Definition AMReX_ParticleContainerBase.cpp:46
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 constexpr GpuTupleElement< I, GpuTuple< Ts... > >::type & get(GpuTuple< Ts... > &tup) noexcept
Definition AMReX_Tuple.H:179
Vector< int > computeNeighborProcs(const ParGDBBase *a_gdb, int ngrow)
Definition AMReX_ParticleUtil.cpp:22
const int[]
Definition AMReX_BLProfiler.cpp:1664
int verbose
Definition AMReX_DistributionMapping.cpp:36
Definition AMReX_MFIter.H:20