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 <memory>
16 #include <unordered_map>
17 
18 namespace amrex
19 {
20  class BoxArray;
21 
23  [[nodiscard]] BoxArray boxComplement (const Box& b1in, const Box& b2);
24 
26  [[nodiscard]] BoxArray complementIn (const Box& b, const BoxArray& ba);
27 
29  [[nodiscard]] BoxArray intersect (const BoxArray& ba, const Box& b, int ng = 0);
30 
31  [[nodiscard]] BoxArray intersect (const BoxArray& ba, const Box& b, const IntVect& ng);
32 
34  [[nodiscard]] BoxArray intersect (const BoxArray& lhs, const BoxArray& rhs);
35 
37  [[nodiscard]] BoxList intersect (const BoxArray& ba, const BoxList& bl);
38 
39  [[nodiscard]] BoxArray convert (const BoxArray& ba, IndexType typ);
40  [[nodiscard]] BoxArray convert (const BoxArray& ba, const IntVect& typ);
41 
42  [[nodiscard]] BoxArray coarsen (const BoxArray& ba, int ratio);
43  [[nodiscard]] BoxArray coarsen (const BoxArray& ba, const IntVect& ratio);
44 
45  [[nodiscard]] BoxArray refine (const BoxArray& ba, int ratio);
46  [[nodiscard]] BoxArray refine (const BoxArray& ba, const IntVect& ratio);
47 
49  [[nodiscard]] BoxList GetBndryCells (const BoxArray& ba, int ngrow);
50 
52  void readBoxArray (BoxArray& ba, std::istream& s, bool b = false);
53 
55  [[nodiscard]] bool match (const BoxArray& x, const BoxArray& y);
56 
71  [[nodiscard]] BoxArray decompose (Box const& domain, int nboxes,
72  Array<bool,AMREX_SPACEDIM> const& decomp
73  = {AMREX_D_DECL(true,true,true)},
74  bool no_overlap = false);
75 
76 struct BARef
77 {
78  BARef ();
79  explicit BARef (size_t size);
80  explicit BARef (const Box& b);
81  explicit BARef (const BoxList& bl);
82  explicit BARef (BoxList&& bl) noexcept;
83  explicit BARef (std::istream& is);
84  BARef (const BARef& rhs);
85  BARef (BARef&& rhs) = delete;
86  BARef& operator= (const BARef& rhs) = delete;
87  BARef& operator= (BARef&& rhs) = delete;
88 
89  ~BARef ();
90 
91  void define (const Box& bx);
92  void define (const BoxList& bl);
93  void define (BoxList&& bl) noexcept;
94  void define (std::istream& is, int& ndims);
96  void resize (Long n);
97 #ifdef AMREX_MEM_PROFILING
98  void updateMemoryUsage_box (int s);
99  void updateMemoryUsage_hash (int s);
100 #endif
101 
102  [[nodiscard]] inline bool HasHashMap () const {
103  bool r;
104 #ifdef AMREX_USE_OMP
105 #pragma omp atomic read
106 #endif
107  r = has_hashmap;
108  return r;
109  }
110 
111  //
114  //
116  mutable Box bbox;
117 
118  mutable IntVect crsn;
119 
120  using HashType = std::unordered_map< IntVect, std::vector<int>, IntVect::shift_hasher > ;
121  //using HashType = std::map< IntVect,std::vector<int> >;
122 
123  mutable HashType hash;
124 
125  mutable bool has_hashmap = false;
126 
127  static int numboxarrays;
128  static int numboxarrays_hwm;
129  static Long total_box_bytes;
130  static Long total_box_bytes_hwm;
131  static Long total_hash_bytes;
132  static Long total_hash_bytes_hwm;
133 
134  static void Initialize ();
135  static void Finalize ();
136  static bool initialized;
137 };
138 
139 struct BATnull
140 {
141  [[nodiscard]] Box operator() (const Box& bx) const noexcept { return bx; }
142  [[nodiscard]] static constexpr Box coarsen (Box const& a_box) { return a_box; }
143  [[nodiscard]] static constexpr IntVect doiLo () { return IntVect::TheZeroVector(); }
144  [[nodiscard]] static constexpr IntVect doiHi () { return IntVect::TheZeroVector(); }
145  [[nodiscard]] static constexpr IndexType index_type () { return IndexType(); }
146  [[nodiscard]] static constexpr IntVect coarsen_ratio () { return IntVect::TheUnitVector(); }
147 };
148 
150 {
151  explicit BATindexType (IndexType a_typ) : m_typ(a_typ) {}
152  [[nodiscard]] Box operator() (const Box& bx) const noexcept { return amrex::convert(bx,m_typ); }
153  [[nodiscard]] static Box coarsen (Box const& a_box) noexcept { return a_box; }
154  [[nodiscard]] static constexpr IntVect doiLo () { return IntVect::TheZeroVector(); }
155  [[nodiscard]] IntVect doiHi () const noexcept { return m_typ.ixType(); }
156  [[nodiscard]] IndexType index_type () const noexcept { return m_typ; }
157  [[nodiscard]] static constexpr IntVect coarsen_ratio () { return IntVect::TheUnitVector(); }
159 };
160 
162 {
163  explicit BATcoarsenRatio (IntVect const& a_crse_ratio) : m_crse_ratio(a_crse_ratio) {}
164  [[nodiscard]] Box operator() (const Box& bx) const noexcept { return amrex::coarsen(bx,m_crse_ratio); }
165  [[nodiscard]] Box coarsen (Box const& a_box) const noexcept { return amrex::coarsen(a_box,m_crse_ratio); }
166  [[nodiscard]] static constexpr IntVect doiLo () { return IntVect::TheZeroVector(); }
167  [[nodiscard]] static constexpr IntVect doiHi () { return IntVect::TheZeroVector(); }
168  [[nodiscard]] static constexpr IndexType index_type () { return IndexType(); }
169  [[nodiscard]] IntVect coarsen_ratio () const noexcept { return m_crse_ratio; }
171 };
172 
174 {
175  BATindexType_coarsenRatio (IndexType a_typ, IntVect const& a_crse_ratio)
176  : m_typ(a_typ), m_crse_ratio(a_crse_ratio) {}
177 
178  [[nodiscard]] Box operator() (const Box& bx) const noexcept {
180  }
181 
182  [[nodiscard]] Box coarsen (Box const& a_box) const noexcept { return amrex::coarsen(a_box,m_crse_ratio); }
183 
184  [[nodiscard]] static constexpr IntVect doiLo () { return IntVect::TheZeroVector(); }
185  [[nodiscard]] IntVect doiHi () const noexcept { return m_typ.ixType(); }
186 
187  [[nodiscard]] IndexType index_type () const noexcept { return m_typ; }
188  [[nodiscard]] IntVect coarsen_ratio () const noexcept { return m_crse_ratio; }
189 
192 };
193 
195 {
197  int a_in_rad, int a_out_rad, int a_extent_rad)
198  : m_face(a_face), m_typ(a_typ), m_crse_ratio(1)
199  {
200  m_loshft = IntVect(-a_extent_rad);
201  m_hishft = IntVect( a_extent_rad);
202  IntVect nodal = a_typ.ixType();
203  m_hishft += nodal;
204  const int d = a_face.coordDir();
205  if (nodal[d]) {
206  // InterpFaceRegister & SyncRegister in IAMR
207  if (m_face.isLow()) {
208  m_loshft[d] = 0;
209  m_hishft[d] = 0;
210  } else {
211  m_loshft[d] = 1;
212  m_hishft[d] = 1;
213  }
214  m_doilo = IntVect(0);
215  m_doihi = nodal;
216  } else {
217  // BndryRegister
218  if (m_face.isLow()) {
219  m_loshft[d] = nodal[d] - a_out_rad;
220  m_hishft[d] = nodal[d] + a_in_rad - 1;
221  } else {
222  m_loshft[d] = 1 - a_in_rad;
223  m_hishft[d] = a_out_rad;
224  }
225  m_doilo = IntVect(a_extent_rad);
226  m_doihi = IntVect(a_extent_rad);
227  m_doihi += nodal;
228  if (m_face.isLow()) { // domain of influence in index space
229  m_doilo[d] = a_out_rad;
230  m_doihi[d] = 0;
231  } else {
232  m_doilo[d] = 0;
233  m_doihi[d] = a_out_rad;
234  }
235  }
236  }
237 
238  [[nodiscard]] Box operator() (const Box& a_bx) const noexcept {
239  IntVect lo = amrex::coarsen(a_bx.smallEnd(), m_crse_ratio);
240  IntVect hi = amrex::coarsen(a_bx.bigEnd(), m_crse_ratio);
241  const int d = m_face.coordDir();
242  if (m_face.isLow()) {
243  hi[d] = lo[d];
244  } else {
245  lo[d] = hi[d];
246  }
247  lo += m_loshft;
248  hi += m_hishft;
249  return Box(lo,hi,m_typ);
250  }
251 
252  [[nodiscard]] Box coarsen (Box const& a_box) const noexcept { return amrex::coarsen(a_box,m_crse_ratio); }
253 
254  [[nodiscard]] IntVect doiLo () const noexcept { return m_doilo; }
255  [[nodiscard]] IntVect doiHi () const noexcept { return m_doihi; }
256 
257  [[nodiscard]] IndexType index_type () const noexcept { return m_typ; }
258  [[nodiscard]] IntVect coarsen_ratio () const noexcept { return m_crse_ratio; }
259 
260  friend bool operator== (BATbndryReg const& a, BATbndryReg const& b) noexcept {
261  return a.m_face == b.m_face && a.m_typ == b.m_typ && a.m_crse_ratio == b.m_crse_ratio
262  && a.m_loshft == b.m_loshft && a.m_hishft == b.m_hishft
263  && a.m_doilo == b.m_doilo && a.m_doihi == b.m_doihi;
264  }
265 
273 };
274 
276 {
278 
279  union BATOp {
280  BATOp () noexcept
281  : m_null() {}
282  BATOp (IndexType t) noexcept
283  : m_indexType(t) {}
284  BATOp (IntVect const& r) noexcept
285  : m_coarsenRatio(r) {}
286  BATOp (IndexType t, IntVect const& r) noexcept
287  : m_indexType_coarsenRatio(t,r) {}
288  BATOp (Orientation f, IndexType t, int in_rad, int out_rad, int extent_rad) noexcept
289  : m_bndryReg(f,t,in_rad,out_rad,extent_rad) {}
295  };
296 
297  BATransformer () = default;
298 
300  : m_bat_type(t.cellCentered() ? BATType::null : BATType::indexType),
301  m_op (t.cellCentered() ? BATOp() : BATOp(t)) {}
302 
303  BATransformer (Orientation f, IndexType t, int in_rad, int out_rad, int extent_rad)
304  : m_bat_type(BATType::bndryReg),
305  m_op(f,t,in_rad,out_rad,extent_rad) {}
306 
307  [[nodiscard]] Box operator() (Box const& ab) const noexcept {
308  switch (m_bat_type)
309  {
310  case BATType::null:
311  return m_op.m_null(ab);
312  case BATType::indexType:
313  return m_op.m_indexType(ab);
315  return m_op.m_coarsenRatio(ab);
317  return m_op.m_indexType_coarsenRatio(ab);
318  default:
319  return m_op.m_bndryReg(ab);
320  }
321  }
322 
323  [[nodiscard]] Box coarsen (Box const& a_box) const noexcept {
324  switch (m_bat_type)
325  {
326  case BATType::null:
327  return amrex::BATnull::coarsen(a_box);
328  case BATType::indexType:
329  return amrex::BATindexType::coarsen(a_box);
331  return m_op.m_coarsenRatio.coarsen(a_box);
333  return m_op.m_indexType_coarsenRatio.coarsen(a_box);
334  default:
335  return m_op.m_bndryReg.coarsen(a_box);
336  }
337  }
338 
339  [[nodiscard]] IntVect doiLo () const noexcept {
340  switch (m_bat_type)
341  {
342  case BATType::null:
343  return amrex::BATnull::doiLo();
344  case BATType::indexType:
350  default:
351  return m_op.m_bndryReg.doiLo();
352  }
353  }
354 
355  [[nodiscard]] IntVect doiHi () const noexcept {
356  switch (m_bat_type)
357  {
358  case BATType::null:
359  return amrex::BATnull::doiHi();
360  case BATType::indexType:
361  return m_op.m_indexType.doiHi();
366  default:
367  return m_op.m_bndryReg.doiHi();
368  }
369  }
370 
371  [[nodiscard]] IndexType index_type () const noexcept {
372  switch (m_bat_type)
373  {
374  case BATType::null:
376  case BATType::indexType:
377  return m_op.m_indexType.index_type();
382  default:
383  return m_op.m_bndryReg.index_type();
384  }
385  }
386 
387  [[nodiscard]] IntVect coarsen_ratio () const noexcept {
388  switch (m_bat_type)
389  {
390  case BATType::null:
392  case BATType::indexType:
398  default:
399  return m_op.m_bndryReg.coarsen_ratio();
400  }
401  }
402 
403  [[nodiscard]] bool is_null () const noexcept {
404  return m_bat_type == BATType::null;
405  }
406 
407  [[nodiscard]] bool is_simple () const noexcept {
408  return m_bat_type != BATType::bndryReg;
409  }
410 
411  void set_coarsen_ratio (IntVect const& a_ratio) noexcept {
412  switch (m_bat_type)
413  {
414  case BATType::null:
415  {
416  if (a_ratio == IntVect::TheUnitVector()) {
417  return;
418  } else {
420  m_op.m_coarsenRatio.m_crse_ratio = a_ratio;
421  return;
422  }
423  }
424  case BATType::indexType:
425  {
426  if (a_ratio == IntVect::TheUnitVector()) {
427  return;
428  } else {
430  auto t = m_op.m_indexType.m_typ;
433  return;
434  }
435  }
437  {
438  if (a_ratio == IntVect::TheUnitVector()) {
440  return;
441  } else {
442  m_op.m_coarsenRatio.m_crse_ratio = a_ratio;
443  return;
444  }
445  }
447  {
448  if (a_ratio == IntVect::TheUnitVector()) {
451  m_op.m_indexType.m_typ = t;
452  return;
453  } else {
455  return;
456  }
457  }
458  default:
459  {
460  m_op.m_bndryReg.m_crse_ratio = a_ratio;
461  return;
462  }
463  }
464  }
465 
466  void set_index_type (IndexType typ) noexcept {
467  switch (m_bat_type)
468  {
469  case BATType::null:
470  {
471  if (typ.cellCentered()) {
472  return;
473  } else {
475  m_op.m_indexType.m_typ = typ;
476  return;
477  }
478  }
479  case BATType::indexType:
480  {
481  if (typ.cellCentered()) {
483  return;
484  } else {
485  m_op.m_indexType.m_typ = typ;
486  return;
487  }
488  }
490  {
491  if (typ.cellCentered()) {
492  return;
493  } else {
498  return;
499  }
500  }
502  {
503  if (typ.cellCentered()) {
507  return;
508  } else {
510  return;
511  }
512  }
513  default:
514  {
515  m_op.m_bndryReg.m_typ = typ;
516  return;
517  }
518  }
519  }
520 
521  friend bool operator== (BATransformer const& a, BATransformer const& b) noexcept {
522  if (a.m_bat_type != BATType::bndryReg && b.m_bat_type != BATType::bndryReg) {
523  return a.index_type() == b.index_type()
524  && a.coarsen_ratio() == b.coarsen_ratio();
525  } else if (a.m_bat_type == BATType::bndryReg && b.m_bat_type == BATType::bndryReg) {
526  return a.m_op.m_bndryReg == b.m_op.m_bndryReg;
527  } else {
528  return false;
529  }
530  }
531 
534 };
535 
536 // for backward compatibility
538 
539 class MFIter;
540 class AmrMesh;
541 class FabArrayBase;
542 
549 class BoxArray
550 {
551 public:
552 
553  BoxArray () noexcept;
554  BoxArray (const BoxArray& rhs) = default;
555  BoxArray (BoxArray&& rhs) noexcept = default;
556  BoxArray& operator= (BoxArray const& rhs) = default;
557  BoxArray& operator= (BoxArray&& rhs) noexcept = default;
558  ~BoxArray() noexcept = default;
559 
561  explicit BoxArray (const Box& bx);
562 
564  explicit BoxArray (size_t n);
565 
567  BoxArray (const Box* bxvec,
568  int nbox);
569 
571  explicit BoxArray (const BoxList& bl);
572  explicit BoxArray (BoxList&& bl) noexcept;
573 
574  BoxArray (const BoxArray& rhs, const BATransformer& trans);
575 
577 
582  void define (const Box& bx);
587  void define (const BoxList& bl);
588  void define (BoxList&& bl) noexcept;
589 
591  void clear ();
592 
594  void resize (Long len);
595 
597  [[nodiscard]] Long size () const noexcept { return m_ref->m_abox.size(); }
598 
600  [[nodiscard]] Long capacity () const noexcept { return static_cast<Long>(m_ref->m_abox.capacity()); }
601 
603  [[nodiscard]] bool empty () const noexcept { return m_ref->m_abox.empty(); }
604 
606  [[nodiscard]] Long numPts() const noexcept;
607 
609  [[nodiscard]] double d_numPts () const noexcept;
616  int readFrom (std::istream& is);
617 
619  std::ostream& writeOn (std::ostream&) const;
620 
622  [[nodiscard]] bool operator== (const BoxArray& rhs) const noexcept;
623 
625  [[nodiscard]] bool operator!= (const BoxArray& rhs) const noexcept;
626 
627  [[nodiscard]] bool operator== (const Vector<Box>& bv) const noexcept;
628  [[nodiscard]] bool operator!= (const Vector<Box>& bv) const noexcept;
629 
631  [[nodiscard]] bool CellEqual (const BoxArray& rhs) const noexcept;
632 
634  BoxArray& maxSize (int block_size);
635 
636  BoxArray& maxSize (const IntVect& block_size);
637 
641  BoxArray& minmaxSize (const IntVect& min_size, const IntVect& max_size);
642 
644  BoxArray& refine (int refinement_ratio);
645 
647  BoxArray& refine (const IntVect& iv);
648 
650  BoxArray& coarsen (int refinement_ratio);
651 
653  [[nodiscard]] bool coarsenable (int refinement_ratio, int min_width=1) const;
654  [[nodiscard]] bool coarsenable (const IntVect& refinement_ratio, int min_width=1) const;
655  [[nodiscard]] bool coarsenable (const IntVect& refinement_ratio, const IntVect& min_width) const;
656 
658  BoxArray& coarsen (const IntVect& iv);
659 
661  BoxArray& growcoarsen (int n, const IntVect& iv);
662  BoxArray& growcoarsen (IntVect const& ngrow, const IntVect& iv);
663 
665  BoxArray& grow (int n);
666 
668  BoxArray& grow (const IntVect& iv);
673  BoxArray& grow (int idir, int n_cell);
678  BoxArray& growLo (int idir, int n_cell);
683  BoxArray& growHi (int idir, int n_cell);
694 
697 
700 
703 
704  BoxArray& convert (const IntVect& iv);
705 
707  BoxArray& convert (Box (*fp)(const Box&));
708 
710  BoxArray& shift (int dir, int nzones);
711 
713  BoxArray& shift (const IntVect &iv);
714 
716  void set (int i, const Box& ibox);
717 
719  [[nodiscard]] Box operator[] (int index) const noexcept {
720  return m_bat(m_ref->m_abox[index]);
721  }
722 
724  [[nodiscard]] Box operator[] (const MFIter& mfi) const noexcept;
725 
727  [[nodiscard]] Box get (int index) const noexcept { return operator[](index); }
728 
730  [[nodiscard]] Box getCellCenteredBox (int index) const noexcept {
731  return m_bat.coarsen(m_ref->m_abox[index]);
732  }
733 
738  [[nodiscard]] bool ok () const;
739 
741  [[nodiscard]] bool isDisjoint () const;
742 
744  [[nodiscard]] BoxList boxList () const;
745 
747  [[nodiscard]] bool contains (const IntVect& v) const;
748 
753  [[nodiscard]]
754  bool contains (const Box& b, bool assume_disjoint_ba = false,
755  const IntVect& ng = IntVect(0)) const;
756 
760  [[nodiscard]]
761  bool contains (const BoxArray& ba, bool assume_disjoint_ba = false,
762  const IntVect& ng = IntVect(0)) const;
763 
771  [[nodiscard]]
772  bool contains (const BoxArray& ba, Periodicity const& period) const;
773 
775  [[nodiscard]] Box minimalBox () const;
776  [[nodiscard]] Box minimalBox (Long& npts_avg_box) const;
777 
782  [[nodiscard]]
783  bool intersects (const Box& b, int ng = 0) const;
784 
785  [[nodiscard]]
786  bool intersects (const Box& b, const IntVect& ng) const;
787 
789  [[nodiscard]]
790  std::vector< std::pair<int,Box> > intersections (const Box& bx) const;
791 
793  [[nodiscard]]
794  std::vector< std::pair<int,Box> > intersections (const Box& bx, bool first_only, int ng) const;
795 
796  [[nodiscard]]
797  std::vector< std::pair<int,Box> > intersections (const Box& bx, bool first_only, const IntVect& ng) const;
798 
800  void intersections (const Box& bx, std::vector< std::pair<int,Box> >& isects) const;
801 
803  void intersections (const Box& bx, std::vector< std::pair<int,Box> >& isects,
804  bool first_only, int ng) const;
805 
806  void intersections (const Box& bx, std::vector< std::pair<int,Box> >& isects,
807  bool first_only, const IntVect& ng) const;
808 
810  [[nodiscard]] BoxList complementIn (const Box& b) const;
811  void complementIn (BoxList& bl, const Box& b) const;
812 
814  void clear_hash_bin () const;
815 
817  void removeOverlap (bool simplify=true);
818 
820  [[nodiscard]] static bool SameRefs (const BoxArray& lhs, const BoxArray& rhs) { return lhs.m_ref == rhs.m_ref; }
821 
822  struct RefID {
823  RefID () noexcept = default;
824  explicit RefID (BARef* data_) noexcept : data(data_) {}
825  bool operator< (const RefID& rhs) const noexcept { return std::less<>()(data,rhs.data); }
826  bool operator== (const RefID& rhs) const noexcept { return data == rhs.data; }
827  bool operator!= (const RefID& rhs) const noexcept { return data != rhs.data; }
828  friend std::ostream& operator<< (std::ostream& os, const RefID& id);
829  private:
830  BARef* data{nullptr};
831  };
832 
834  [[nodiscard]] RefID getRefID () const noexcept { return RefID { m_ref.get() }; }
835 
837  [[nodiscard]] IndexType ixType () const noexcept { return m_bat.index_type(); }
838 
840  [[nodiscard]] IntVect crseRatio () const noexcept { return m_bat.coarsen_ratio(); }
841 
842  static void Initialize ();
843  static void Finalize ();
844  static bool initialized;
845 
847  void uniqify ();
848 
849  [[nodiscard]] BoxList const& simplified_list () const; // For regular AMR grids only!!!
850  [[nodiscard]] BoxArray simplified () const; // For regular AMR grids only!!!
851 
852  [[nodiscard]] BATransformer const& transformer () const;
853 
854  [[nodiscard]] std::weak_ptr<BARef> getWeakRef () const;
855  [[nodiscard]] std::shared_ptr<BARef> const& getSharedRef () const;
856  std::shared_ptr<BARef>& getSharedRef ();
857 
858  friend class AmrMesh;
859  friend class FabArrayBase;
860 
861 private:
863  void type_update ();
864 
865  [[nodiscard]] BARef::HashType& getHashMap () const;
866 
867  [[nodiscard]] IntVect getDoiLo () const noexcept;
868  [[nodiscard]] IntVect getDoiHi () const noexcept;
869 
872  std::shared_ptr<BARef> m_ref;
873  mutable std::shared_ptr<BoxList> m_simplified_list;
874 };
875 
877 std::ostream& operator<< (std::ostream& os, const BoxArray& ba);
878 
879 std::ostream& operator<< (std::ostream& os, const BoxArray::RefID& id);
880 
881 }
882 
883 #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:550
IndexType ixType() const noexcept
Return index type of this BoxArray.
Definition: AMReX_BoxArray.H:837
Box operator[](int index) const noexcept
Return element index of this BoxArray.
Definition: AMReX_BoxArray.H:719
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:600
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:834
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:870
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:840
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:873
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.
std::shared_ptr< BARef > & getSharedRef()
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:730
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:844
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:820
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:872
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:597
Box get(int index) const noexcept
Return element index of this BoxArray.
Definition: AMReX_BoxArray.H:727
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:603
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::shared_ptr< BARef > const & getSharedRef() const
std::vector< std::pair< int, Box > > intersections(const Box &bx, bool first_only, int ng) const
Return intersections of Box and BoxArray(+ghostcells).
std::weak_ptr< BARef > getWeakRef() const
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:77
static Long total_box_bytes_hwm
Definition: AMReX_BoxArray.H:130
BARef & operator=(const BARef &rhs)=delete
static int numboxarrays
Definition: AMReX_BoxArray.H:127
static Long total_hash_bytes
Definition: AMReX_BoxArray.H:131
BARef(BARef &&rhs)=delete
static void Initialize()
IntVect crsn
Definition: AMReX_BoxArray.H:118
std::unordered_map< IntVect, std::vector< int >, IntVect::shift_hasher > HashType
Definition: AMReX_BoxArray.H:120
static Long total_hash_bytes_hwm
Definition: AMReX_BoxArray.H:132
void resize(Long n)
static bool initialized
Definition: AMReX_BoxArray.H:136
HashType hash
Definition: AMReX_BoxArray.H:123
BARef(const Box &b)
BARef(size_t size)
BARef(const BoxList &bl)
Box bbox
Box hash stuff.
Definition: AMReX_BoxArray.H:116
void define(const BoxList &bl)
BARef(BoxList &&bl) noexcept
BARef(std::istream &is)
static int numboxarrays_hwm
Definition: AMReX_BoxArray.H:128
bool has_hashmap
Definition: AMReX_BoxArray.H:125
Vector< Box > m_abox
The data.
Definition: AMReX_BoxArray.H:113
static Long total_box_bytes
Definition: AMReX_BoxArray.H:129
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:102
Definition: AMReX_BoxArray.H:195
IndexType index_type() const noexcept
Definition: AMReX_BoxArray.H:257
Box operator()(const Box &a_bx) const noexcept
Definition: AMReX_BoxArray.H:238
BATbndryReg(Orientation a_face, IndexType a_typ, int a_in_rad, int a_out_rad, int a_extent_rad)
Definition: AMReX_BoxArray.H:196
friend bool operator==(BATbndryReg const &a, BATbndryReg const &b) noexcept
Definition: AMReX_BoxArray.H:260
IntVect doiHi() const noexcept
Definition: AMReX_BoxArray.H:255
IntVect coarsen_ratio() const noexcept
Definition: AMReX_BoxArray.H:258
IntVect doiLo() const noexcept
Definition: AMReX_BoxArray.H:254
Orientation m_face
Definition: AMReX_BoxArray.H:266
IntVect m_loshft
Definition: AMReX_BoxArray.H:269
IndexType m_typ
Definition: AMReX_BoxArray.H:267
Box coarsen(Box const &a_box) const noexcept
Definition: AMReX_BoxArray.H:252
IntVect m_doihi
Definition: AMReX_BoxArray.H:272
IntVect m_hishft
Definition: AMReX_BoxArray.H:270
IntVect m_doilo
Definition: AMReX_BoxArray.H:271
IntVect m_crse_ratio
Definition: AMReX_BoxArray.H:268
Definition: AMReX_BoxArray.H:162
BATcoarsenRatio(IntVect const &a_crse_ratio)
Definition: AMReX_BoxArray.H:163
static constexpr IntVect doiHi()
Definition: AMReX_BoxArray.H:167
Box coarsen(Box const &a_box) const noexcept
Definition: AMReX_BoxArray.H:165
static constexpr IndexType index_type()
Definition: AMReX_BoxArray.H:168
static constexpr IntVect doiLo()
Definition: AMReX_BoxArray.H:166
IntVect m_crse_ratio
Definition: AMReX_BoxArray.H:170
Box operator()(const Box &bx) const noexcept
Definition: AMReX_BoxArray.H:164
IntVect coarsen_ratio() const noexcept
Definition: AMReX_BoxArray.H:169
Definition: AMReX_BoxArray.H:174
BATindexType_coarsenRatio(IndexType a_typ, IntVect const &a_crse_ratio)
Definition: AMReX_BoxArray.H:175
static constexpr IntVect doiLo()
Definition: AMReX_BoxArray.H:184
IntVect doiHi() const noexcept
Definition: AMReX_BoxArray.H:185
Box operator()(const Box &bx) const noexcept
Definition: AMReX_BoxArray.H:178
IntVect m_crse_ratio
Definition: AMReX_BoxArray.H:191
IntVect coarsen_ratio() const noexcept
Definition: AMReX_BoxArray.H:188
IndexType index_type() const noexcept
Definition: AMReX_BoxArray.H:187
Box coarsen(Box const &a_box) const noexcept
Definition: AMReX_BoxArray.H:182
IndexType m_typ
Definition: AMReX_BoxArray.H:190
Definition: AMReX_BoxArray.H:150
IntVect doiHi() const noexcept
Definition: AMReX_BoxArray.H:155
Box operator()(const Box &bx) const noexcept
Definition: AMReX_BoxArray.H:152
IndexType index_type() const noexcept
Definition: AMReX_BoxArray.H:156
static Box coarsen(Box const &a_box) noexcept
Definition: AMReX_BoxArray.H:153
BATindexType(IndexType a_typ)
Definition: AMReX_BoxArray.H:151
IndexType m_typ
Definition: AMReX_BoxArray.H:158
static constexpr IntVect coarsen_ratio()
Definition: AMReX_BoxArray.H:157
static constexpr IntVect doiLo()
Definition: AMReX_BoxArray.H:154
Definition: AMReX_BoxArray.H:140
static constexpr Box coarsen(Box const &a_box)
Definition: AMReX_BoxArray.H:142
static constexpr IndexType index_type()
Definition: AMReX_BoxArray.H:145
static constexpr IntVect doiHi()
Definition: AMReX_BoxArray.H:144
Box operator()(const Box &bx) const noexcept
Definition: AMReX_BoxArray.H:141
static constexpr IntVect coarsen_ratio()
Definition: AMReX_BoxArray.H:146
static constexpr IntVect doiLo()
Definition: AMReX_BoxArray.H:143
Definition: AMReX_BoxArray.H:276
BATransformer(IndexType t)
Definition: AMReX_BoxArray.H:299
void set_index_type(IndexType typ) noexcept
Definition: AMReX_BoxArray.H:466
BATransformer(Orientation f, IndexType t, int in_rad, int out_rad, int extent_rad)
Definition: AMReX_BoxArray.H:303
void set_coarsen_ratio(IntVect const &a_ratio) noexcept
Definition: AMReX_BoxArray.H:411
IntVect doiLo() const noexcept
Definition: AMReX_BoxArray.H:339
IntVect coarsen_ratio() const noexcept
Definition: AMReX_BoxArray.H:387
bool is_null() const noexcept
Definition: AMReX_BoxArray.H:403
Box operator()(Box const &ab) const noexcept
Definition: AMReX_BoxArray.H:307
BATType
Definition: AMReX_BoxArray.H:277
Box coarsen(Box const &a_box) const noexcept
Definition: AMReX_BoxArray.H:323
friend bool operator==(BATransformer const &a, BATransformer const &b) noexcept
Definition: AMReX_BoxArray.H:521
bool is_simple() const noexcept
Definition: AMReX_BoxArray.H:407
IndexType index_type() const noexcept
Definition: AMReX_BoxArray.H:371
IntVect doiHi() const noexcept
Definition: AMReX_BoxArray.H:355
BATType m_bat_type
Definition: AMReX_BoxArray.H:532
BATOp m_op
Definition: AMReX_BoxArray.H:533
Definition: AMReX_BoxArray.H:822
friend std::ostream & operator<<(std::ostream &os, const RefID &id)
RefID() noexcept=default
bool operator<(const RefID &rhs) const noexcept
Definition: AMReX_BoxArray.H:825
bool operator!=(const RefID &rhs) const noexcept
Definition: AMReX_BoxArray.H:827
bool operator==(const RefID &rhs) const noexcept
Definition: AMReX_BoxArray.H:826
BARef * data
Definition: AMReX_BoxArray.H:830
Definition: AMReX_BoxArray.H:279
BATindexType m_indexType
Definition: AMReX_BoxArray.H:291
BATOp(IndexType t) noexcept
Definition: AMReX_BoxArray.H:282
BATOp() noexcept
Definition: AMReX_BoxArray.H:280
BATOp(IntVect const &r) noexcept
Definition: AMReX_BoxArray.H:284
BATnull m_null
Definition: AMReX_BoxArray.H:290
BATbndryReg m_bndryReg
Definition: AMReX_BoxArray.H:294
BATindexType_coarsenRatio m_indexType_coarsenRatio
Definition: AMReX_BoxArray.H:293
BATOp(Orientation f, IndexType t, int in_rad, int out_rad, int extent_rad) noexcept
Definition: AMReX_BoxArray.H:288
BATcoarsenRatio m_coarsenRatio
Definition: AMReX_BoxArray.H:292
BATOp(IndexType t, IntVect const &r) noexcept
Definition: AMReX_BoxArray.H:286