00001 #ifdef CH_LANG_CC
00002
00003
00004
00005
00006
00007
00008
00009 #endif
00010
00011 #ifndef _PROBLEMDOMAIN_H_
00012 #define _PROBLEMDOMAIN_H_
00013
00014 #ifndef WRAPPER
00015 #include <iostream>
00016
00017 #include "Vector.H"
00018 #include "IntVect.H"
00019 #include "Box.H"
00020 #include "Misc.H"
00021 #endif
00022
00023 #include "SPACE.H"
00024 #include "NamespaceHeader.H"
00025
00027
00033 class ShiftIterator
00034 {
00035 public:
00037 ShiftIterator();
00038
00040
00044 ShiftIterator(const bool* a_isPeriodic);
00045
00047 ShiftIterator(const ShiftIterator& a_shiftIt);
00048
00050 ~ShiftIterator();
00051
00053 ShiftIterator& operator=(const ShiftIterator& a_src);
00054
00056 void computeShifts(const bool* a_isPeriodic);
00057
00059 inline IntVect operator()() const;
00060
00062 IntVect i() const {return this->operator()();}
00064 inline void operator++();
00065
00067 void incr(){++(*this);}
00068
00070 inline bool ok() const;
00071
00073 int index() const { return m_index;}
00074
00075 const IntVect& operator[](int index) const {return m_shift_vectors[index];}
00077 inline void reset();
00078
00080 inline void begin();
00081
00083
00086 void end();
00087
00088 private:
00089
00090 int m_index;
00091
00092 Vector<IntVect> m_shift_vectors;
00093
00094 };
00095
00097
00118 class ProblemDomain
00119 {
00120 public:
00121
00122
00123
00125
00127 ProblemDomain ();
00128
00130
00133 ProblemDomain(const Box& a_domBox);
00134
00136
00141 ProblemDomain(const Box& a_domBox, const bool* a_isPeriodic);
00142
00144
00148 ProblemDomain (const IntVect& small,
00149 const IntVect& big);
00150
00152
00158 ProblemDomain (const IntVect& small,
00159 const IntVect& big,
00160 const bool* a_isPeriodic);
00161
00163
00167 ProblemDomain (const IntVect& small,
00168 const int* vec_len);
00169
00171
00177 ProblemDomain (const IntVect& small,
00178 const int* vec_len,
00179 const bool* a_isPeriodic);
00180
00182
00184 ProblemDomain (const ProblemDomain& a_src);
00185
00187
00190 void define(const Box& a_domBox);
00191
00193
00198 void define(const Box& a_domBox, const bool* a_isPeriodic);
00199
00201
00205 void define (const IntVect& small,
00206 const IntVect& big);
00207
00209
00215 void define (const IntVect& small,
00216 const IntVect& big,
00217 const bool* a_isPeriodic);
00218
00220
00224 void define (const IntVect& small,
00225 const int* vec_len);
00226
00228
00234 void define (const IntVect& small,
00235 const int* vec_len,
00236 const bool* a_isPeriodic);
00237
00239
00241 void define (const ProblemDomain& a_src);
00242
00243
00244
00246
00248 const Box& domainBox() const;
00249
00251
00253 bool isPeriodic(int a_dir) const;
00254
00256
00258 bool isPeriodic() const;
00259
00261
00265 ShiftIterator shiftIterator() const;
00266
00268
00270 bool isEmpty () const;
00271
00273
00275 int size (const int& a_idir) const;
00276
00278
00280 IntVect size () const;
00281
00283
00289 bool contains (const IntVect& p) const;
00290
00292
00296 bool image (IntVect& p) const;
00297
00299
00303 bool contains (const Box& b) const;
00304
00306 bool contains_box(const Box& b) const {return contains(b);}
00307
00309
00316 bool intersects (const Box& a_box) const;
00317
00319
00327 bool intersectsNotEmpty (const Box& a_box) const;
00328
00330
00332 bool periodicAdjacent(const Box& a_box) const;
00333
00335
00337 void insertImages(std::list<Box>& a_list, const Box& a_box) const;
00338
00339
00341
00344 bool intersects(const Box& box1, const Box& box2) const;
00345
00347
00349 bool operator== (const ProblemDomain& a_otherDomain) const;
00350
00352
00354 bool operator!= (const ProblemDomain& a_otherDomain) const;
00355
00357
00359 friend void operator &=(Box& a_box, const ProblemDomain& a_probdomain);
00360
00362
00364 friend Box operator & (const Box& a_box, const ProblemDomain& a_probdomain);
00365
00366
00367
00369
00371 ProblemDomain& operator= (const ProblemDomain& b);
00372
00374
00376 void setPeriodic(int a_dir, bool a_isPeriodic);
00377
00379
00381 inline ProblemDomain& grow (int i);
00382
00384
00386 friend inline ProblemDomain grow(const ProblemDomain& pd,
00387 int i);
00388
00390
00394 inline ProblemDomain& grow(const IntVect& v);
00395
00397
00401 friend inline ProblemDomain grow (const ProblemDomain& pd,
00402 const IntVect& v);
00403
00405
00409 inline ProblemDomain& grow(int idir, int n_cell);
00410
00412
00416 inline ProblemDomain& growLo(int idir, int n_cell=1);
00417
00419
00422 inline ProblemDomain& growHi(int idir, int n_cell=1);
00423
00425
00433 friend Box bdryLo (const ProblemDomain& a_pd,
00434 int a_dir,
00435 int a_len=1);
00436
00438
00446 friend Box bdryHi (const ProblemDomain& a_pd,
00447 int a_dir,
00448 int a_len=1);
00449
00451
00476 friend Box adjCellLo (const ProblemDomain& a_pd,
00477 int a_dir,
00478 int a_len=1);
00479
00481
00507 friend Box adjCellHi (const ProblemDomain& a_pd,
00508 int a_dir,
00509 int a_len=1);
00510
00511
00512
00514
00520 Box operator& (const Box& a_b) const;
00521
00522
00523
00525
00529 ProblemDomain& refine (int a_refinement_ratio);
00530
00532
00537 friend ProblemDomain refine (const ProblemDomain& a_probdomain,
00538 int a_refinement_ratio);
00539
00541
00545 ProblemDomain& refine (const IntVect& a_refinement_ratio);
00546
00548
00554 friend ProblemDomain refine (const ProblemDomain& a_probdomain,
00555 const IntVect& a_refinement_ratio);
00556
00557
00558
00560
00566 ProblemDomain& coarsen (int a_refinement_ratio);
00567
00569
00574 friend ProblemDomain coarsen (const ProblemDomain& a_probdomain,
00575 int a_refinement_ratio);
00576
00578
00582 ProblemDomain& coarsen (const IntVect& refinement_ratio);
00583
00585
00591 friend ProblemDomain coarsen (const ProblemDomain& a_probdomain,
00592 const IntVect& a_refinement_ratio);
00593
00594
00596 void shift(const IntVect& a_shift)
00597 {
00598 m_domainBox.shift(a_shift);
00599 }
00600
00601
00603
00605 friend std::ostream& operator<< (std::ostream& os,
00606 const ProblemDomain& bx);
00607
00609
00611 friend std::istream& operator>> (std::istream& is,
00612 ProblemDomain& bx);
00613
00614 void shiftIt(Box& a_box, int shiftIndex) const ;
00615 void unshiftIt(Box& a_box, int shiftIndex) const ;
00617
00620 void dumpOn (std::ostream& strm) const;
00621
00622
00623 protected:
00624 friend class HDF5Handle;
00625
00629 bool m_isPeriodic[SpaceDim];
00630
00634 Box m_domainBox;
00635
00639 ShiftIterator m_shiftIt;
00640 };
00641
00642 class ImageIterator
00643 {
00644 public:
00645 ImageIterator(const ProblemDomain& a_domain){ define(a_domain);}
00646 void define(const ProblemDomain& a_domain);
00647 void begin(const Box& a_box)
00648 {
00649 m_counter=-1;
00650 m_box = a_box;
00651 this->operator++();
00652 }
00653 void operator++();
00654 bool ok() { return m_shifter[m_counter] != IntVect::Zero;}
00655 const Box& box() const {return m_current;}
00656 const ProblemDomain& domain() const {return m_domain;}
00657 void checkDefine(const ProblemDomain& a_domain)
00658 {
00659 if(!(m_domain == a_domain)) define(a_domain);
00660 }
00661 protected:
00662 ProblemDomain m_domain;
00663 Box m_quadrant[D_TERM6(3,*3,*3,*3,*3,*3)];
00664 IntVect m_shifter[D_TERM6(3,*3,*3,*3,*3,*3)];
00665 Box m_box;
00666 Box m_current;
00667 int m_counter;
00668 };
00669
00670
00671
00672
00673
00674
00675
00676
00677 #ifndef WRAPPER
00678
00679 inline
00680 ShiftIterator::ShiftIterator()
00681 : m_index(100), m_shift_vectors()
00682 {
00683 }
00684
00685 inline
00686 ShiftIterator::ShiftIterator(const ShiftIterator& a_src)
00687 {
00688 m_index = a_src.m_index;
00689 m_shift_vectors = a_src.m_shift_vectors;
00690 }
00691
00692 inline
00693 ShiftIterator&
00694 ShiftIterator::operator=(const ShiftIterator& a_src)
00695 {
00696 m_index = a_src.m_index;
00697 m_shift_vectors = a_src.m_shift_vectors;
00698 return *this;
00699 }
00700
00701 inline
00702 IntVect
00703 ShiftIterator::operator()() const
00704 {
00705 CH_assert(ok());
00706 return m_shift_vectors[m_index];
00707 }
00708
00709 inline
00710 void
00711 ShiftIterator::operator++()
00712 {
00713 m_index++;
00714 }
00715
00716 inline
00717 bool
00718 ShiftIterator::ok() const
00719 {
00720 return (m_index < m_shift_vectors.size());
00721 }
00722
00723 inline
00724 void
00725 ShiftIterator::reset()
00726 {
00727 m_index = 0;
00728 }
00729
00730 inline
00731 void
00732 ShiftIterator::begin()
00733 {
00734 m_index = 0;
00735 }
00736
00737 inline
00738 void
00739 ShiftIterator::end()
00740 {
00741 m_index = m_shift_vectors.size();
00742 }
00743
00744 inline
00745 ProblemDomain::ProblemDomain ()
00746 {
00747
00748 for (int dir=0; dir<SpaceDim; dir++)
00749 {
00750 m_isPeriodic[dir] = false;
00751 }
00752
00753 }
00754
00755 inline
00756 ProblemDomain::ProblemDomain (const ProblemDomain& b)
00757 : m_domainBox(b.m_domainBox), m_shiftIt(b.m_shiftIt)
00758 {
00759 for (int dir=0; dir<SpaceDim; dir++)
00760 {
00761 m_isPeriodic[dir] = b.m_isPeriodic[dir];
00762 }
00763 }
00764
00765 inline
00766 bool
00767 ProblemDomain::operator== (const ProblemDomain& a_otherDomain) const
00768 {
00769 bool result = true;
00770
00771 if (m_domainBox != a_otherDomain.m_domainBox)
00772 {
00773 result = false;
00774 }
00775 else
00776 {
00777 for (int dir=0; dir<SpaceDim; dir++)
00778 {
00779 if (m_isPeriodic[dir] != a_otherDomain.m_isPeriodic[dir])
00780 {
00781 result = false;
00782 break;
00783 }
00784 }
00785 }
00786
00787 return result;
00788 }
00789
00790 inline
00791 bool
00792 ProblemDomain::operator!= (const ProblemDomain& a_otherDomain) const
00793 {
00794 return !(*this == a_otherDomain);
00795 }
00796
00797 inline
00798 ProblemDomain&
00799 ProblemDomain::operator= (const ProblemDomain& b)
00800 {
00801 m_domainBox = b.m_domainBox;
00802 for (int dir=0; dir<SpaceDim; dir++)
00803 {
00804 m_isPeriodic[dir] = b.m_isPeriodic[dir];
00805 }
00806 m_shiftIt = b.m_shiftIt;
00807 return *this;
00808 }
00809
00810 inline void
00811 ProblemDomain::shiftIt(Box& a_box, int a_shiftIndex) const
00812 {
00813 a_box.shift(m_shiftIt[a_shiftIndex]*m_domainBox.size());
00814 }
00815
00816 inline void
00817 ProblemDomain::unshiftIt(Box& a_box, int a_shiftIndex) const
00818 {
00819 a_box.shift(- m_shiftIt[a_shiftIndex]*m_domainBox.size());
00820 }
00821
00822
00823 inline
00824 const Box&
00825 ProblemDomain::domainBox() const
00826 {
00827 return m_domainBox;
00828 }
00829
00830 inline
00831 bool
00832 ProblemDomain::isPeriodic(int a_dir) const
00833 {
00834 return m_isPeriodic[a_dir];
00835 }
00836
00837 inline
00838 bool
00839 ProblemDomain::isPeriodic() const
00840 {
00841 return D_TERM6(m_isPeriodic[0], ||
00842 m_isPeriodic[1], ||
00843 m_isPeriodic[2], ||
00844 m_isPeriodic[3], ||
00845 m_isPeriodic[4], ||
00846 m_isPeriodic[5]);
00847 }
00848
00849 inline
00850 ProblemDomain&
00851 ProblemDomain::grow(int i)
00852 {
00853 m_domainBox.grow(i);
00854 return *this;
00855 }
00856
00857 inline
00858 ProblemDomain
00859 grow(const ProblemDomain& pd, int i)
00860 {
00861 ProblemDomain newPd(pd);
00862 newPd.grow(i);
00863 return newPd;
00864 }
00865
00866 inline
00867 ProblemDomain&
00868 ProblemDomain::grow(const IntVect& v)
00869 {
00870 m_domainBox.grow(v);
00871 return *this;
00872 }
00873
00874 inline
00875 ProblemDomain
00876 grow(const ProblemDomain& pd, const IntVect& v)
00877 {
00878 ProblemDomain newPd(pd);
00879 newPd.grow(v);
00880 return newPd;
00881 }
00882
00883 inline
00884 ProblemDomain&
00885 ProblemDomain::grow(int idir, int n_cell)
00886 {
00887 m_domainBox.grow(idir, n_cell);
00888 return *this;
00889 }
00890
00891 inline
00892 ProblemDomain&
00893 ProblemDomain::growLo(int idir, int n_cell)
00894 {
00895 m_domainBox.growLo(idir, n_cell);
00896 return *this;
00897 }
00898
00899 inline
00900 ProblemDomain&
00901 ProblemDomain::growHi(int idir, int n_cell)
00902 {
00903 m_domainBox.growHi(idir, n_cell);
00904 return *this;
00905 }
00906
00907 inline
00908 ShiftIterator
00909 ProblemDomain::shiftIterator() const
00910 {
00911 return m_shiftIt;
00912 }
00913
00914 inline
00915 bool
00916 ProblemDomain::isEmpty () const
00917 {
00918 return (m_domainBox.isEmpty());
00919 }
00920
00921 inline
00922 int
00923 ProblemDomain::size (const int& a_idir) const
00924 {
00925 return (m_domainBox.size(a_idir));
00926 }
00927
00928 inline
00929 IntVect
00930 ProblemDomain::size () const
00931 {
00932 return (m_domainBox.size());
00933 }
00934
00935 inline
00936 bool
00937 ProblemDomain::contains (const IntVect& p) const
00938 {
00939
00940
00941 return ( !isEmpty()
00942 && (D_TERM6((m_isPeriodic[0]
00943 || ( p[0] >= m_domainBox.smallEnd(0)
00944 && p[0] <= m_domainBox.bigEnd (0))),
00945 && (m_isPeriodic[1]
00946 || ( p[1] >= m_domainBox.smallEnd(1)
00947 && p[1] <= m_domainBox.bigEnd (1))),
00948 && (m_isPeriodic[2]
00949 || ( p[2] >= m_domainBox.smallEnd(2)
00950 && p[2] <= m_domainBox.bigEnd (2))),
00951 && (m_isPeriodic[3]
00952 || ( p[3] >= m_domainBox.smallEnd(3)
00953 && p[3] <= m_domainBox.bigEnd (3))),
00954 && (m_isPeriodic[4]
00955 || ( p[4] >= m_domainBox.smallEnd(4)
00956 && p[4] <= m_domainBox.bigEnd (4))),
00957 && (m_isPeriodic[5]
00958 || ( p[5] >= m_domainBox.smallEnd(5)
00959 && p[5] <= m_domainBox.bigEnd (5))))));
00960 }
00961
00962 inline
00963 bool
00964 ProblemDomain::image(IntVect& p) const
00965 {
00966 if(m_domainBox.contains(p)) return true;
00967 if(!contains(p)) return false;
00968
00969 D_TERM6(
00970 if(m_isPeriodic[0])
00971 {
00972 if(p[0]<m_domainBox.smallEnd(0)) p[0]+= m_domainBox.size(0);
00973 else if(p[0]>m_domainBox.bigEnd(0)) p[0]-= m_domainBox.size(0);
00974 },
00975 if(m_isPeriodic[1])
00976 {
00977 if(p[1]<m_domainBox.smallEnd(1)) p[1]+= m_domainBox.size(1);
00978 else if(p[1]>m_domainBox.bigEnd(1)) p[1]-= m_domainBox.size(1);
00979 },
00980 if(m_isPeriodic[2])
00981 {
00982 if(p[2]<m_domainBox.smallEnd(2)) p[2]+= m_domainBox.size(2);
00983 else if(p[2]>m_domainBox.bigEnd(2)) p[2]-= m_domainBox.size(2);
00984 },
00985 if(m_isPeriodic[3])
00986 {
00987 if(p[3]<m_domainBox.smallEnd(3)) p[3]+= m_domainBox.size(3);
00988 else if(p[3]>m_domainBox.bigEnd(3)) p[3]-= m_domainBox.size(3);
00989 },
00990 if(m_isPeriodic[4])
00991 {
00992 if(p[4]<m_domainBox.smallEnd(4)) p[4]+= m_domainBox.size(4);
00993 else if(p[4]>m_domainBox.bigEnd(4)) p[4]-= m_domainBox.size(4);
00994 },
00995 if(m_isPeriodic[5])
00996 {
00997 if(p[5]<m_domainBox.smallEnd(5)) p[5]+= m_domainBox.size(5);
00998 else if(p[5]>m_domainBox.bigEnd(5)) p[5]-= m_domainBox.size(5);
00999 });
01000
01001 return true;
01002 }
01003
01004 inline
01005 bool
01006 ProblemDomain::contains (const Box& b) const
01007 {
01008
01009 if (b.type() != m_domainBox.type())
01010 {
01011 return ( !isEmpty()
01012 && (D_TERM6((m_isPeriodic[0]
01013 || ( b.smallEnd(0) >= m_domainBox.smallEnd(0)
01014 && b.bigEnd (0) <= m_domainBox.bigEnd (0))), &&
01015 (m_isPeriodic[1]
01016 || ( b.smallEnd(1) >= m_domainBox.smallEnd(1)
01017 && b.bigEnd (1) <= m_domainBox.bigEnd (1))), &&
01018 (m_isPeriodic[2]
01019 || ( b.smallEnd(2) >= m_domainBox.smallEnd(2)
01020 && b.bigEnd (2) <= m_domainBox.bigEnd (2))), &&
01021 (m_isPeriodic[3]
01022 || ( b.smallEnd(3) >= m_domainBox.smallEnd(3)
01023 && b.bigEnd (3) <= m_domainBox.bigEnd (3))), &&
01024 (m_isPeriodic[4]
01025 || ( b.smallEnd(4) >= m_domainBox.smallEnd(4)
01026 && b.bigEnd (4) <= m_domainBox.bigEnd (4))), &&
01027 (m_isPeriodic[5]
01028 || ( b.smallEnd(5) >= m_domainBox.smallEnd(5)
01029 && b.bigEnd (5) <= m_domainBox.bigEnd (5))))));
01030 }
01031 else
01032 {
01033 Box domainBox = m_domainBox;
01034
01035
01036 for (int dir = 0; dir < SpaceDim; dir++)
01037 {
01038 if (b.type(dir) != domainBox.type(dir))
01039 {
01040 if (b.type(dir) == IndexType::NODE)
01041 {
01042 domainBox.surroundingNodes(dir);
01043 }
01044 else
01045 {
01046 domainBox.enclosedCells(dir);
01047 }
01048 }
01049 }
01050
01051 return ( !isEmpty()
01052 && (D_TERM6((m_isPeriodic[0]
01053 || ( b.smallEnd(0) >= domainBox.smallEnd(0)
01054 && b.bigEnd (0) <= domainBox.bigEnd (0))), &&
01055 (m_isPeriodic[1]
01056 || ( b.smallEnd(1) >= domainBox.smallEnd(1)
01057 && b.bigEnd (1) <= domainBox.bigEnd (1))), &&
01058 (m_isPeriodic[2]
01059 || ( b.smallEnd(2) >= domainBox.smallEnd(2)
01060 && b.bigEnd (2) <= domainBox.bigEnd (2))), &&
01061 (m_isPeriodic[3]
01062 || ( b.smallEnd(3) >= domainBox.smallEnd(3)
01063 && b.bigEnd (3) <= domainBox.bigEnd (3))), &&
01064 (m_isPeriodic[4]
01065 || ( b.smallEnd(4) >= domainBox.smallEnd(4)
01066 && b.bigEnd (4) <= domainBox.bigEnd (4))), &&
01067 (m_isPeriodic[5]
01068 || ( b.smallEnd(5) >= domainBox.smallEnd(5)
01069 && b.bigEnd (5) <= domainBox.bigEnd (5))))));
01070 }
01071 }
01072
01073 #endif
01074
01075 #include "NamespaceFooter.H"
01076 #endif