Block-Structured AMR Software Framework
AMReX_BoxArray.H
Go to the documentation of this file.
1 
2 #ifndef BL_BOXARRAY_H
3 #define BL_BOXARRAY_H
4 #include <AMReX_Config.H>
5 
6 #include <AMReX_IndexType.H>
7 #include <AMReX_BoxList.H>
8 #include <AMReX_Array.H>
9 #include <AMReX_Periodicity.H>
10 #include <AMReX_Vector.H>
11 
12 #include <iosfwd>
13 #include <cstddef>
14 #include <map>
15 #include <unordered_map>
16 
17 namespace amrex
18 {
19  class BoxArray;
20 
22  [[nodiscard]] BoxArray boxComplement (const Box& b1in, const Box& b2);
23 
25  [[nodiscard]] BoxArray complementIn (const Box& b, const BoxArray& ba);
26 
28  [[nodiscard]] BoxArray intersect (const BoxArray& ba, const Box& b, int ng = 0);
29 
30  [[nodiscard]] BoxArray intersect (const BoxArray& ba, const Box& b, const IntVect& ng);
31 
33  [[nodiscard]] BoxArray intersect (const BoxArray& lhs, const BoxArray& rhs);
34 
36  [[nodiscard]] BoxList intersect (const BoxArray& ba, const BoxList& bl);
37 
38  [[nodiscard]] BoxArray convert (const BoxArray& ba, IndexType typ);
39  [[nodiscard]] BoxArray convert (const BoxArray& ba, const IntVect& typ);
40 
41  [[nodiscard]] BoxArray coarsen (const BoxArray& ba, int ratio);
42  [[nodiscard]] BoxArray coarsen (const BoxArray& ba, const IntVect& ratio);
43 
44  [[nodiscard]] BoxArray refine (const BoxArray& ba, int ratio);
45  [[nodiscard]] BoxArray refine (const BoxArray& ba, const IntVect& ratio);
46 
48  [[nodiscard]] BoxList GetBndryCells (const BoxArray& ba, int ngrow);
49 
51  void readBoxArray (BoxArray& ba, std::istream& s, bool b = false);
52 
54  [[nodiscard]] bool match (const BoxArray& x, const BoxArray& y);
55 
70  [[nodiscard]] BoxArray decompose (Box const& domain, int nboxes,
71  Array<bool,AMREX_SPACEDIM> const& decomp
72  = {AMREX_D_DECL(true,true,true)},
73  bool no_overlap = false);
74 
75 struct BARef
76 {
77  BARef ();
78  explicit BARef (size_t size);
79  explicit BARef (const Box& b);
80  explicit BARef (const BoxList& bl);
81  explicit BARef (BoxList&& bl) noexcept;
82  explicit BARef (std::istream& is);
83  BARef (const BARef& rhs);
84  BARef (BARef&& rhs) = delete;
85  BARef& operator= (const BARef& rhs) = delete;
86  BARef& operator= (BARef&& rhs) = delete;
87 
88  ~BARef ();
89 
90  void define (const Box& bx);
91  void define (const BoxList& bl);
92  void define (BoxList&& bl) noexcept;
93  void define (std::istream& is, int& ndims);
95  void resize (Long n);
96 #ifdef AMREX_MEM_PROFILING
97  void updateMemoryUsage_box (int s);
98  void updateMemoryUsage_hash (int s);
99 #endif
100 
101  [[nodiscard]] inline bool HasHashMap () const {
102  bool r;
103 #ifdef AMREX_USE_OMP
104 #pragma omp atomic read
105 #endif
106  r = has_hashmap;
107  return r;
108  }
109 
110  //
113  //
115  mutable Box bbox;
116 
117  mutable IntVect crsn;
118 
119  using HashType = std::unordered_map< IntVect, std::vector<int>, IntVect::shift_hasher > ;
120  //using HashType = std::map< IntVect,std::vector<int> >;
121 
122  mutable HashType hash;
123 
124  mutable bool has_hashmap = false;
125 
126  static int numboxarrays;
127  static int numboxarrays_hwm;
128  static Long total_box_bytes;
129  static Long total_box_bytes_hwm;
130  static Long total_hash_bytes;
131  static Long total_hash_bytes_hwm;
132 
133  static void Initialize ();
134  static void Finalize ();
135  static bool initialized;
136 };
137 
138 struct BATnull
139 {
140  [[nodiscard]] Box operator() (const Box& bx) const noexcept { return bx; }
141  [[nodiscard]] static constexpr Box coarsen (Box const& a_box) { return a_box; }
142  [[nodiscard]] static constexpr IntVect doiLo () { return IntVect::TheZeroVector(); }
143  [[nodiscard]] static constexpr IntVect doiHi () { return IntVect::TheZeroVector(); }
144  [[nodiscard]] static constexpr IndexType index_type () { return IndexType(); }
145  [[nodiscard]] static constexpr IntVect coarsen_ratio () { return IntVect::TheUnitVector(); }
146 };
147 
149 {
150  explicit BATindexType (IndexType a_typ) : m_typ(a_typ) {}
151  [[nodiscard]] Box operator() (const Box& bx) const noexcept { return amrex::convert(bx,m_typ); }
152  [[nodiscard]] static Box coarsen (Box const& a_box) noexcept { return a_box; }
153  [[nodiscard]] static constexpr IntVect doiLo () { return IntVect::TheZeroVector(); }
154  [[nodiscard]] IntVect doiHi () const noexcept { return m_typ.ixType(); }
155  [[nodiscard]] IndexType index_type () const noexcept { return m_typ; }
156  [[nodiscard]] static constexpr IntVect coarsen_ratio () { return IntVect::TheUnitVector(); }
158 };
159 
161 {
162  explicit BATcoarsenRatio (IntVect const& a_crse_ratio) : m_crse_ratio(a_crse_ratio) {}
163  [[nodiscard]] Box operator() (const Box& bx) const noexcept { return amrex::coarsen(bx,m_crse_ratio); }
164  [[nodiscard]] Box coarsen (Box const& a_box) const noexcept { return amrex::coarsen(a_box,m_crse_ratio); }
165  [[nodiscard]] static constexpr IntVect doiLo () { return IntVect::TheZeroVector(); }
166  [[nodiscard]] static constexpr IntVect doiHi () { return IntVect::TheZeroVector(); }
167  [[nodiscard]] static constexpr IndexType index_type () { return IndexType(); }
168  [[nodiscard]] IntVect coarsen_ratio () const noexcept { return m_crse_ratio; }
170 };
171 
173 {
174  BATindexType_coarsenRatio (IndexType a_typ, IntVect const& a_crse_ratio)
175  : m_typ(a_typ), m_crse_ratio(a_crse_ratio) {}
176 
177  [[nodiscard]] Box operator() (const Box& bx) const noexcept {
179  }
180 
181  [[nodiscard]] Box coarsen (Box const& a_box) const noexcept { return amrex::coarsen(a_box,m_crse_ratio); }
182 
183  [[nodiscard]] static constexpr IntVect doiLo () { return IntVect::TheZeroVector(); }
184  [[nodiscard]] IntVect doiHi () const noexcept { return m_typ.ixType(); }
185 
186  [[nodiscard]] IndexType index_type () const noexcept { return m_typ; }
187  [[nodiscard]] IntVect coarsen_ratio () const noexcept { return m_crse_ratio; }
188 
191 };
192 
194 {
196  int a_in_rad, int a_out_rad, int a_extent_rad)
197  : m_face(a_face), m_typ(a_typ), m_crse_ratio(1)
198  {
199  m_loshft = IntVect(-a_extent_rad);
200  m_hishft = IntVect( a_extent_rad);
201  IntVect nodal = a_typ.ixType();
202  m_hishft += nodal;
203  const int d = a_face.coordDir();
204  if (nodal[d]) {
205  // InterpFaceRegister & SyncRegister in IAMR
206  if (m_face.isLow()) {
207  m_loshft[d] = 0;
208  m_hishft[d] = 0;
209  } else {
210  m_loshft[d] = 1;
211  m_hishft[d] = 1;
212  }
213  m_doilo = IntVect(0);
214  m_doihi = nodal;
215  } else {
216  // BndryRegister
217  if (m_face.isLow()) {
218  m_loshft[d] = nodal[d] - a_out_rad;
219  m_hishft[d] = nodal[d] + a_in_rad - 1;
220  } else {
221  m_loshft[d] = 1 - a_in_rad;
222  m_hishft[d] = a_out_rad;
223  }
224  m_doilo = IntVect(a_extent_rad);
225  m_doihi = IntVect(a_extent_rad);
226  m_doihi += nodal;
227  if (m_face.isLow()) { // domain of influence in index space
228  m_doilo[d] = a_out_rad;
229  m_doihi[d] = 0;
230  } else {
231  m_doilo[d] = 0;
232  m_doihi[d] = a_out_rad;
233  }
234  }
235  }
236 
237  [[nodiscard]] Box operator() (const Box& a_bx) const noexcept {
238  IntVect lo = amrex::coarsen(a_bx.smallEnd(), m_crse_ratio);
239  IntVect hi = amrex::coarsen(a_bx.bigEnd(), m_crse_ratio);
240  const int d = m_face.coordDir();
241  if (m_face.isLow()) {
242  hi[d] = lo[d];
243  } else {
244  lo[d] = hi[d];
245  }
246  lo += m_loshft;
247  hi += m_hishft;
248  return Box(lo,hi,m_typ);
249  }
250 
251  [[nodiscard]] Box coarsen (Box const& a_box) const noexcept { return amrex::coarsen(a_box,m_crse_ratio); }
252 
253  [[nodiscard]] IntVect doiLo () const noexcept { return m_doilo; }
254  [[nodiscard]] IntVect doiHi () const noexcept { return m_doihi; }
255 
256  [[nodiscard]] IndexType index_type () const noexcept { return m_typ; }
257  [[nodiscard]] IntVect coarsen_ratio () const noexcept { return m_crse_ratio; }
258 
259  friend bool operator== (BATbndryReg const& a, BATbndryReg const& b) noexcept {
260  return a.m_face == b.m_face && a.m_typ == b.m_typ && a.m_crse_ratio == b.m_crse_ratio
261  && a.m_loshft == b.m_loshft && a.m_hishft == b.m_hishft
262  && a.m_doilo == b.m_doilo && a.m_doihi == b.m_doihi;
263  }
264 
272 };
273 
275 {
277 
278  union BATOp {
279  BATOp () noexcept
280  : m_null() {}
281  BATOp (IndexType t) noexcept
282  : m_indexType(t) {}
283  BATOp (IntVect const& r) noexcept
284  : m_coarsenRatio(r) {}
285  BATOp (IndexType t, IntVect const& r) noexcept
286  : m_indexType_coarsenRatio(t,r) {}
287  BATOp (Orientation f, IndexType t, int in_rad, int out_rad, int extent_rad) noexcept
288  : m_bndryReg(f,t,in_rad,out_rad,extent_rad) {}
294  };
295 
296  BATransformer () = default;
297 
299  : m_bat_type(t.cellCentered() ? BATType::null : BATType::indexType),
300  m_op (t.cellCentered() ? BATOp() : BATOp(t)) {}
301 
302  BATransformer (Orientation f, IndexType t, int in_rad, int out_rad, int extent_rad)
303  : m_bat_type(BATType::bndryReg),
304  m_op(f,t,in_rad,out_rad,extent_rad) {}
305 
306  [[nodiscard]] Box operator() (Box const& ab) const noexcept {
307  switch (m_bat_type)
308  {
309  case BATType::null:
310  return m_op.m_null(ab);
311  case BATType::indexType:
312  return m_op.m_indexType(ab);
314  return m_op.m_coarsenRatio(ab);
316  return m_op.m_indexType_coarsenRatio(ab);
317  default:
318  return m_op.m_bndryReg(ab);
319  }
320  }
321 
322  [[nodiscard]] Box coarsen (Box const& a_box) const noexcept {
323  switch (m_bat_type)
324  {
325  case BATType::null:
326  return amrex::BATnull::coarsen(a_box);
327  case BATType::indexType:
328  return amrex::BATindexType::coarsen(a_box);
330  return m_op.m_coarsenRatio.coarsen(a_box);
332  return m_op.m_indexType_coarsenRatio.coarsen(a_box);
333  default:
334  return m_op.m_bndryReg.coarsen(a_box);
335  }
336  }
337 
338  [[nodiscard]] IntVect doiLo () const noexcept {
339  switch (m_bat_type)
340  {
341  case BATType::null:
342  return amrex::BATnull::doiLo();
343  case BATType::indexType:
349  default:
350  return m_op.m_bndryReg.doiLo();
351  }
352  }
353 
354  [[nodiscard]] IntVect doiHi () const noexcept {
355  switch (m_bat_type)
356  {
357  case BATType::null:
358  return amrex::BATnull::doiHi();
359  case BATType::indexType:
360  return m_op.m_indexType.doiHi();
365  default:
366  return m_op.m_bndryReg.doiHi();
367  }
368  }
369 
370  [[nodiscard]] IndexType index_type () const noexcept {
371  switch (m_bat_type)
372  {
373  case BATType::null:
375  case BATType::indexType:
376  return m_op.m_indexType.index_type();
381  default:
382  return m_op.m_bndryReg.index_type();
383  }
384  }
385 
386  [[nodiscard]] IntVect coarsen_ratio () const noexcept {
387  switch (m_bat_type)
388  {
389  case BATType::null:
391  case BATType::indexType:
397  default:
398  return m_op.m_bndryReg.coarsen_ratio();
399  }
400  }
401 
402  [[nodiscard]] bool is_null () const noexcept {
403  return m_bat_type == BATType::null;
404  }
405 
406  [[nodiscard]] bool is_simple () const noexcept {
407  return m_bat_type != BATType::bndryReg;
408  }
409 
410  void set_coarsen_ratio (IntVect const& a_ratio) noexcept {
411  switch (m_bat_type)
412  {
413  case BATType::null:
414  {
415  if (a_ratio == IntVect::TheUnitVector()) {
416  return;
417  } else {
419  m_op.m_coarsenRatio.m_crse_ratio = a_ratio;
420  return;
421  }
422  }
423  case BATType::indexType:
424  {
425  if (a_ratio == IntVect::TheUnitVector()) {
426  return;
427  } else {
429  auto t = m_op.m_indexType.m_typ;
432  return;
433  }
434  }
436  {
437  if (a_ratio == IntVect::TheUnitVector()) {
439  return;
440  } else {
441  m_op.m_coarsenRatio.m_crse_ratio = a_ratio;
442  return;
443  }
444  }
446  {
447  if (a_ratio == IntVect::TheUnitVector()) {
450  m_op.m_indexType.m_typ = t;
451  return;
452  } else {
454  return;
455  }
456  }
457  default:
458  {
459  m_op.m_bndryReg.m_crse_ratio = a_ratio;
460  return;
461  }
462  }
463  }
464 
465  void set_index_type (IndexType typ) noexcept {
466  switch (m_bat_type)
467  {
468  case BATType::null:
469  {
470  if (typ.cellCentered()) {
471  return;
472  } else {
474  m_op.m_indexType.m_typ = typ;
475  return;
476  }
477  }
478  case BATType::indexType:
479  {
480  if (typ.cellCentered()) {
482  return;
483  } else {
484  m_op.m_indexType.m_typ = typ;
485  return;
486  }
487  }
489  {
490  if (typ.cellCentered()) {
491  return;
492  } else {
497  return;
498  }
499  }
501  {
502  if (typ.cellCentered()) {
506  return;
507  } else {
509  return;
510  }
511  }
512  default:
513  {
514  m_op.m_bndryReg.m_typ = typ;
515  return;
516  }
517  }
518  }
519 
520  friend bool operator== (BATransformer const& a, BATransformer const& b) noexcept {
521  if (a.m_bat_type != BATType::bndryReg && b.m_bat_type != BATType::bndryReg) {
522  return a.index_type() == b.index_type()
523  && a.coarsen_ratio() == b.coarsen_ratio();
524  } else if (a.m_bat_type == BATType::bndryReg && b.m_bat_type == BATType::bndryReg) {
525  return a.m_op.m_bndryReg == b.m_op.m_bndryReg;
526  } else {
527  return false;
528  }
529  }
530 
533 };
534 
535 // for backward compatibility
537 
538 class MFIter;
539 class AmrMesh;
540 class FabArrayBase;
541 
548 class BoxArray
549 {
550 public:
551 
552  BoxArray () noexcept;
553  BoxArray (const BoxArray& rhs) = default;
554  BoxArray (BoxArray&& rhs) noexcept = default;
555  BoxArray& operator= (BoxArray const& rhs) = default;
556  BoxArray& operator= (BoxArray&& rhs) noexcept = default;
557  ~BoxArray() noexcept = default;
558 
560  explicit BoxArray (const Box& bx);
561 
563  explicit BoxArray (size_t n);
564 
566  BoxArray (const Box* bxvec,
567  int nbox);
568 
570  explicit BoxArray (const BoxList& bl);
571  explicit BoxArray (BoxList&& bl) noexcept;
572 
573  BoxArray (const BoxArray& rhs, const BATransformer& trans);
574 
576 
581  void define (const Box& bx);
586  void define (const BoxList& bl);
587  void define (BoxList&& bl) noexcept;
588 
590  void clear ();
591 
593  void resize (Long len);
594 
596  [[nodiscard]] Long size () const noexcept { return m_ref->m_abox.size(); }
597 
599  [[nodiscard]] Long capacity () const noexcept { return static_cast<Long>(m_ref->m_abox.capacity()); }
600 
602  [[nodiscard]] bool empty () const noexcept { return m_ref->m_abox.empty(); }
603 
605  [[nodiscard]] Long numPts() const noexcept;
606 
608  [[nodiscard]] double d_numPts () const noexcept;
615  int readFrom (std::istream& is);
616 
618  std::ostream& writeOn (std::ostream&) const;
619 
621  [[nodiscard]] bool operator== (const BoxArray& rhs) const noexcept;
622 
624  [[nodiscard]] bool operator!= (const BoxArray& rhs) const noexcept;
625 
626  [[nodiscard]] bool operator== (const Vector<Box>& bv) const noexcept;
627  [[nodiscard]] bool operator!= (const Vector<Box>& bv) const noexcept;
628 
630  [[nodiscard]] bool CellEqual (const BoxArray& rhs) const noexcept;
631 
633  BoxArray& maxSize (int block_size);
634 
635  BoxArray& maxSize (const IntVect& block_size);
636 
640  BoxArray& minmaxSize (const IntVect& min_size, const IntVect& max_size);
641 
643  BoxArray& refine (int refinement_ratio);
644 
646  BoxArray& refine (const IntVect& iv);
647 
649  BoxArray& coarsen (int refinement_ratio);
650 
652  [[nodiscard]] bool coarsenable (int refinement_ratio, int min_width=1) const;
653  [[nodiscard]] bool coarsenable (const IntVect& refinement_ratio, int min_width=1) const;
654  [[nodiscard]] bool coarsenable (const IntVect& refinement_ratio, const IntVect& min_width) const;
655 
657  BoxArray& coarsen (const IntVect& iv);
658 
660  BoxArray& growcoarsen (int n, const IntVect& iv);
661  BoxArray& growcoarsen (IntVect const& ngrow, const IntVect& iv);
662 
664  BoxArray& grow (int n);
665 
667  BoxArray& grow (const IntVect& iv);
672  BoxArray& grow (int idir, int n_cell);
677  BoxArray& growLo (int idir, int n_cell);
682  BoxArray& growHi (int idir, int n_cell);
693 
696 
699 
702 
703  BoxArray& convert (const IntVect& iv);
704 
706  BoxArray& convert (Box (*fp)(const Box&));
707 
709  BoxArray& shift (int dir, int nzones);
710 
712  BoxArray& shift (const IntVect &iv);
713 
715  void set (int i, const Box& ibox);
716 
718  [[nodiscard]] Box operator[] (int index) const noexcept {
719  return m_bat(m_ref->m_abox[index]);
720  }
721 
723  [[nodiscard]] Box operator[] (const MFIter& mfi) const noexcept;
724 
726  [[nodiscard]] Box get (int index) const noexcept { return operator[](index); }
727 
729  [[nodiscard]] Box getCellCenteredBox (int index) const noexcept {
730  return m_bat.coarsen(m_ref->m_abox[index]);
731  }
732 
737  [[nodiscard]] bool ok () const;
738 
740  [[nodiscard]] bool isDisjoint () const;
741 
743  [[nodiscard]] BoxList boxList () const;
744 
746  [[nodiscard]] bool contains (const IntVect& v) const;
747 
752  [[nodiscard]]
753  bool contains (const Box& b, bool assume_disjoint_ba = false,
754  const IntVect& ng = IntVect(0)) const;
755 
759  [[nodiscard]]
760  bool contains (const BoxArray& ba, bool assume_disjoint_ba = false,
761  const IntVect& ng = IntVect(0)) const;
762 
770  [[nodiscard]]
771  bool contains (const BoxArray& ba, Periodicity const& period) const;
772 
774  [[nodiscard]] Box minimalBox () const;
775  [[nodiscard]] Box minimalBox (Long& npts_avg_box) const;
776 
781  [[nodiscard]]
782  bool intersects (const Box& b, int ng = 0) const;
783 
784  [[nodiscard]]
785  bool intersects (const Box& b, const IntVect& ng) const;
786 
788  [[nodiscard]]
789  std::vector< std::pair<int,Box> > intersections (const Box& bx) const;
790 
792  [[nodiscard]]
793  std::vector< std::pair<int,Box> > intersections (const Box& bx, bool first_only, int ng) const;
794 
795  [[nodiscard]]
796  std::vector< std::pair<int,Box> > intersections (const Box& bx, bool first_only, const IntVect& ng) const;
797 
799  void intersections (const Box& bx, std::vector< std::pair<int,Box> >& isects) const;
800 
802  void intersections (const Box& bx, std::vector< std::pair<int,Box> >& isects,
803  bool first_only, int ng) const;
804 
805  void intersections (const Box& bx, std::vector< std::pair<int,Box> >& isects,
806  bool first_only, const IntVect& ng) const;
807 
809  [[nodiscard]] BoxList complementIn (const Box& b) const;
810  void complementIn (BoxList& bl, const Box& b) const;
811 
813  void clear_hash_bin () const;
814 
816  void removeOverlap (bool simplify=true);
817 
819  [[nodiscard]] static bool SameRefs (const BoxArray& lhs, const BoxArray& rhs) { return lhs.m_ref == rhs.m_ref; }
820 
821  struct RefID {
822  RefID () noexcept = default;
823  explicit RefID (BARef* data_) noexcept : data(data_) {}
824  bool operator< (const RefID& rhs) const noexcept { return std::less<>()(data,rhs.data); }
825  bool operator== (const RefID& rhs) const noexcept { return data == rhs.data; }
826  bool operator!= (const RefID& rhs) const noexcept { return data != rhs.data; }
827  friend std::ostream& operator<< (std::ostream& os, const RefID& id);
828  private:
829  BARef* data{nullptr};
830  };
831 
833  [[nodiscard]] RefID getRefID () const noexcept { return RefID { m_ref.get() }; }
834 
836  [[nodiscard]] IndexType ixType () const noexcept { return m_bat.index_type(); }
837 
839  [[nodiscard]] IntVect crseRatio () const noexcept { return m_bat.coarsen_ratio(); }
840 
841  static void Initialize ();
842  static void Finalize ();
843  static bool initialized;
844 
846  void uniqify ();
847 
848  [[nodiscard]] BoxList const& simplified_list () const; // For regular AMR grids only!!!
849  [[nodiscard]] BoxArray simplified () const; // For regular AMR grids only!!!
850 
851  [[nodiscard]] BATransformer const& transformer () const;
852 
853  friend class AmrMesh;
854  friend class FabArrayBase;
855 
856 private:
858  void type_update ();
859 
860  [[nodiscard]] BARef::HashType& getHashMap () const;
861 
862  [[nodiscard]] IntVect getDoiLo () const noexcept;
863  [[nodiscard]] IntVect getDoiHi () const noexcept;
864 
867  std::shared_ptr<BARef> m_ref;
868  mutable std::shared_ptr<BoxList> m_simplified_list;
869 };
870 
872 std::ostream& operator<< (std::ostream& os, const BoxArray& ba);
873 
874 std::ostream& operator<< (std::ostream& os, const BoxArray::RefID& id);
875 
876 }
877 
878 #endif /*BL_BOXARRAY_H*/
int idir
Definition: AMReX_HypreMLABecLap.cpp:1093
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
Definition: AMReX_AmrMesh.H:62
A collection of Boxes stored in an Array.
Definition: AMReX_BoxArray.H:549
IndexType ixType() const noexcept
Return index type of this BoxArray.
Definition: AMReX_BoxArray.H:836
Box operator[](int index) const noexcept
Return element index of this BoxArray.
Definition: AMReX_BoxArray.H:718
BoxArray & minmaxSize(const IntVect &min_size, const IntVect &max_size)
Long capacity() const noexcept
Return the number of boxes that can be held in the current allocated storage.
Definition: AMReX_BoxArray.H:599
BoxArray & growcoarsen(int n, const IntVect &iv)
Grow and then coarsen each Box in the BoxArray.
RefID getRefID() const noexcept
Return a unique ID of the reference.
Definition: AMReX_BoxArray.H:833
bool intersects(const Box &b, int ng=0) const
True if the Box intersects with this BoxArray(+ghostcells). The Box must have the same IndexType as t...
bool intersects(const Box &b, const IntVect &ng) const
BoxArray() noexcept
BoxArray & coarsen(int refinement_ratio)
Coarsen each Box in the BoxArray to the specified ratio.
bool contains(const Box &b, bool assume_disjoint_ba=false, const IntVect &ng=IntVect(0)) const
True if the Box is contained in this BoxArray(+ng). The Box must also have the same IndexType as thos...
std::vector< std::pair< int, Box > > intersections(const Box &bx) const
Return intersections of Box and BoxArray.
bool contains(const IntVect &v) const
True if the IntVect is within any of the Boxes in this BoxArray.
bool contains(const BoxArray &ba, bool assume_disjoint_ba=false, const IntVect &ng=IntVect(0)) const
True if all Boxes in ba are contained in this BoxArray(+ng).
static void Finalize()
BATransformer m_bat
Definition: AMReX_BoxArray.H:865
BoxList boxList() const
Create a BoxList from this BoxArray.
std::vector< std::pair< int, Box > > intersections(const Box &bx, bool first_only, const IntVect &ng) const
BoxArray simplified() const
IntVect crseRatio() const noexcept
Return crse ratio of this BoxArray.
Definition: AMReX_BoxArray.H:839
int readFrom(std::istream &is)
Initialize the BoxArray from the supplied istream. It is an error if the BoxArray has already been in...
IntVect getDoiLo() const noexcept
bool contains(const BoxArray &ba, Periodicity const &period) const
True if all cells in ba are periodically contained in this BoxArray.
void define(const Box &bx)
Initialize the BoxArray from a single box. It is an error if the BoxArray has already been initialize...
BoxArray & growLo(int idir, int n_cell)
Grow each Box in the BoxArray on the low end by n_cell cells in the idir direction.
Long numPts() const noexcept
Returns the total number of cells contained in all boxes in the BoxArray.
std::shared_ptr< BoxList > m_simplified_list
Definition: AMReX_BoxArray.H:868
static void Initialize()
IntVect getDoiHi() const noexcept
BoxList const & simplified_list() const
void clear_hash_bin() const
Clear out the internal hash table used by intersections.
Box minimalBox() const
Return smallest Box that contains all Boxes in this BoxArray.
BoxArray & enclosedCells()
Apply Box::enclosedCells() to each Box in the BoxArray.
bool coarsenable(int refinement_ratio, int min_width=1) const
Coarsen each Box in the BoxArray to the specified ratio.
BoxArray & refine(int refinement_ratio)
Refine each Box in the BoxArray to the specified ratio.
void intersections(const Box &bx, std::vector< std::pair< int, Box > > &isects, bool first_only, int ng) const
intersect Box and BoxArray(+ghostcells), then store the result in isects
void resize(Long len)
Resize the BoxArray. See Vector<T>::resize() for the gory details.
void intersections(const Box &bx, std::vector< std::pair< int, Box > > &isects) const
intersect Box and BoxArray, then store the result in isects
Box getCellCenteredBox(int index) const noexcept
Return cell-centered box at element index of this BoxArray.
Definition: AMReX_BoxArray.H:729
double d_numPts() const noexcept
Returns the total number of cells (in double type) contained in all boxes in the BoxArray.
static bool initialized
Definition: AMReX_BoxArray.H:843
BoxArray & maxSize(int block_size)
Forces each Box in BoxArray to have sides <= block_size.
static bool SameRefs(const BoxArray &lhs, const BoxArray &rhs)
whether two BoxArrays share the same data
Definition: AMReX_BoxArray.H:819
BoxArray & surroundingNodes()
Apply surroundingNodes(Box) to each Box in BoxArray. See the documentation of Box for details.
BoxList complementIn(const Box &b) const
Return box - boxarray.
BATransformer const & transformer() const
void type_update()
Update BoxArray index type according the box type, and then convert boxes to cell-centered.
BoxArray & convert(IndexType typ)
Apply Box::convert(IndexType) to each Box in the BoxArray.
void complementIn(BoxList &bl, const Box &b) const
std::shared_ptr< BARef > m_ref
The data – a reference-counted pointer to a Ref.
Definition: AMReX_BoxArray.H:867
void removeOverlap(bool simplify=true)
Change the BoxArray to one with no overlap and then simplify it (see the simplify function in BoxList...
BARef::HashType & getHashMap() const
std::ostream & writeOn(std::ostream &) const
Output this BoxArray to a checkpoint file.
void intersections(const Box &bx, std::vector< std::pair< int, Box > > &isects, bool first_only, const IntVect &ng) const
Long size() const noexcept
Return the number of boxes in the BoxArray.
Definition: AMReX_BoxArray.H:596
Box get(int index) const noexcept
Return element index of this BoxArray.
Definition: AMReX_BoxArray.H:726
bool CellEqual(const BoxArray &rhs) const noexcept
Are the BoxArrays equal after conversion to cell-centered.
void set(int i, const Box &ibox)
Set element i in this BoxArray to Box ibox.
BoxArray & growHi(int idir, int n_cell)
Grow each Box in the BoxArray on the high end by n_cell cells in the idir direction.
bool ok() const
Return true if Box is valid and they all have the same IndexType. Is true by default if the BoxArray ...
Box minimalBox(Long &npts_avg_box) const
void uniqify()
Make ourselves unique.
bool empty() const noexcept
Return whether the BoxArray is empty.
Definition: AMReX_BoxArray.H:602
bool isDisjoint() const
Return true if set of intersecting Boxes in BoxArray is null.
BoxArray & grow(int n)
Grow each Box in the BoxArray by the specified amount.
std::vector< std::pair< int, Box > > intersections(const Box &bx, bool first_only, int ng) const
Return intersections of Box and BoxArray(+ghostcells).
BoxArray & shift(int dir, int nzones)
Apply Box::shift(int,int) to each Box in the BoxArray.
void clear()
Remove all Boxes from the BoxArray.
A class for managing a List of Boxes that share a common IndexType. This class implements operations ...
Definition: AMReX_BoxList.H:52
Base class for FabArray.
Definition: AMReX_FabArrayBase.H:41
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE CellIndex ixType(int dir) const noexcept
Returns the CellIndex in direction dir.
Definition: AMReX_IndexType.H:116
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE 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:682
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
Definition: AMReX_MFIter.H:57
Encapsulation of the Orientation of the Faces of a Box.
Definition: AMReX_Orientation.H:29
AMREX_GPU_HOST_DEVICE int coordDir() const noexcept
Returns the coordinate direction.
Definition: AMReX_Orientation.H:83
AMREX_GPU_HOST_DEVICE bool isLow() const noexcept
Returns true if Orientation is low.
Definition: AMReX_Orientation.H:89
This provides length of period for periodic domains. 0 means it is not periodic in that direction....
Definition: AMReX_Periodicity.H:17
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition: AMReX_Vector.H:27
AMREX_EXPORT int max_grid_size
Definition: AMReX_EB2.cpp:23
AMREX_GPU_HOST_DEVICE Long size(T const &b) noexcept
integer version
Definition: AMReX_GpuRange.H:26
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
Definition: AMReX_Amr.cpp:49
BoxND< AMREX_SPACEDIM > Box
Definition: AMReX_BaseFwd.H:27
BoxArray intersect(const BoxArray &ba, const Box &b, int ng=0)
Make a BoxArray from the intersection of Box b and BoxArray(+ghostcells).
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
BoxList GetBndryCells(const BoxArray &ba, int ngrow)
Find the ghost cells of a given BoxArray.
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
bool match(const BoxArray &x, const BoxArray &y)
Note that two BoxArrays that match are not necessarily equal.
BoxArray complementIn(const Box &b, const BoxArray &ba)
Make a BoxArray from the complement of BoxArray ba in Box b.
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > refine(const BoxND< dim > &b, int ref_ratio) noexcept
Definition: AMReX_Box.H:1342
BoxArray decompose(Box const &domain, int nboxes, Array< bool, AMREX_SPACEDIM > const &decomp={AMREX_D_DECL(true, true, true)}, bool no_overlap=false)
Decompose domain box into BoxArray.
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
IndexTypeND< AMREX_SPACEDIM > IndexType
Definition: AMReX_BaseFwd.H:33
void readBoxArray(BoxArray &ba, std::istream &s, bool b=false)
Read a BoxArray from a stream. If b is true, read in a special way.
BoxArray boxComplement(const Box &b1in, const Box &b2)
Make a BoxArray from the the complement of b2 in b1in.
std::array< T, N > Array
Definition: AMReX_Array.H:24
Definition: AMReX_BoxArray.H:76
static Long total_box_bytes_hwm
Definition: AMReX_BoxArray.H:129
BARef & operator=(const BARef &rhs)=delete
static int numboxarrays
Definition: AMReX_BoxArray.H:126
static Long total_hash_bytes
Definition: AMReX_BoxArray.H:130
BARef(BARef &&rhs)=delete
static void Initialize()
IntVect crsn
Definition: AMReX_BoxArray.H:117
std::unordered_map< IntVect, std::vector< int >, IntVect::shift_hasher > HashType
Definition: AMReX_BoxArray.H:119
static Long total_hash_bytes_hwm
Definition: AMReX_BoxArray.H:131
void resize(Long n)
static bool initialized
Definition: AMReX_BoxArray.H:135
HashType hash
Definition: AMReX_BoxArray.H:122
BARef(const Box &b)
BARef(size_t size)
BARef(const BoxList &bl)
Box bbox
Box hash stuff.
Definition: AMReX_BoxArray.H:115
void define(const BoxList &bl)
BARef(BoxList &&bl) noexcept
BARef(std::istream &is)
static int numboxarrays_hwm
Definition: AMReX_BoxArray.H:127
bool has_hashmap
Definition: AMReX_BoxArray.H:124
Vector< Box > m_abox
The data.
Definition: AMReX_BoxArray.H:112
static Long total_box_bytes
Definition: AMReX_BoxArray.H:128
void define(const Box &bx)
void define(std::istream &is, int &ndims)
void define(BoxList &&bl) noexcept
BARef(const BARef &rhs)
static void Finalize()
bool HasHashMap() const
Definition: AMReX_BoxArray.H:101
Definition: AMReX_BoxArray.H:194
IndexType index_type() const noexcept
Definition: AMReX_BoxArray.H:256
Box operator()(const Box &a_bx) const noexcept
Definition: AMReX_BoxArray.H:237
BATbndryReg(Orientation a_face, IndexType a_typ, int a_in_rad, int a_out_rad, int a_extent_rad)
Definition: AMReX_BoxArray.H:195
friend bool operator==(BATbndryReg const &a, BATbndryReg const &b) noexcept
Definition: AMReX_BoxArray.H:259
IntVect doiHi() const noexcept
Definition: AMReX_BoxArray.H:254
IntVect coarsen_ratio() const noexcept
Definition: AMReX_BoxArray.H:257
IntVect doiLo() const noexcept
Definition: AMReX_BoxArray.H:253
Orientation m_face
Definition: AMReX_BoxArray.H:265
IntVect m_loshft
Definition: AMReX_BoxArray.H:268
IndexType m_typ
Definition: AMReX_BoxArray.H:266
Box coarsen(Box const &a_box) const noexcept
Definition: AMReX_BoxArray.H:251
IntVect m_doihi
Definition: AMReX_BoxArray.H:271
IntVect m_hishft
Definition: AMReX_BoxArray.H:269
IntVect m_doilo
Definition: AMReX_BoxArray.H:270
IntVect m_crse_ratio
Definition: AMReX_BoxArray.H:267
Definition: AMReX_BoxArray.H:161
BATcoarsenRatio(IntVect const &a_crse_ratio)
Definition: AMReX_BoxArray.H:162
static constexpr IntVect doiHi()
Definition: AMReX_BoxArray.H:166
Box coarsen(Box const &a_box) const noexcept
Definition: AMReX_BoxArray.H:164
static constexpr IndexType index_type()
Definition: AMReX_BoxArray.H:167
static constexpr IntVect doiLo()
Definition: AMReX_BoxArray.H:165
IntVect m_crse_ratio
Definition: AMReX_BoxArray.H:169
Box operator()(const Box &bx) const noexcept
Definition: AMReX_BoxArray.H:163
IntVect coarsen_ratio() const noexcept
Definition: AMReX_BoxArray.H:168
Definition: AMReX_BoxArray.H:173
BATindexType_coarsenRatio(IndexType a_typ, IntVect const &a_crse_ratio)
Definition: AMReX_BoxArray.H:174
static constexpr IntVect doiLo()
Definition: AMReX_BoxArray.H:183
IntVect doiHi() const noexcept
Definition: AMReX_BoxArray.H:184
Box operator()(const Box &bx) const noexcept
Definition: AMReX_BoxArray.H:177
IntVect m_crse_ratio
Definition: AMReX_BoxArray.H:190
IntVect coarsen_ratio() const noexcept
Definition: AMReX_BoxArray.H:187
IndexType index_type() const noexcept
Definition: AMReX_BoxArray.H:186
Box coarsen(Box const &a_box) const noexcept
Definition: AMReX_BoxArray.H:181
IndexType m_typ
Definition: AMReX_BoxArray.H:189
Definition: AMReX_BoxArray.H:149
IntVect doiHi() const noexcept
Definition: AMReX_BoxArray.H:154
Box operator()(const Box &bx) const noexcept
Definition: AMReX_BoxArray.H:151
IndexType index_type() const noexcept
Definition: AMReX_BoxArray.H:155
static Box coarsen(Box const &a_box) noexcept
Definition: AMReX_BoxArray.H:152
BATindexType(IndexType a_typ)
Definition: AMReX_BoxArray.H:150
IndexType m_typ
Definition: AMReX_BoxArray.H:157
static constexpr IntVect coarsen_ratio()
Definition: AMReX_BoxArray.H:156
static constexpr IntVect doiLo()
Definition: AMReX_BoxArray.H:153
Definition: AMReX_BoxArray.H:139
static constexpr Box coarsen(Box const &a_box)
Definition: AMReX_BoxArray.H:141
static constexpr IndexType index_type()
Definition: AMReX_BoxArray.H:144
static constexpr IntVect doiHi()
Definition: AMReX_BoxArray.H:143
Box operator()(const Box &bx) const noexcept
Definition: AMReX_BoxArray.H:140
static constexpr IntVect coarsen_ratio()
Definition: AMReX_BoxArray.H:145
static constexpr IntVect doiLo()
Definition: AMReX_BoxArray.H:142
Definition: AMReX_BoxArray.H:275
BATransformer(IndexType t)
Definition: AMReX_BoxArray.H:298
void set_index_type(IndexType typ) noexcept
Definition: AMReX_BoxArray.H:465
BATransformer(Orientation f, IndexType t, int in_rad, int out_rad, int extent_rad)
Definition: AMReX_BoxArray.H:302
void set_coarsen_ratio(IntVect const &a_ratio) noexcept
Definition: AMReX_BoxArray.H:410
IntVect doiLo() const noexcept
Definition: AMReX_BoxArray.H:338
IntVect coarsen_ratio() const noexcept
Definition: AMReX_BoxArray.H:386
bool is_null() const noexcept
Definition: AMReX_BoxArray.H:402
Box operator()(Box const &ab) const noexcept
Definition: AMReX_BoxArray.H:306
BATType
Definition: AMReX_BoxArray.H:276
Box coarsen(Box const &a_box) const noexcept
Definition: AMReX_BoxArray.H:322
friend bool operator==(BATransformer const &a, BATransformer const &b) noexcept
Definition: AMReX_BoxArray.H:520
bool is_simple() const noexcept
Definition: AMReX_BoxArray.H:406
IndexType index_type() const noexcept
Definition: AMReX_BoxArray.H:370
IntVect doiHi() const noexcept
Definition: AMReX_BoxArray.H:354
BATType m_bat_type
Definition: AMReX_BoxArray.H:531
BATOp m_op
Definition: AMReX_BoxArray.H:532
Definition: AMReX_BoxArray.H:821
friend std::ostream & operator<<(std::ostream &os, const RefID &id)
RefID() noexcept=default
bool operator<(const RefID &rhs) const noexcept
Definition: AMReX_BoxArray.H:824
bool operator!=(const RefID &rhs) const noexcept
Definition: AMReX_BoxArray.H:826
bool operator==(const RefID &rhs) const noexcept
Definition: AMReX_BoxArray.H:825
BARef * data
Definition: AMReX_BoxArray.H:829
Definition: AMReX_BoxArray.H:278
BATindexType m_indexType
Definition: AMReX_BoxArray.H:290
BATOp(IndexType t) noexcept
Definition: AMReX_BoxArray.H:281
BATOp() noexcept
Definition: AMReX_BoxArray.H:279
BATOp(IntVect const &r) noexcept
Definition: AMReX_BoxArray.H:283
BATnull m_null
Definition: AMReX_BoxArray.H:289
BATbndryReg m_bndryReg
Definition: AMReX_BoxArray.H:293
BATindexType_coarsenRatio m_indexType_coarsenRatio
Definition: AMReX_BoxArray.H:292
BATOp(Orientation f, IndexType t, int in_rad, int out_rad, int extent_rad) noexcept
Definition: AMReX_BoxArray.H:287
BATcoarsenRatio m_coarsenRatio
Definition: AMReX_BoxArray.H:291
BATOp(IndexType t, IntVect const &r) noexcept
Definition: AMReX_BoxArray.H:285