Chombo + EB  3.0
ProblemDomain.H
Go to the documentation of this file.
1 #ifdef CH_LANG_CC
2 /*
3  * _______ __
4  * / ___/ / ___ __ _ / / ___
5  * / /__/ _ \/ _ \/ V \/ _ \/ _ \
6  * \___/_//_/\___/_/_/_/_.__/\___/
7  * Please refer to Copyright.txt, in Chombo's root directory.
8  */
9 #endif
10 
11 #ifndef _PROBLEMDOMAIN_H_
12 #define _PROBLEMDOMAIN_H_
13 
14 #ifndef WRAPPER
15 #include <iostream>
16 
17 #include "Vector.H"
18 #include "IntVect.H"
19 #include "Box.H"
20 #include "Misc.H"
21 #endif
22 
23 #include "SPACE.H"
24 #include "NamespaceHeader.H"
25 
26 /// Class to manage box-shifting used to enforce periodic BC's
27 /** The ShiftIterator class contains a list of shift vectors necessary
28  to enforce periodic boundary conditions (enforced by shifting the
29  Box of each Fab, then copying any overlapping valid cells with
30  ghost cells of the shifted Fab. Depending on the number of periodic
31  directions, the list of shift vectors over which to iterate differs.
32 */
34 {
35 public:
36  /// Basic constructor
37  ShiftIterator();
38 
39  /// Defining Constructor
40  /**
41  Builds a ShiftIterator based on the periodicity in
42  \em a_isPeriodic
43  */
44  ShiftIterator(const bool* a_isPeriodic);
45 
46  /// Copy constructor
47  ShiftIterator(const ShiftIterator& a_shiftIt);
48 
49  /// Destructor
51 
52  /// Assignment operator
53  ShiftIterator& operator=(const ShiftIterator& a_src);
54 
55  /// Recompute shift vectors based on periodic directions
56  void computeShifts(const bool* a_isPeriodic);
57 
58  /// Returns the current shift unit vector
59  inline IntVect operator()() const;
60 
61  /// Equivalent to the () operator
62  IntVect i() const
63  {
64  return this->operator()();
65  }
66  /// Increment to the next shift unit vector
67  inline void operator++();
68 
69  /// Equivalent to the ++ operator
70  void incr()
71  {
72  ++(*this);
73  }
74 
75  /// Is the iterator still within its range of shift vectors?
76  inline bool ok() const;
77 
78  ///which periodic image is this
79  int index() const
80  {
81  return m_index;
82  }
83 
84  const IntVect& operator[](int index) const
85  {
86  return m_shift_vectors[index];
87  }
88  /// Reset to first shift unit vector
89  inline void reset();
90 
91  /// Equivalent to reset()
92  inline void begin();
93 
94  /// Skip the iterator to the end.
95  /**
96  ok() will return false after this method is called.
97  */
98  void end();
99 
100 private:
101 
102  int m_index;
103 
105 
106 };
107 
108 /// A class to facilitate interaction with physical boundary conditions
109 /**
110 
111  ProblemDomain is a class which facilitates the application of physical
112  boundary conditions, both periodic and non-periodic. This class contains
113  much of the functionality of the Box class, since logically the
114  computational domain is generally a Box.
115 
116  Intersection with a ProblemDomain object will result in only removing
117  regions which are outside the physical domain in non-periodic directions.
118  Regions outside the logical computational domain in periodic directions will
119  be treated as ghost cells which can be filled with an exchange() function
120  or through suitable interpolation from a coarser domain.
121 
122  Since ProblemDomain will contain a Box, it is a dimension dependent class,
123  so SpaceDim must be defined as either 1, 2, or 3 when compiling.
124 
125  Note that this implementation of ProblemDomain is inherently
126  cell-centered.
127 
128 */
129 
131 {
132 public:
133 
134  // Constructors
135 
136  /// The default constructor. The constructed domain box is empty.
137  /**
138  */
139  ProblemDomain ();
140 
141  /// Construct ProblemDomain with \em a_domBox as computational domain
142  /**
143  This constructor defaults to non-periodic domain
144  */
145  ProblemDomain(const Box& a_domBox);
146 
147  /// Construct ProblemDomain with a_domBox as computational domain.
148  /**
149  \em a_isPeriodic is a SpaceDim array of bools; if true, the
150  physical boundary condition is periodic in the coordinate direction.
151  False means a non-periodic BC.
152  */
153  ProblemDomain(const Box& a_domBox, const bool* a_isPeriodic);
154 
155  /// Construct a ProblemDomain.
156  /**
157  It is an error if small is greater than big.
158  Defaults to non-periodic domain.
159  */
160  ProblemDomain (const IntVect& small,
161  const IntVect& big);
162 
163  /// Construct a ProblemDomain.
164  /**
165  It is an error if small is greater than big.
166  \em a_isPeriodic is a SpaceDim array of bools; if true, the
167  physical boundary condition is periodic in the coordinate direction.
168  False means a non-periodic BC.
169  */
170  ProblemDomain (const IntVect& small,
171  const IntVect& big,
172  const bool* a_isPeriodic);
173 
174  /// Construct ProblemDomain with specified lengths.
175  /**
176  It is an error if the lengths are negative.
177  Defaults to non-periodic domain.
178  */
179  ProblemDomain (const IntVect& small,
180  const int* vec_len);
181 
182  /// Construct ProblemDomain with specified lengths.
183  /**
184  It is an error if the lengths are negative.
185  \em a_isPeriodic is a SpaceDim array of bools; if true, the
186  physical boundary condition is periodic in the coordinate direction.
187  False means a non-periodic BC.
188  */
189  ProblemDomain (const IntVect& small,
190  const int* vec_len,
191  const bool* a_isPeriodic);
192 
193  /// The copy constructor.
194  /**
195  */
196  ProblemDomain (const ProblemDomain& a_src);
197 
198  /// Construct ProblemDomain with \em a_domBox as computational domain
199  /**
200  This constructor defaults to non-periodic domain
201  */
202  void define(const Box& a_domBox);
203 
204  /// Construct ProblemDomain with a_domBox as computational domain.
205  /**
206  \em a_isPeriodic is a SpaceDim array of bools; if true, the
207  physical boundary condition is periodic in the coordinate direction.
208  False means a non-periodic BC.
209  */
210  void define(const Box& a_domBox, const bool* a_isPeriodic);
211 
212  /// Construct a ProblemDomain.
213  /**
214  It is an error if small is greater than big.
215  Defaults to non-periodic domain.
216  */
217  void define (const IntVect& small,
218  const IntVect& big);
219 
220  /// Construct a ProblemDomain.
221  /**
222  It is an error if small is greater than big.
223  \em a_isPeriodic is a SpaceDim array of bools; if true, the
224  physical boundary condition is periodic in the coordinate direction.
225  False means a non-periodic BC.
226  */
227  void define (const IntVect& small,
228  const IntVect& big,
229  const bool* a_isPeriodic);
230 
231  /// Construct ProblemDomain with specified lengths.
232  /**
233  It is an error if the lengths are negative.
234  Defaults to non-periodic domain.
235  */
236  void define (const IntVect& small,
237  const int* vec_len);
238 
239  /// Construct ProblemDomain with specified lengths.
240  /**
241  It is an error if the lengths are negative.
242  \em a_isPeriodic is a SpaceDim array of bools; if true, the
243  physical boundary condition is periodic in the coordinate direction.
244  False means a non-periodic BC.
245  */
246  void define (const IntVect& small,
247  const int* vec_len,
248  const bool* a_isPeriodic);
249 
250  /// The copy constructor.
251  /**
252  */
253  void define (const ProblemDomain& a_src);
254 
255  // Accessors
256 
257  /// Returns the logical computational domain
258  /**
259  */
260  const Box& domainBox() const;
261 
262  /// Returns true if BC is periodic in direction a_dir
263  /**
264  */
265  bool isPeriodic(int a_dir) const;
266 
267  /// Returns true is BC is periodic in _any_ direction
268  /**
269  */
270  bool isPeriodic() const;
271 
272  /// Returns the shiftIterator for this ProblemDomain
273  /**
274  The shiftIterator is defined based on the periodicity of this
275  ProblemDomain.
276  */
278 
279  /// Returns true if this ProblemDomain is empty or undefined.
280  /**
281  */
282  bool isEmpty () const;
283 
284  /// Return the size of the domainBox in the specified coordinate direction.
285  /**
286  */
287  int size (const int& a_idir) const;
288 
289  /// Return the size of the domainBox.
290  /**
291  */
292  IntVect size () const;
293 
294  /// Returns true if argument is contained within this ProblemDomain.
295  /**
296  An empty ProblemDomain does not contain and is not contained by
297  any ProblemDomain, including itself. Note also that in a periodic
298  Direction, any index is valid, since the domain is infinite in that
299  direction.
300  */
301  bool contains (const IntVect& p) const;
302 
303  /// Returns the periodic image of this IntVect inside of the ProblemDomain.
304  /**
305  Return true if the domain contains this IntVect, returning the image
306  in 'p', otherwise returns false
307  */
308  bool image (IntVect& p) const;
309 
310  /// Returns true if argument is contained within this ProblemDomain.
311  /**
312  An empty ProblemDomain does not contain any Box. An entirely periodic
313  domain contains any Box.
314  */
315  bool contains (const Box& b) const;
316 
317  /// Equivalent to \em contains() function.
318  bool contains_box(const Box& b) const
319  {
320  return contains(b);
321  }
322 
323  /// Returns true if this ProblemDomain and the argument intersect.
324  /**
325  It is an error if a_box is not cell-centered.
326  An empty ProblemDomain does not intersect any Box.
327  This will do nothing in periodic directions (since a periodic domain
328  is infinite in the periodic direction. If periodic in all dimensions,
329  this will always return true.
330  */
331  bool intersects (const Box& a_box) const;
332 
333  /// Returns true if this ProblemDomain and the argument intersect.
334  /**
335  It is an error if a_box is not cell-centered.
336  This routine does not perform the check to see if *this or b are
337  empty Boxes. It is the callers responsibility to ensure that
338  this never happens. If you are unsure, the use the .intersects(..)
339  routine. In periodic directions, will always return true.
340 
341  */
342  bool intersectsNotEmpty (const Box& a_box) const;
343 
344  /// Returns true if this argument is adjacent to a periodic boundary
345  /**
346  */
347  bool periodicAdjacent(const Box& a_box) const;
348 
349  ///
350  /**
351  */
352  void insertImages(std::list<Box>& a_list, const Box& a_box) const;
353 
354 
355  /// Returns true if box1 and box2 or any of their periodic images intersect
356  /**
357  (useful for checking disjointness)
358  */
359  bool intersects(const Box& box1, const Box& box2) const;
360 
361  /// Returns true if the two domain are equivalent
362  /**
363  */
364  bool operator== (const ProblemDomain& a_otherDomain) const;
365 
366  /// Returns true if the two domain are not equivalent
367  /**
368  */
369  bool operator!= (const ProblemDomain& a_otherDomain) const;
370 
371  /// Modifies a_box to be the intersection of a_box and \em a_probdomain
372  /**
373  */
374  friend void operator &=(Box& a_box, const ProblemDomain& a_probdomain);
375 
376  /// Returns a Box which is the interesection of a_box and \em a_probdomain
377  /**
378  */
379  friend Box operator & (const Box& a_box, const ProblemDomain& a_probdomain);
380 
381  // Modification Functions
382 
383  /// The assignment operator.
384  /**
385  */
387 
388  /// Sets whether BC is periodic in direction a_dir (true == periodic)
389  /**
390  */
391  void setPeriodic(int a_dir, bool a_isPeriodic);
392 
393  /// Grows (or shrinks) the domain Box by i in all directions
394  /**
395  */
396  inline ProblemDomain& grow (int i);
397 
398  /// Returns a ProblemDomain with a domainBox grown by the given amount
399  /**
400  */
401  friend inline ProblemDomain grow(const ProblemDomain& pd,
402  int i);
403 
404  /// Grow this ProblemDomain
405  /**
406  Modifies this ProblemDomain by growing the domainBox in each
407  direction by the specified amount
408  */
409  inline ProblemDomain& grow(const IntVect& v);
410 
411  /// Returns a grown version of \em pd.
412  /**
413  Returns a ProblemDomain that is the argument ProblemDomain
414  with a DomainBox grown by the given amount.
415  */
416  friend inline ProblemDomain grow (const ProblemDomain& pd,
417  const IntVect& v);
418 
419  /// Grow this ProblemDomain
420  /**
421  Modifies this ProblemDomain by growing it on the low and high end
422  by n_cell cells in direction idir.
423  */
424  inline ProblemDomain& grow(int idir, int n_cell);
425 
426  /// Grow this ProblemDomain on the low side
427  /**
428  Modifies this ProblemDomain by growing it on the low end by n_cell
429  cells in direction idir.
430  */
431  inline ProblemDomain& growLo(int idir, int n_cell=1);
432 
433  /// Grow this ProblemDomain on the high side
434  /** Modifies this ProblemDomain by growing it on the high end by n_Cell
435  cells in direction idir
436  */
437  inline ProblemDomain& growHi(int idir, int n_cell=1);
438 
439  /// Returns a face-centered Box at the low side of \em a_pd
440  /**
441  Returns the edge-centered Box (in direction \em a_dir) defining
442  the low side of the argument ProblemDomain. The output Box will
443  have the given length in the given direction. Directions are
444  zero-based. It is an error if not 0 <= dir < SpaceDim. The neighbor
445  of an Empty ProblemDomain is an Empty Box of the appropriate type.
446  If \em a_dir is a periodic direction, will return an empty Box.
447  */
448  friend Box bdryLo (const ProblemDomain& a_pd,
449  int a_dir,
450  int a_len=1);
451 
452  /// Returns a face-centered Box at the high side of \em a_pd
453  /**
454  Returns the edge-centered Box (in direction \em a_dir) defining the
455  high side of the argument ProblemDomain. The return Box will have the
456  given length in the given direction. Directions are zero-based.
457  It is an error if not 0 <= dir < SpaceDim. The neighbor of an
458  Empty ProblemDomain is an Empty Box of the appropriate type.
459  If \em a_dir is a periodic direction, will return an empty Box.
460  */
461  friend Box bdryHi (const ProblemDomain& a_pd,
462  int a_dir,
463  int a_len=1);
464 
465  /// Returns the cell-centered Box adjacent to the low side of \em a_pd.
466  /**
467  Returns the cell centered Box of the given length adjacent to the
468  argument ProblemDomain on the low end along the given coordinate direction.
469  The return Box is identical to the argument ProblemDomain in the other
470  directions. The return ProblemDomain and the argument ProblemDomain
471  have an empty intersection.
472 
473  \b NOTES:
474  - len >= 1.
475  - Box retval = adjCellLo(b,dir,len) is equivalent to the following set of
476  operations:
477  \code
478  Box retval(b);
479 
480  retval.convert(dir,ProblemDomain::CELL);
481 
482  retval.setrange(dir,retval.smallEnd(dir)-len,len);
483  \endcode
484 
485  Directions are zero-based. It is an error if not 0 <= dir <
486  SpaceDim. The neighbor of an Empty ProblemDomain is an Empty Box of the
487  appropriate type. If \em a_dir is a periodic direction, will return
488  an empty Box as well.
489 
490  */
491  friend Box adjCellLo (const ProblemDomain& a_pd,
492  int a_dir,
493  int a_len=1);
494 
495  /// Returns the cell-centered Box adjacent to the low side of \em a_pd.
496  /**
497  Returns the cell centered Box of the given length adjacent to the
498  argument ProblemDomain on the high end along the given coordinate
499  direction. The return Box is identical to the argument ProblemDomain in
500  the other directions. The return Box and the argument ProblemDomain have
501  an empty intersection.
502 
503  \b NOTES:
504  - len >= 1.
505 
506  - Box retval = adjCellHi(b,dir,len) is equivalent to the
507  following set of operations:
508  \code
509  Box retval(b);
510 
511  retval.convert(dir,ProblemDomain::CELL);
512 
513  retval.setrange(dir,retval.bigEnd(dir)+1,len);
514  \endcode
515 
516  Directions are zero-based. It is an error if not 0 <= dir <
517  SpaceDim. The neighbor of an Empty ProblemDomain is an Empty Box of the
518  appropriate type. If \em a_dir is a periodic direction, will return
519  an empty Box.
520 
521  */
522  friend Box adjCellHi (const ProblemDomain& a_pd,
523  int a_dir,
524  int a_len=1);
525 
526  // Intersection functions
527 
528  /// Returns the Box intersection of this ProblemDomain and \em a_b
529  /**
530  The intersection of the Empty ProblemDomain and any Box is the Empty
531  Box. This operator does nothing in periodic directions (since
532  a periodic domain is an infinite domain).
533 
534  */
535  Box operator& (const Box& a_b) const;
536 
537  // refinement
538 
539  /// Refine this problem domain.
540  /**
541  Modifies this ProblemDomain by refining it by given (positive) refinement
542  ratio. The Empty ProblemDomain is not modified by this function.
543  */
544  ProblemDomain& refine (int a_refinement_ratio);
545 
546  /// Return a ProblemDomain which is a refinement of \em a_probdomain
547  /**
548  Returns a ProblemDomain that is the argument ProblemDomain refined by
549  given (positive) refinement ratio. The Empty ProblemDomain is not
550  modified by this function.
551  */
552  friend ProblemDomain refine (const ProblemDomain& a_probdomain,
553  int a_refinement_ratio);
554 
555  /// Refine this ProblemDomain
556  /**
557  Modifies this ProblemDomain by refining it by given (positive) refinement
558  ratio. The Empty ProblemDomain is not modified by this function.
559  */
560  ProblemDomain& refine (const IntVect& a_refinement_ratio);
561 
562  /// Refinement function
563  /**
564  Returns a ProblemDomain that is the argument ProblemDomain refined
565  by given (positive) refinement ratio. The Empty ProblemDomain is
566  not modified by this function.
567 
568  */
569  friend ProblemDomain refine (const ProblemDomain& a_probdomain,
570  const IntVect& a_refinement_ratio);
571 
572  // coarsening
573 
574  /// Coarsen this ProblemDomain
575  /**
576  Modifies this ProblemDomain by coarsening it by given (positive)
577  refinement ratio. The Empty ProblemDomain is not modified by
578  this function.
579  */
580 
581  ProblemDomain& coarsen (int a_refinement_ratio);
582 
583  /// Coarsening function
584  /**
585  Returns a ProblemDomain that is the argument ProblemDomain coarsened
586  by given (positive) refinement ratio. The Empty ProblemDomain is not
587  modified by this function.
588  */
589  friend ProblemDomain coarsen (const ProblemDomain& a_probdomain,
590  int a_refinement_ratio);
591 
592  /// Coarsen this ProblemDomain
593  /**
594  Modifies this ProblemDomain by coarsening by given (positive) refinement
595  ratio. The Empty ProblemDomain is not modified by this function.
596  */
597  ProblemDomain& coarsen (const IntVect& refinement_ratio);
598 
599  /// Coarsening function
600  /**
601  Returns a ProblemDomain that is the argument ProblemDomain coarsened
602  by given (positive) refinement ratio. The Empty ProblemDomain is
603  not modified by this function.
604 
605  */
606  friend ProblemDomain coarsen (const ProblemDomain& a_probdomain,
607  const IntVect& a_refinement_ratio);
608 
609 
610  ///
611  void shift(const IntVect& a_shift)
612  {
613  m_domainBox.shift(a_shift);
614  }
615  // I/O Functions
616 
617  /// Write an ASCII representation to the ostream.
618  /**
619  */
620  friend std::ostream& operator<< (std::ostream& os,
621  const ProblemDomain& bx);
622 
623  /// Read from istream.
624  /**
625  */
626  friend std::istream& operator>> (std::istream& is,
627  ProblemDomain& bx);
628 
629  void shiftIt(Box& a_box, int shiftIndex) const ;
630  void unshiftIt(Box& a_box, int shiftIndex) const ;
631  /// Gives more detail than printOn.
632  /**
633  Useful for exiting due to an error.
634  */
635  void dumpOn (std::ostream& strm) const;
636 
637 
638 protected:
639  friend class HDF5Handle;
640 
641  /**
642  Periodicity info
643  */
645 
646  /**
647  Domain index extents
648  */
650 
651  /**
652  Shift iterator for this ProblemDomain
653  */
655 };
656 
658 {
659 public:
660  ImageIterator(const ProblemDomain& a_domain)
661  {
662  define(a_domain);
663  }
664 
665  void define(const ProblemDomain& a_domain);
666 
667  void begin(const Box& a_box)
668  {
669  m_counter=-1;
670  m_box = a_box;
671  this->operator++();
672  }
673 
674  void operator++();
675 
676  bool ok()
677  {
678  return m_shifter[m_counter] != IntVect::Zero;
679  }
680 
681  const Box& box() const
682  {
683  return m_current;
684  }
685 
686  const ProblemDomain& domain() const
687  {
688  return m_domain;
689  }
690 
691  void checkDefine(const ProblemDomain& a_domain)
692  {
693  if (!(m_domain == a_domain)) define(a_domain);
694  }
695 
696 protected:
698  Box m_quadrant[D_TERM6(3,*3,*3,*3,*3,*3)];
699  IntVect m_shifter[D_TERM6(3,*3,*3,*3,*3,*3)];
703 };
704 
705 //
706 // Inlines.
707 //
708 #ifndef WRAPPER
709 
710 inline
712  : m_index(100), m_shift_vectors()
713 {
714 }
715 
716 inline
718 {
719  m_index = a_src.m_index;
721 }
722 
723 inline
726 {
727  m_index = a_src.m_index;
729  return *this;
730 }
731 
732 inline
733 IntVect
735 {
736  CH_assert(ok());
737  return m_shift_vectors[m_index];
738 }
739 
740 inline
741 void
743 {
744  m_index++;
745 }
746 
747 inline
748 bool
750 {
751  return (m_index < m_shift_vectors.size());
752 }
753 
754 inline
755 void
757 {
758  m_index = 0;
759 }
760 
761 inline
762 void
764 {
765  m_index = 0;
766 }
767 
768 inline
769 void
771 {
773 }
774 
775 inline
777 {
778  // default is a non-periodic domain
779  for (int dir=0; dir<SpaceDim; dir++)
780  {
781  m_isPeriodic[dir] = false;
782  }
783 
784 }
785 
786 inline
788  : m_domainBox(b.m_domainBox), m_shiftIt(b.m_shiftIt)
789 {
790  for (int dir=0; dir<SpaceDim; dir++)
791  {
792  m_isPeriodic[dir] = b.m_isPeriodic[dir];
793  }
794 }
795 
796 inline void
798 
799 {
801  m_shiftIt = (b.m_shiftIt);
802  for (int dir=0; dir<SpaceDim; dir++)
803  {
804  m_isPeriodic[dir] = b.m_isPeriodic[dir];
805  }
806 }
807 
808 inline
809 bool
810 ProblemDomain::operator== (const ProblemDomain& a_otherDomain) const
811 {
812  bool result = true;
813 
814  if (m_domainBox != a_otherDomain.m_domainBox)
815  {
816  result = false;
817  }
818  else
819  {
820  for (int dir=0; dir<SpaceDim; dir++)
821  {
822  if (m_isPeriodic[dir] != a_otherDomain.m_isPeriodic[dir])
823  {
824  result = false;
825  break;
826  }
827  }
828  }
829 
830  return result;
831 }
832 
833 inline
834 bool
835 ProblemDomain::operator!= (const ProblemDomain& a_otherDomain) const
836 {
837  return !(*this == a_otherDomain);
838 }
839 
840 inline
843 {
845  for (int dir=0; dir<SpaceDim; dir++)
846  {
847  m_isPeriodic[dir] = b.m_isPeriodic[dir];
848  }
849  m_shiftIt = b.m_shiftIt;
850  return *this;
851 }
852 
853 inline void
854 ProblemDomain::shiftIt(Box& a_box, int a_shiftIndex) const
855 {
856  a_box.shift(m_shiftIt[a_shiftIndex]*m_domainBox.size());
857 }
858 
859 inline void
860 ProblemDomain::unshiftIt(Box& a_box, int a_shiftIndex) const
861 {
862  a_box.shift(- m_shiftIt[a_shiftIndex]*m_domainBox.size());
863 }
864 
865 
866 inline
867 const Box&
869 {
870  return m_domainBox;
871 }
872 
873 inline
874 bool
876 {
877  return m_isPeriodic[a_dir];
878 }
879 
880 inline
881 bool
883 {
884  return D_TERM6(m_isPeriodic[0], ||
885  m_isPeriodic[1], ||
886  m_isPeriodic[2], ||
887  m_isPeriodic[3], ||
888  m_isPeriodic[4], ||
889  m_isPeriodic[5]);
890 }
891 
892 inline
895 {
896  m_domainBox.grow(i);
897  return *this;
898 }
899 
900 inline
902 grow(const ProblemDomain& pd, int i)
903 {
904  ProblemDomain newPd(pd);
905  newPd.grow(i);
906  return newPd;
907 }
908 
909 inline
912 {
913  m_domainBox.grow(v);
914  return *this;
915 }
916 
917 inline
919 grow(const ProblemDomain& pd, const IntVect& v)
920 {
921  ProblemDomain newPd(pd);
922  newPd.grow(v);
923  return newPd;
924 }
925 
926 inline
928 ProblemDomain::grow(int idir, int n_cell)
929 {
930  m_domainBox.grow(idir, n_cell);
931  return *this;
932 }
933 
934 inline
936 ProblemDomain::growLo(int idir, int n_cell)
937 {
938  m_domainBox.growLo(idir, n_cell);
939  return *this;
940 }
941 
942 inline
944 ProblemDomain::growHi(int idir, int n_cell)
945 {
946  m_domainBox.growHi(idir, n_cell);
947  return *this;
948 }
949 
950 inline
953 {
954  return m_shiftIt;
955 }
956 
957 inline
958 bool
960 {
961  return (m_domainBox.isEmpty());
962 }
963 
964 inline
965 int
966 ProblemDomain::size (const int& a_idir) const
967 {
968  return (m_domainBox.size(a_idir));
969 }
970 
971 inline
972 IntVect
974 {
975  return (m_domainBox.size());
976 }
977 
978 inline
979 bool
981 {
982 
983  // boy is this ugly!
984  return ( !isEmpty()
985  && (D_TERM6((m_isPeriodic[0]
986  || ( p[0] >= m_domainBox.smallEnd(0)
987  && p[0] <= m_domainBox.bigEnd (0))),
988  && (m_isPeriodic[1]
989  || ( p[1] >= m_domainBox.smallEnd(1)
990  && p[1] <= m_domainBox.bigEnd (1))),
991  && (m_isPeriodic[2]
992  || ( p[2] >= m_domainBox.smallEnd(2)
993  && p[2] <= m_domainBox.bigEnd (2))),
994  && (m_isPeriodic[3]
995  || ( p[3] >= m_domainBox.smallEnd(3)
996  && p[3] <= m_domainBox.bigEnd (3))),
997  && (m_isPeriodic[4]
998  || ( p[4] >= m_domainBox.smallEnd(4)
999  && p[4] <= m_domainBox.bigEnd (4))),
1000  && (m_isPeriodic[5]
1001  || ( p[5] >= m_domainBox.smallEnd(5)
1002  && p[5] <= m_domainBox.bigEnd (5))))));
1003 }
1004 
1005 inline
1006 bool
1008 {
1009  if (m_domainBox.contains(p)) return true;
1010  if (!contains(p)) return false;
1011 
1012  D_TERM6(
1013  if (m_isPeriodic[0])
1014  {
1015  if (p[0]<m_domainBox.smallEnd(0)) p[0]+= m_domainBox.size(0);
1016  else if (p[0]>m_domainBox.bigEnd(0)) p[0]-= m_domainBox.size(0);
1017  },
1018  if (m_isPeriodic[1])
1019  {
1020  if (p[1]<m_domainBox.smallEnd(1)) p[1]+= m_domainBox.size(1);
1021  else if (p[1]>m_domainBox.bigEnd(1)) p[1]-= m_domainBox.size(1);
1022  },
1023  if (m_isPeriodic[2])
1024  {
1025  if (p[2]<m_domainBox.smallEnd(2)) p[2]+= m_domainBox.size(2);
1026  else if (p[2]>m_domainBox.bigEnd(2)) p[2]-= m_domainBox.size(2);
1027  },
1028  if (m_isPeriodic[3])
1029  {
1030  if (p[3]<m_domainBox.smallEnd(3)) p[3]+= m_domainBox.size(3);
1031  else if (p[3]>m_domainBox.bigEnd(3)) p[3]-= m_domainBox.size(3);
1032  },
1033  if (m_isPeriodic[4])
1034  {
1035  if (p[4]<m_domainBox.smallEnd(4)) p[4]+= m_domainBox.size(4);
1036  else if (p[4]>m_domainBox.bigEnd(4)) p[4]-= m_domainBox.size(4);
1037  },
1038  if (m_isPeriodic[5])
1039  {
1040  if (p[5]<m_domainBox.smallEnd(5)) p[5]+= m_domainBox.size(5);
1041  else if (p[5]>m_domainBox.bigEnd(5)) p[5]-= m_domainBox.size(5);
1042  });
1043 
1044  return true;
1045 }
1046 
1047 inline
1048 bool
1050 {
1051  // boy is this ugly!
1052  if (b.type() == m_domainBox.type())
1053  {
1054  return ( !isEmpty()
1055  && (D_TERM6((m_isPeriodic[0]
1056  || ( b.smallEnd(0) >= m_domainBox.smallEnd(0)
1057  && b.bigEnd (0) <= m_domainBox.bigEnd (0))), &&
1058  (m_isPeriodic[1]
1059  || ( b.smallEnd(1) >= m_domainBox.smallEnd(1)
1060  && b.bigEnd (1) <= m_domainBox.bigEnd (1))), &&
1061  (m_isPeriodic[2]
1062  || ( b.smallEnd(2) >= m_domainBox.smallEnd(2)
1063  && b.bigEnd (2) <= m_domainBox.bigEnd (2))), &&
1064  (m_isPeriodic[3]
1065  || ( b.smallEnd(3) >= m_domainBox.smallEnd(3)
1066  && b.bigEnd (3) <= m_domainBox.bigEnd (3))), &&
1067  (m_isPeriodic[4]
1068  || ( b.smallEnd(4) >= m_domainBox.smallEnd(4)
1069  && b.bigEnd (4) <= m_domainBox.bigEnd (4))), &&
1070  (m_isPeriodic[5]
1071  || ( b.smallEnd(5) >= m_domainBox.smallEnd(5)
1072  && b.bigEnd (5) <= m_domainBox.bigEnd (5))))));
1073  }
1074  else
1075  {
1077 
1078  // Check b's centering and adjust intersectBox as needed
1079  for (int dir = 0; dir < SpaceDim; dir++)
1080  {
1081  if (b.type(dir) != domainBox.type(dir))
1082  {
1083  if (b.type(dir) == IndexType::NODE)
1084  {
1085  domainBox.surroundingNodes(dir);
1086  }
1087  else
1088  {
1089  domainBox.enclosedCells(dir);
1090  }
1091  }
1092  }
1093 
1094  return ( !isEmpty()
1095  && (D_TERM6((m_isPeriodic[0]
1096  || ( b.smallEnd(0) >= domainBox.smallEnd(0)
1097  && b.bigEnd (0) <= domainBox.bigEnd (0))), &&
1098  (m_isPeriodic[1]
1099  || ( b.smallEnd(1) >= domainBox.smallEnd(1)
1100  && b.bigEnd (1) <= domainBox.bigEnd (1))), &&
1101  (m_isPeriodic[2]
1102  || ( b.smallEnd(2) >= domainBox.smallEnd(2)
1103  && b.bigEnd (2) <= domainBox.bigEnd (2))), &&
1104  (m_isPeriodic[3]
1105  || ( b.smallEnd(3) >= domainBox.smallEnd(3)
1106  && b.bigEnd (3) <= domainBox.bigEnd (3))), &&
1107  (m_isPeriodic[4]
1108  || ( b.smallEnd(4) >= domainBox.smallEnd(4)
1109  && b.bigEnd (4) <= domainBox.bigEnd (4))), &&
1110  (m_isPeriodic[5]
1111  || ( b.smallEnd(5) >= domainBox.smallEnd(5)
1112  && b.bigEnd (5) <= domainBox.bigEnd (5))))));
1113  }
1114 }
1115 
1116 #endif /*WRAPPER*/
1117 
1118 #include "NamespaceFooter.H"
1119 #endif
Box m_domainBox
Definition: ProblemDomain.H:649
Definition: Box.H:52
Box m_quadrant[D_TERM6(3,*3,*3,*3,*3,*3)]
Definition: ProblemDomain.H:698
bool image(IntVect &p) const
Returns the periodic image of this IntVect inside of the ProblemDomain.
Definition: ProblemDomain.H:1007
friend Box bdryHi(const ProblemDomain &a_pd, int a_dir, int a_len=1)
Returns a face-centered Box at the high side of a_pd.
#define D_TERM6(a, b, c, d, e, f)
Definition: CHArray.H:40
#define CH_assert(cond)
Definition: CHArray.H:37
A class to facilitate interaction with physical boundary conditions.
Definition: ProblemDomain.H:130
ProblemDomain & growLo(int idir, int n_cell=1)
Grow this ProblemDomain on the low side.
Definition: ProblemDomain.H:936
int m_index
Definition: ProblemDomain.H:102
void shift(const IntVect &a_shift)
Definition: ProblemDomain.H:611
bool intersectsNotEmpty(const Box &a_box) const
Returns true if this ProblemDomain and the argument intersect.
bool contains_box(const Box &b) const
Equivalent to contains() function.
Definition: ProblemDomain.H:318
friend std::istream & operator>>(std::istream &is, ProblemDomain &bx)
Read from istream.
friend void operator&=(Box &a_box, const ProblemDomain &a_probdomain)
Modifies a_box to be the intersection of a_box and a_probdomain.
friend Box adjCellLo(const ProblemDomain &a_pd, int a_dir, int a_len=1)
Returns the cell-centered Box adjacent to the low side of a_pd.
friend Box operator&(const Box &a_box, const ProblemDomain &a_probdomain)
Returns a Box which is the interesection of a_box and a_probdomain.
IntVect size() const
size functions
Definition: Box.H:1814
bool isEmpty() const
Returns true if this ProblemDomain is empty or undefined.
Definition: ProblemDomain.H:959
ProblemDomain & operator=(const ProblemDomain &b)
The assignment operator.
Definition: ProblemDomain.H:842
Box & surroundingNodes()
ShiftIterator shiftIterator() const
Returns the shiftIterator for this ProblemDomain.
Definition: ProblemDomain.H:952
bool m_isPeriodic[SpaceDim]
Definition: ProblemDomain.H:644
ProblemDomain m_domain
Definition: ProblemDomain.H:697
const int SpaceDim
Definition: SPACE.H:39
const Box & box() const
Definition: ProblemDomain.H:681
bool ok()
Definition: ProblemDomain.H:676
void incr()
Equivalent to the ++ operator.
Definition: ProblemDomain.H:70
int index() const
which periodic image is this
Definition: ProblemDomain.H:79
void setPeriodic(int a_dir, bool a_isPeriodic)
Sets whether BC is periodic in direction a_dir (true == periodic)
Box & enclosedCells()
ShiftIterator()
Basic constructor.
Definition: ProblemDomain.H:711
friend Box adjCellHi(const ProblemDomain &a_pd, int a_dir, int a_len=1)
Returns the cell-centered Box adjacent to the low side of a_pd.
void begin(const Box &a_box)
Definition: ProblemDomain.H:667
bool isPeriodic() const
Returns true is BC is periodic in any direction.
Definition: ProblemDomain.H:882
IntVect i() const
Equivalent to the () operator.
Definition: ProblemDomain.H:62
Class to manage box-shifting used to enforce periodic BC's.
Definition: ProblemDomain.H:33
void computeShifts(const bool *a_isPeriodic)
Recompute shift vectors based on periodic directions.
void define(const ProblemDomain &a_domain)
bool intersects(const Box &a_box) const
Returns true if this ProblemDomain and the argument intersect.
void begin()
Equivalent to reset()
Definition: ProblemDomain.H:763
bool ok() const
Is the iterator still within its range of shift vectors?
Definition: ProblemDomain.H:749
friend std::ostream & operator<<(std::ostream &os, const ProblemDomain &bx)
Write an ASCII representation to the ostream.
const IntVect & bigEnd() const
Definition: Box.H:1779
ShiftIterator m_shiftIt
Definition: ProblemDomain.H:654
const Box & domainBox() const
Returns the logical computational domain.
Definition: ProblemDomain.H:868
Vector< IntVect > m_shift_vectors
Definition: ProblemDomain.H:104
void insertImages(std::list< Box > &a_list, const Box &a_box) const
int m_counter
Definition: ProblemDomain.H:702
bool isEmpty() const
{ Comparison Functions}
Definition: Box.H:1863
IntVect type() const
Definition: Box.H:1800
static const IntVect Zero
Definition: IntVect.H:627
void define(const Box &a_domBox)
Construct ProblemDomain with a_domBox as computational domain.
void checkDefine(const ProblemDomain &a_domain)
Definition: ProblemDomain.H:691
ProblemDomain & refine(int a_refinement_ratio)
Refine this problem domain.
bool operator==(const ProblemDomain &a_otherDomain) const
Returns true if the two domain are equivalent.
Definition: ProblemDomain.H:810
A Rectangular Domain on an Integer Lattice.
Definition: Box.H:465
ProblemDomain & growHi(int idir, int n_cell=1)
Grow this ProblemDomain on the high side.
Definition: ProblemDomain.H:944
void unshiftIt(Box &a_box, int shiftIndex) const
Definition: ProblemDomain.H:860
ImageIterator(const ProblemDomain &a_domain)
Definition: ProblemDomain.H:660
void operator++()
Increment to the next shift unit vector.
Definition: ProblemDomain.H:742
bool contains(const IntVect &p) const
Definition: Box.H:1888
Handle to a particular group in an HDF file.
Definition: CH_HDF5.H:267
IntVect m_shifter[D_TERM6(3,*3,*3,*3,*3,*3)]
Definition: ProblemDomain.H:699
friend Box bdryLo(const ProblemDomain &a_pd, int a_dir, int a_len=1)
Returns a face-centered Box at the low side of a_pd.
An integer Vector in SpaceDim-dimensional space.
Definition: CHArray.H:42
ProblemDomain & coarsen(int a_refinement_ratio)
Coarsen this ProblemDomain.
bool periodicAdjacent(const Box &a_box) const
Returns true if this argument is adjacent to a periodic boundary.
IntVect operator()() const
Returns the current shift unit vector.
Definition: ProblemDomain.H:734
ProblemDomain()
The default constructor. The constructed domain box is empty.
Definition: ProblemDomain.H:776
bool contains(const IntVect &p) const
Returns true if argument is contained within this ProblemDomain.
Definition: ProblemDomain.H:980
ProblemDomain grow(const ProblemDomain &pd, int i)
Definition: ProblemDomain.H:902
size_t size() const
Definition: Vector.H:177
ShiftIterator & operator=(const ShiftIterator &a_src)
Assignment operator.
Definition: ProblemDomain.H:725
void reset()
Reset to first shift unit vector.
Definition: ProblemDomain.H:756
Box & shift(int dir, int nzones)
shift functions
Definition: Box.H:2088
Box m_box
Definition: ProblemDomain.H:700
Box & grow(int i)
grow functions
Definition: Box.H:2268
Box & growLo(int idir, int n_cell=1)
Definition: Box.H:2354
void operator++()
ProblemDomain & grow(int i)
Grows (or shrinks) the domain Box by i in all directions.
Definition: ProblemDomain.H:894
Definition: ProblemDomain.H:657
IntVect size() const
Return the size of the domainBox.
Definition: ProblemDomain.H:973
Box m_current
Definition: ProblemDomain.H:701
~ShiftIterator()
Destructor.
const IntVect & smallEnd() const
{ Accessors}
Definition: Box.H:1765
bool operator!=(const ProblemDomain &a_otherDomain) const
Returns true if the two domain are not equivalent.
Definition: ProblemDomain.H:835
void end()
Skip the iterator to the end.
Definition: ProblemDomain.H:770
Box & growHi(int idir, int n_cell=1)
Definition: Box.H:2385
void shiftIt(Box &a_box, int shiftIndex) const
Definition: ProblemDomain.H:854
const ProblemDomain & domain() const
Definition: ProblemDomain.H:686
void dumpOn(std::ostream &strm) const
Gives more detail than printOn.
const IntVect & operator[](int index) const
Definition: ProblemDomain.H:84