dune-grid
2.2.0
|
00001 #ifndef DUNE_SGRID_HH 00002 #define DUNE_SGRID_HH 00003 00004 #include <limits> 00005 #include <vector> 00006 #include <stack> 00007 00008 #include <dune/common/fvector.hh> 00009 #include <dune/common/fmatrix.hh> 00010 #include <dune/common/bigunsignedint.hh> 00011 #include <dune/common/collectivecommunication.hh> 00012 #include <dune/common/reservedvector.hh> 00013 #include <dune/geometry/genericgeometry/topologytypes.hh> 00014 #include <dune/grid/common/capabilities.hh> 00015 #include <dune/grid/common/grid.hh> 00016 #include <dune/grid/sgrid/numbering.hh> 00017 #include <dune/grid/common/indexidset.hh> 00018 00024 namespace Dune { 00025 00026 //************************************************************************ 00030 typedef double sgrid_ctype; 00031 00032 // globally define the persistent index type 00033 const int sgrid_dim_bits = 24; // bits for encoding each dimension 00034 const int sgrid_level_bits = 6; // bits for encoding level number 00035 const int sgrid_codim_bits = 4; // bits for encoding codimension 00036 00037 //************************************************************************ 00038 // forward declaration of templates 00039 00040 template<int dim, int dimworld, class GridImp> class SGeometry; 00041 template<int codim, int dim, class GridImp> class SEntity; 00042 template<int codim, class GridImp> class SEntityPointer; 00043 template<int codim, class GridImp> class SEntitySeed; 00044 template<int codim, PartitionIteratorType, class GridImp> class SLevelIterator; 00045 template<int dim, int dimworld, class ctype> class SGrid; 00046 template<class GridImp> class SIntersection; 00047 template<class GridImp> class SIntersectionIterator; 00048 template<class GridImp> class SHierarchicIterator; 00049 00050 namespace FacadeOptions 00051 { 00052 00053 template<int dim, int dimworld, class ctype, int mydim, int cdim> 00054 struct StoreGeometryReference<mydim, cdim, 00055 SGrid<dim,dimworld,ctype>, SGeometry> 00056 { 00057 static const bool v = false; 00058 }; 00059 00060 template<int dim, int dimworld, class ctype, int mydim, int cdim> 00061 struct StoreGeometryReference<mydim, cdim, 00062 const SGrid<dim,dimworld,ctype>, SGeometry> 00063 { 00064 static const bool v = false; 00065 }; 00066 00067 } 00068 00069 //************************************************************************ 00100 template<int mydim, int cdim, class GridImp> 00101 class SGeometry 00102 : public GeometryDefaultImplementation<mydim,cdim,GridImp,SGeometry> 00103 { 00104 public: 00106 typedef typename GridImp::ctype ctype; 00107 00109 GeometryType type () const 00110 { 00111 static const GeometryType cubeType(GeometryType::cube,mydim); 00112 return cubeType; 00113 } 00114 00116 bool affine() const { return true ; } 00117 00119 int corners () const 00120 { 00121 return 1<<mydim; 00122 } 00123 00125 FieldVector< ctype, cdim > corner ( const int i ) const 00126 { 00127 return c[i]; 00128 } 00129 00131 FieldVector<ctype, cdim > center ( ) const 00132 { 00133 return centroid; 00134 } 00135 00137 FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const; 00138 00140 FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& global) const; 00141 00161 ctype integrationElement (const FieldVector<ctype, mydim>& local) const 00162 { 00163 return volume(); 00164 } 00165 00167 ctype volume() const; 00168 00169 const FieldMatrix<ctype, mydim, cdim > &jacobianTransposed ( const FieldVector< ctype, mydim > &local ) const; 00170 const FieldMatrix<ctype,cdim,mydim>& jacobianInverseTransposed (const FieldVector<ctype, mydim>& local) const; 00171 00173 void print (std::ostream& ss, int indent) const; 00174 00179 void make (FieldMatrix<ctype,mydim+1,cdim>& __As); 00180 00182 SGeometry () : builtinverse(false) {}; 00183 00184 private: 00185 FieldVector<ctype, cdim> s; 00186 FieldVector<ctype, cdim> centroid; 00187 FieldMatrix<ctype,mydim,cdim> A; 00188 array<FieldVector<ctype, cdim>, 1<<mydim> c; 00189 mutable FieldMatrix<ctype,cdim,mydim> Jinv; 00190 mutable bool builtinverse; 00191 }; 00192 00194 template<int cdim, class GridImp> 00195 class SGeometry<0,cdim,GridImp> 00196 : public GeometryDefaultImplementation<0,cdim,GridImp,SGeometry> 00197 { 00198 public: 00200 typedef typename GridImp::ctype ctype; 00201 00203 GeometryType type () const 00204 { 00205 static const GeometryType cubeType(GeometryType::cube,0); 00206 return cubeType; 00207 } 00208 00210 bool affine() const { return true ; } 00211 00213 int corners () const 00214 { 00215 return 1; 00216 } 00217 00219 FieldVector<ctype, cdim > corner ( const int i ) const 00220 { 00221 return s; 00222 } 00223 00225 FieldVector<ctype, cdim > center ( ) const 00226 { 00227 return s; 00228 } 00229 00231 void print (std::ostream& ss, int indent) const; 00232 00234 void make (FieldMatrix<ctype,1,cdim>& __As); 00235 00237 FieldVector<ctype, cdim> global (const FieldVector<ctype, 0>& local) const { return corner(0); } 00238 00240 FieldVector<ctype, 0> local (const FieldVector<ctype, cdim>& global) const { return FieldVector<ctype,0> (0.0); } 00241 00263 00264 SGeometry () {}; 00265 00271 ctype volume() const 00272 { 00273 return 1; 00274 } 00275 00281 ctype integrationElement(const FieldVector<ctype, 0>& local) const { 00282 return volume(); 00283 } 00284 00285 const FieldMatrix<ctype, 0, cdim > &jacobianTransposed ( const FieldVector<ctype, 0 > &local ) const 00286 { 00287 static const FieldMatrix<ctype, 0, cdim > dummy ( ctype( 0 ) ); 00288 return dummy; 00289 } 00290 00291 const FieldMatrix<ctype,cdim,0>& jacobianInverseTransposed (const FieldVector<ctype, 0>& local) const 00292 { 00293 static const FieldMatrix<ctype,cdim,0> dummy( ctype(0) ); 00294 return dummy; 00295 } 00296 00297 protected: 00298 FieldVector<ctype, cdim> s; 00299 }; 00300 00301 template <int mydim, int cdim, class GridImp> 00302 inline std::ostream& operator<< (std::ostream& s, SGeometry<mydim,cdim,GridImp>& e) 00303 { 00304 e.print(s,0); 00305 return s; 00306 } 00307 00308 //************************************************************************ 00313 template<int codim, int dim, class GridImp, template<int,int,class> class EntityImp> 00314 class SEntityBase : 00315 public EntityDefaultImplementation<codim,dim,GridImp,EntityImp> 00316 { 00317 friend class SEntityPointer<codim,GridImp>; 00318 friend class SIntersectionIterator<GridImp>; 00319 enum { dimworld = GridImp::dimensionworld }; 00320 00321 typedef typename GridImp::Traits::template Codim< codim >::GeometryImpl GeometryImpl; 00322 00323 public: 00324 typedef typename GridImp::ctype ctype; 00325 typedef typename GridImp::template Codim<codim>::Geometry Geometry; 00326 typedef typename GridImp::PersistentIndexType PersistentIndexType; 00327 00329 int level () const 00330 { 00331 return l; 00332 } 00333 00335 int globalIndex() const; 00336 00340 SEntitySeed<codim, GridImp> seed () const { 00341 return SEntitySeed<codim, GridImp>(l, index); 00342 } 00343 00345 GeometryType type () const 00346 { 00347 static const GeometryType cubeType(GeometryType::cube,dim-codim); 00348 return cubeType; 00349 } 00350 00352 Geometry geometry () const 00353 { 00354 if (!builtgeometry) makegeometry(); 00355 00356 // return result 00357 return Geometry( geo ); 00358 } 00359 00360 PartitionType partitionType () const { return InteriorEntity; } 00361 00363 SEntityBase (GridImp* _grid, int _l, int _index) : 00364 grid(_grid), 00365 l(_l), 00366 index(_index), 00367 z(grid->z(l,index,codim)), 00368 builtgeometry(false) {} 00369 00371 SEntityBase () : 00372 builtgeometry(false) // mark geometry as not built 00373 {} 00374 00376 SEntityBase ( const SEntityBase& other ) : 00377 grid(other.grid), 00378 l(other.l), 00379 index(other.index), 00380 z(other.z), 00381 geo(), // do not copy geometry 00382 builtgeometry(false) // mark geometry as not built 00383 {} 00384 00386 void make (GridImp* _grid, int _l, int _id); 00387 00389 void make (int _l, int _id); 00390 00392 void makegeometry () const; 00393 00395 PersistentIndexType persistentIndex () const 00396 { 00397 return grid->persistentIndex(l, codim, z); 00398 } 00399 00401 int compressedIndex () const 00402 { 00403 return index; 00404 } 00405 00407 int compressedLeafIndex () const 00408 { 00409 // codim != dim -> there are no copies of entities 00410 // maxlevel -> ids are fine 00411 if (codim<dim || l==grid->maxLevel()) 00412 return compressedIndex(); 00413 00414 // this is a vertex which is not on the finest level 00415 // move coordinates up to maxlevel (multiply by 2 for each level 00416 array<int,dim> coord; 00417 for (int k=0; k<dim; k++) 00418 coord[k] = z[k]*(1<<(grid->maxLevel()-l)); 00419 00420 // compute number with respect to maxLevel 00421 return grid->n(grid->maxLevel(),coord); 00422 } 00423 00424 protected: 00425 // this is how we implement our elements 00426 GridImp* grid; 00427 int l; 00428 int index; 00429 array<int,dim> z; 00430 mutable GeometryImpl geo; 00431 mutable bool builtgeometry; 00432 }; 00433 00434 00440 template<int codim, int dim, class GridImp> 00441 class SEntity : public SEntityBase<codim,dim,GridImp,SEntity> 00442 { 00443 typedef Dune::SEntityBase<codim,dim,GridImp,Dune::SEntity> SEntityBase; 00444 friend class SEntityPointer<codim,GridImp>; 00445 friend class SIntersectionIterator<GridImp>; 00446 public: 00448 SEntity (GridImp* _grid, int _l, int _id) : 00449 SEntityBase(_grid,_l,_id) {}; 00450 }; 00451 00478 template<int dim, class GridImp> 00479 class SEntity<0,dim,GridImp> : public SEntityBase<0,dim,GridImp,SEntity> 00480 { 00481 enum { dimworld = GridImp::dimensionworld }; 00482 typedef Dune::SEntityBase<0,dim,GridImp,Dune::SEntity> SEntityBase; 00483 using SEntityBase::grid; 00484 using SEntityBase::l; 00485 using SEntityBase::index; 00486 using SEntityBase::z; 00487 00488 typedef typename GridImp::Traits::template Codim< 0 >::GeometryImpl GeometryImpl; 00489 typedef typename GridImp::Traits::template Codim< 0 >::LocalGeometryImpl LocalGeometryImpl; 00490 00491 friend class SEntityPointer<0,GridImp>; 00492 friend class SIntersectionIterator<GridImp>; 00493 00494 public: 00495 typedef typename GridImp::ctype ctype; 00496 typedef typename GridImp::template Codim<0>::Geometry Geometry; 00497 typedef typename GridImp::template Codim<0>::LocalGeometry LocalGeometry; 00498 template <int cd> 00499 struct Codim 00500 { 00501 typedef typename GridImp::template Codim<cd>::EntityPointer EntityPointer; 00502 }; 00503 typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer; 00504 typedef typename GridImp::LeafIntersectionIterator IntersectionIterator; 00505 typedef typename GridImp::HierarchicIterator HierarchicIterator; 00506 typedef typename GridImp::PersistentIndexType PersistentIndexType; 00507 00509 friend class SHierarchicIterator<GridImp>; 00510 00515 template<int cc> int count () const; 00516 00521 template<int cc> typename Codim<cc>::EntityPointer subEntity (int i) const; 00522 00524 int subCompressedIndex (int codim, int i) const 00525 { 00526 if (codim==0) return this->compressedIndex(); 00527 // compute subIndex 00528 return (this->grid)->n(this->l, this->grid->subz(this->z,i,codim)); 00529 } 00530 00534 int subCompressedLeafIndex (int codim, int i) const 00535 { 00536 if (codim==0) return this->compressedLeafIndex(); 00537 00538 assert(this->l == this->grid->maxLevel()); 00539 // compute subIndex 00540 return (this->grid)->n(this->l, this->grid->subz(this->z,i,codim)); 00541 } 00542 00544 PersistentIndexType subPersistentIndex (int codim, int i) const 00545 { 00546 if (codim==0) return this->persistentIndex(); 00547 // compute subId 00548 return this->grid->persistentIndex(this->l, codim, this->grid->subz(this->z,i,codim)); 00549 } 00550 00558 IntersectionIterator ibegin () const; 00559 IntersectionIterator ileafbegin () const; 00560 IntersectionIterator ilevelbegin () const; 00562 IntersectionIterator iend () const; 00563 IntersectionIterator ileafend () const; 00564 IntersectionIterator ilevelend () const; 00565 00571 EntityPointer father () const; 00572 00574 bool hasFather () const 00575 { 00576 return (this->level()>0); 00577 } 00578 00580 bool isLeaf () const 00581 { 00582 return ( this->grid->maxLevel() == this->level() ); 00583 } 00584 00596 LocalGeometry geometryInFather () const; 00597 00604 HierarchicIterator hbegin (int maxLevel) const; 00605 00607 HierarchicIterator hend (int maxLevel) const; 00608 00609 // members specific to SEntity 00611 SEntity (GridImp* _grid, int _l, int _index) : 00612 SEntityBase(_grid,_l,_index), 00613 built_father(false) 00614 {} 00615 00616 SEntity (const SEntity& other ) : 00617 SEntityBase(other.grid, other.l, other.index ), 00618 built_father(false) 00619 {} 00620 00622 void make (GridImp* _grid, int _l, int _id) 00623 { 00624 SEntityBase::make(_grid,_l,_id); 00625 built_father = false; 00626 } 00627 00629 void make (int _l, int _id) 00630 { 00631 SEntityBase::make(_l,_id); 00632 built_father = false; 00633 } 00634 00635 private: 00636 00637 SEntity(); 00638 00639 mutable bool built_father; 00640 mutable int father_index; 00641 mutable LocalGeometryImpl in_father_local; 00642 void make_father() const; 00643 }; 00644 00645 00646 //************************************************************************ 00654 struct SHierarchicStackElem { 00655 int l; 00656 int index; 00657 SHierarchicStackElem () : l(-1), index(-1) {} 00658 SHierarchicStackElem (int _l, int _index) {l=_l; index=_index;} 00659 bool operator== (const SHierarchicStackElem& s) const {return !operator!=(s);} 00660 bool operator!= (const SHierarchicStackElem& s) const {return l!=s.l || index!=s.index;} 00661 }; 00662 00663 template<class GridImp> 00664 class SHierarchicIterator : 00665 public Dune::SEntityPointer <0,GridImp> 00666 { 00667 friend class SHierarchicIterator<const GridImp>; 00668 enum { dim = GridImp::dimension }; 00669 enum { dimworld = GridImp::dimensionworld }; 00670 typedef Dune::SEntityPointer<0,GridImp> SEntityPointer; 00671 using SEntityPointer::realEntity; 00672 using SEntityPointer::grid; 00673 using SEntityPointer::l; 00674 using SEntityPointer::index; 00675 public: 00676 typedef typename GridImp::template Codim<0>::Entity Entity; 00677 typedef typename GridImp::ctype ctype; 00678 00680 void increment(); 00681 00688 SHierarchicIterator (GridImp* _grid, 00689 const Dune::SEntity<0,GridImp::dimension,GridImp>& _e, 00690 int _maxLevel, bool makeend) : 00691 SEntityPointer(_grid,_e.level(),_e.compressedIndex()) 00692 { 00693 // without sons, we are done 00694 // (the end iterator is equal to the calling iterator) 00695 if (makeend) return; 00696 00697 // remember element where begin has been called 00698 orig_l = this->entity().level(); 00699 orig_index = _grid->getRealImplementation(this->entity()).compressedIndex(); 00700 00701 // push original element on stack 00702 SHierarchicStackElem originalElement(orig_l, orig_index); 00703 stack.push(originalElement); 00704 00705 // compute maxLevel 00706 maxLevel = std::min(_maxLevel,this->grid->maxLevel()); 00707 00708 // ok, push all the sons as well 00709 push_sons(orig_l,orig_index); 00710 00711 // and pop the first son 00712 increment(); 00713 } 00714 00715 private: 00716 int maxLevel; 00717 int orig_l, orig_index; 00718 00720 std::stack<SHierarchicStackElem, Dune::ReservedVector<SHierarchicStackElem,GridImp::MAXL> > stack; 00721 00722 void push_sons (int level, int fatherid); 00723 }; 00724 00725 //************************************************************************ 00732 template<class GridImp> 00733 class SIntersectionIterator 00734 { 00735 enum { dim=GridImp::dimension }; 00736 enum { dimworld=GridImp::dimensionworld }; 00737 00738 typedef typename GridImp::Traits::template Codim< 1 >::GeometryImpl GeometryImpl; 00739 typedef typename GridImp::Traits::template Codim< 1 >::LocalGeometryImpl LocalGeometryImpl; 00740 00741 public: 00742 typedef typename GridImp::template Codim<0>::Entity Entity; 00743 typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer; 00744 typedef typename GridImp::template Codim<1>::Geometry Geometry; 00745 typedef typename GridImp::template Codim<1>::LocalGeometry LocalGeometry; 00746 typedef Dune::SIntersection<GridImp> IntersectionImp; 00747 typedef Dune::Intersection<const GridImp, Dune::SIntersection> Intersection; 00749 enum { dimension=dim }; 00751 enum { dimensionworld=dimworld }; 00753 typedef typename GridImp::ctype ctype; 00754 00756 bool equals(const SIntersectionIterator<GridImp>& i) const; 00758 void increment(); 00759 00761 const Intersection & dereference() const 00762 { 00763 return intersection; 00764 } 00765 00768 EntityPointer inside() const; 00769 00772 EntityPointer outside() const; 00773 00775 bool boundary () const; 00776 00778 bool conforming () const; 00779 00780 int boundaryId () const { 00781 if (boundary()) return count + 1; 00782 return 0; 00783 }; 00784 00785 int boundarySegmentIndex () const { 00786 if (boundary()) 00787 return grid->boundarySegmentIndex(self.level(), count, zred); 00788 return -1; 00789 }; 00790 00792 bool neighbor () const; 00793 00795 FieldVector<ctype, GridImp::dimensionworld> outerNormal (const FieldVector<ctype, GridImp::dimension-1>& local) const 00796 { 00797 return unitOuterNormal(local); 00798 } 00800 FieldVector<ctype, GridImp::dimensionworld> unitOuterNormal (const FieldVector<ctype, GridImp::dimension-1>& local) const 00801 { 00802 return centerUnitOuterNormal(); 00803 } 00805 FieldVector<ctype, GridImp::dimensionworld> centerUnitOuterNormal () const 00806 { 00807 // while we are at it, compute normal direction 00808 FieldVector<ctype, dimworld> normal(0.0); 00809 if (count%2) 00810 normal[count/2] = 1.0; // odd 00811 else 00812 normal[count/2] = -1.0; // even 00813 00814 return normal; 00815 } 00817 FieldVector<ctype, GridImp::dimensionworld> integrationOuterNormal (const FieldVector<ctype, GridImp::dimension-1>& local) const 00818 { 00819 FieldVector<ctype, dimworld> n = unitOuterNormal(local); 00820 n *= geometry().integrationElement(local); 00821 return n; 00822 } 00823 00827 LocalGeometry geometryInInside () const; 00831 LocalGeometry geometryInOutside () const; 00835 Geometry geometry () const; 00836 00838 GeometryType type () const 00839 { 00840 static const GeometryType cubeType(GeometryType::cube,dim-1); 00841 return cubeType; 00842 } 00843 00845 int indexInInside () const; 00847 int indexInOutside () const; 00848 00850 SIntersectionIterator (GridImp* _grid, const SEntity<0,dim,GridImp >* _self, int _count) : 00851 self(*_self), ne(self), grid(_grid), 00852 partition(_grid->partition(grid->getRealImplementation(ne).l,_self->z)), 00853 zred(_grid->compress(grid->getRealImplementation(ne).l,_self->z)), 00854 intersection(IntersectionImp(*this)) 00855 { 00856 // make neighbor 00857 make(_count); 00858 } 00859 00860 SIntersectionIterator (const SIntersectionIterator & other) : 00861 self(other.self), ne(other.ne), grid(other.grid), 00862 partition(other.partition), zred(other.zred), 00863 count(other.count), valid_count(other.valid_count), 00864 valid_nb(other.valid_nb), is_on_boundary(other.is_on_boundary), 00865 built_intersections(false), 00866 intersection(IntersectionImp(*this)) 00867 { 00868 } 00869 00871 SIntersectionIterator& operator = (const SIntersectionIterator& other) 00872 { 00873 /* We can't assign the grid */ 00874 assert(grid == other.grid); 00875 00876 /* Assign data from other */ 00877 self = other.self; 00878 ne = other.ne; 00879 partition = other.partition; 00880 zred = other.zred; 00881 count = other.count; 00882 valid_count = other.valid_count; 00883 valid_nb = other.valid_nb; 00884 is_on_boundary = other.is_on_boundary; 00885 00886 /* mark cached data as invalid */ 00887 built_intersections = false; 00888 00889 return *this; 00890 } 00891 00892 private: 00893 void make (int _count) const; 00894 void makeintersections () const; 00895 EntityPointer self; 00896 mutable EntityPointer ne; 00897 const GridImp * grid; 00898 int partition; 00899 array<int,dim> zred; 00900 mutable int count; 00901 mutable bool valid_count; 00902 mutable bool valid_nb; 00903 mutable bool is_on_boundary; 00904 mutable bool built_intersections; 00905 mutable LocalGeometryImpl is_self_local; 00906 mutable GeometryImpl is_global; 00907 mutable LocalGeometryImpl is_nb_local; 00908 Intersection intersection; 00909 }; 00910 00911 template<class GridImp> 00912 class SIntersection 00913 { 00914 enum { dim=GridImp::dimension }; 00915 enum { dimworld=GridImp::dimensionworld }; 00916 public: 00917 typedef typename GridImp::template Codim<0>::Entity Entity; 00918 typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer; 00919 typedef typename GridImp::template Codim<1>::Geometry Geometry; 00920 typedef typename Geometry::LocalCoordinate LocalCoordinate; 00921 typedef typename Geometry::GlobalCoordinate GlobalCoordinate; 00922 typedef typename GridImp::template Codim<1>::LocalGeometry LocalGeometry; 00923 typedef Dune::Intersection<const GridImp, Dune::SIntersectionIterator> Intersection; 00925 enum { dimension=dim }; 00927 enum { dimensionworld=dimworld }; 00929 typedef typename GridImp::ctype ctype; 00930 00931 bool boundary () const 00932 { 00933 return is.boundary(); 00934 } 00935 00937 int boundaryId () const 00938 { 00939 return is.boundaryId(); 00940 } 00941 00943 size_t boundarySegmentIndex () const 00944 { 00945 return is.boundarySegmentIndex(); 00946 } 00947 00949 bool neighbor () const 00950 { 00951 return is.neighbor(); 00952 } 00953 00955 EntityPointer inside() const 00956 { 00957 return is.inside(); 00958 } 00959 00961 EntityPointer outside() const 00962 { 00963 return is.outside(); 00964 } 00965 00967 bool conforming () const 00968 { 00969 return is.conforming(); 00970 } 00971 00973 LocalGeometry geometryInInside () const 00974 { 00975 return is.geometryInInside(); 00976 } 00977 00979 LocalGeometry geometryInOutside () const 00980 { 00981 return is.geometryInOutside(); 00982 } 00983 00985 Geometry geometry () const 00986 { 00987 return is.geometry(); 00988 } 00989 00991 GeometryType type () const 00992 { 00993 return is.type(); 00994 } 00995 00997 int indexInInside () const 00998 { 00999 return is.indexInInside(); 01000 } 01001 01003 int indexInOutside () const 01004 { 01005 return is.indexInOutside(); 01006 } 01007 01009 GlobalCoordinate outerNormal (const LocalCoordinate& local) const 01010 { 01011 return is.outerNormal(local); 01012 } 01013 01015 GlobalCoordinate integrationOuterNormal (const LocalCoordinate& local) const 01016 { 01017 return is.integrationOuterNormal(local); 01018 } 01019 01021 GlobalCoordinate unitOuterNormal (const LocalCoordinate& local) const 01022 { 01023 return is.unitOuterNormal(local); 01024 } 01025 01027 GlobalCoordinate centerUnitOuterNormal () const 01028 { 01029 return is.centerUnitOuterNormal(); 01030 } 01031 01033 SIntersection (const SIntersectionIterator<GridImp> & is_) : is(is_) {} 01034 01035 private: 01036 #ifndef DOXYGEN // doxygen can't handle this recursive usage 01037 const SIntersectionIterator<GridImp> & is; 01038 #endif 01039 }; 01040 01041 //************************************************************************ 01042 01046 template <class T> 01047 class AutoPtrStack : public std::stack<T*> 01048 { 01049 public: 01050 ~AutoPtrStack() 01051 { 01052 while(! this->empty()) 01053 { 01054 T* e = this->top(); 01055 delete e; 01056 this->pop(); 01057 } 01058 } 01059 }; 01060 01063 template<int codim, class GridImp> 01064 class SEntityPointer 01065 { 01066 enum { dim = GridImp::dimension }; 01067 friend class SIntersectionIterator<GridImp>; 01068 public: 01069 typedef SEntityPointer<codim,GridImp> EntityPointerImp; 01070 typedef typename GridImp::template Codim<codim>::Entity Entity; 01072 enum { codimension = codim }; 01073 01075 bool equals(const SEntityPointer<codim,GridImp>& i) const; 01077 Entity& dereference() const; 01079 int level () const; 01080 01082 SEntityPointer (GridImp * _grid, int _l, int _index) : 01083 grid(_grid), l(_l), index(_index), 01084 e(0) 01085 { 01086 } 01087 01089 SEntityPointer (const SEntity<codim,dim,GridImp> & _e) : 01090 grid(_e.grid), l(_e.l), index(_e.index), 01091 e(0) 01092 { 01093 } 01094 01096 SEntityPointer (const SEntityPointer<codim,GridImp>& other) : 01097 grid(other.grid), l(other.l), index(other.index), 01098 e( 0 ) 01099 { 01100 } 01101 01103 ~SEntityPointer() 01104 { 01105 if( e ) 01106 enStack().push( e ); 01107 #ifndef NDEBUG 01108 index = -1; 01109 #endif 01110 } 01111 01113 SEntityPointer& operator = (const SEntityPointer& other) 01114 { 01115 grid = other.grid; 01116 l = other.l; 01117 index = other.index; 01118 01119 // free current entity 01120 if( e ) 01121 enStack().push( e ); 01122 e = 0; 01123 01124 return *this; 01125 } 01126 01127 protected: 01128 inline SEntity<codim,dim,GridImp>& realEntity() const 01129 { 01130 return grid->getRealImplementation(entity()); 01131 } 01132 01133 inline Entity& entity() const 01134 { 01135 if( ! e ) 01136 { 01137 e = getEntity( grid, l, index ); 01138 } 01139 return *e; 01140 } 01141 01142 typedef AutoPtrStack< Entity > EntityStackType; 01143 static inline EntityStackType& enStack() 01144 { 01145 static EntityStackType eStack; 01146 return eStack; 01147 } 01148 01149 inline Entity* getEntity(GridImp* _grid, int _l, int _id ) const 01150 { 01151 // get stack reference 01152 EntityStackType& enSt = enStack(); 01153 01154 if( enSt.empty() ) 01155 { 01156 return (new Entity(SEntity<codim,dim,GridImp>(_grid, _l, _id))); 01157 } 01158 else 01159 { 01160 Entity* e = enSt.top(); 01161 enSt.pop(); 01162 grid->getRealImplementation(*e).make(_grid, _l,_id); 01163 return e; 01164 } 01165 } 01166 01167 GridImp* grid; 01168 int l; 01169 mutable int index; 01170 mutable Entity* e; 01171 }; 01172 01175 template<int codim, class GridImp> 01176 class SEntitySeed 01177 { 01178 enum { dim = GridImp::dimension }; 01179 public: 01180 enum { codimension = codim }; 01181 01183 SEntitySeed (int l, int index) : 01184 _l(l), _index(index) 01185 {} 01186 01187 int level () const { return this->_l; } 01188 int index () const { return this->_index; } 01189 01190 private: 01191 int _l; 01192 int _index; 01193 }; 01194 01195 //************************************************************************ 01196 01197 01200 template<int codim, PartitionIteratorType pitype, class GridImp> 01201 class SLevelIterator : 01202 public Dune::SEntityPointer <codim,GridImp> 01203 { 01204 friend class SLevelIterator<codim, pitype,const GridImp>; 01205 enum { dim = GridImp::dimension }; 01206 typedef Dune::SEntityPointer<codim,GridImp> SEntityPointer; 01207 using SEntityPointer::realEntity; 01208 using SEntityPointer::l; 01209 using SEntityPointer::index; 01210 public: 01211 typedef typename GridImp::template Codim<codim>::Entity Entity; 01212 01214 void increment(); 01215 01217 SLevelIterator (GridImp * _grid, int _l, int _id) : 01218 SEntityPointer(_grid,_l,_id) {} 01219 }; 01220 01221 01222 //======================================================================== 01227 //======================================================================== 01228 01229 template<class GridImp> 01230 class SGridLevelIndexSet : public IndexSet<GridImp,SGridLevelIndexSet<GridImp> > 01231 { 01232 typedef SGridLevelIndexSet< GridImp > This; 01233 typedef IndexSet< GridImp, This > Base; 01234 01235 enum { dim = GridImp::dimension }; 01236 01237 public: 01238 01240 SGridLevelIndexSet ( const GridImp &g, int l ) 01241 : grid( g ), 01242 level( l ) 01243 { 01244 // TODO move list of geometrytypes to grid, can be computed static (singleton) 01245 // contains a single element type; 01246 for (int codim=0; codim<=GridImp::dimension; codim++) 01247 mytypes[codim].push_back(GeometryType(GeometryType::cube,GridImp::dimension-codim)); 01248 } 01249 01251 template<int cd> 01252 int index (const typename GridImp::Traits::template Codim<cd>::Entity& e) const 01253 { 01254 return grid.getRealImplementation(e).compressedIndex(); 01255 } 01256 01257 template< int cc > 01258 int subIndex ( const typename GridImp::Traits::template Codim< cc >::Entity &e, 01259 int i, unsigned int codim ) const 01260 { 01261 if( cc == 0 ) 01262 return grid.getRealImplementation(e).subCompressedIndex(codim, i); 01263 else 01264 DUNE_THROW( NotImplemented, "subIndex for higher codimension entity not implemented for SGrid." ); 01265 } 01266 01267 // return true if the given entity is contained in \f$E\f$. 01268 template< class EntityType > 01269 bool contains ( const EntityType &e ) const 01270 { 01271 return (e.level() == level); 01272 } 01273 01275 int size (GeometryType type) const 01276 { 01277 return grid.size( level, type ); 01278 } 01279 01281 int size (int codim) const 01282 { 01283 return grid.size( level, codim ); 01284 } 01285 01287 const std::vector<GeometryType>& geomTypes (int codim) const 01288 { 01289 return mytypes[codim]; 01290 } 01291 01292 private: 01293 const GridImp& grid; 01294 int level; 01295 std::vector<GeometryType> mytypes[GridImp::dimension+1]; 01296 }; 01297 01298 // Leaf Index Set 01299 01300 template<class GridImp> 01301 class SGridLeafIndexSet : public IndexSet<GridImp,SGridLeafIndexSet<GridImp> > 01302 { 01303 typedef SGridLeafIndexSet< GridImp > This; 01304 typedef IndexSet< GridImp, This > Base; 01305 01306 enum { dim = GridImp::dimension }; 01307 01308 public: 01309 01311 explicit SGridLeafIndexSet ( const GridImp &g ) 01312 : grid( g ) 01313 { 01314 // TODO move list of geometrytypes to grid, can be computed static (singleton) 01315 // contains a single element type; 01316 for (int codim=0; codim<=dim; codim++) 01317 mytypes[codim].push_back(GeometryType(GeometryType::cube,dim-codim)); 01318 } 01319 01321 /* 01322 We use the remove_const to extract the Type from the mutable class, 01323 because the const class is not instantiated yet. 01324 */ 01325 template<int cd> 01326 int index (const typename remove_const<GridImp>::type::Traits::template Codim<cd>::Entity& e) const 01327 { 01328 return grid.getRealImplementation(e).compressedLeafIndex(); 01329 } 01330 01331 template< int cc > 01332 int subIndex ( const typename GridImp::Traits::template Codim< cc >::Entity &e, 01333 int i, unsigned int codim ) const 01334 { 01335 if( cc == 0 ) 01336 return grid.getRealImplementation(e).subCompressedIndex(codim, i); 01337 else 01338 DUNE_THROW( NotImplemented, "subIndex for higher codimension entity not implemented for SGrid." ); 01339 } 01340 01342 int size (GeometryType type) const 01343 { 01344 return grid.size( grid.maxLevel(), type ); 01345 } 01346 01348 int size (int codim) const 01349 { 01350 return grid.size( grid.maxLevel(), codim ); 01351 } 01352 01353 // return true if the given entity is contained in \f$E\f$. 01354 template< class EntityType > 01355 bool contains ( const EntityType &e ) const 01356 { 01357 return (e.level() == grid.maxLevel()); 01358 } 01359 01361 const std::vector<GeometryType>& geomTypes (int codim) const 01362 { 01363 return mytypes[codim]; 01364 } 01365 01366 private: 01367 const GridImp& grid; 01368 std::vector<GeometryType> mytypes[dim+1]; 01369 }; 01370 01371 01372 //======================================================================== 01377 //======================================================================== 01378 01379 template<class GridImp> 01380 class SGridGlobalIdSet : 01381 public IdSet<GridImp,SGridGlobalIdSet<GridImp>, typename remove_const<GridImp>::type::PersistentIndexType> 01382 /* 01383 We used the remove_const to extract the Type from the mutable class, 01384 because the const class is not instantiated yet. 01385 */ 01386 { 01387 typedef SGridGlobalIdSet< GridImp > This; 01388 01389 public: 01390 01393 using IdSet<GridImp,SGridGlobalIdSet<GridImp>, typename remove_const<GridImp>::type::PersistentIndexType>::subId; 01394 01396 /* 01397 We use the remove_const to extract the Type from the mutable class, 01398 because the const class is not instantiated yet. 01399 */ 01400 typedef typename remove_const<GridImp>::type::PersistentIndexType IdType; 01401 01403 explicit SGridGlobalIdSet ( const GridImp &g ) 01404 : grid( g ) 01405 {} 01406 01408 /* 01409 We use the remove_const to extract the Type from the mutable class, 01410 because the const class is not instantiated yet. 01411 */ 01412 template<int cd> 01413 IdType id (const typename remove_const<GridImp>::type::Traits::template Codim<cd>::Entity& e) const 01414 { 01415 return grid.getRealImplementation(e).persistentIndex(); 01416 } 01417 01419 /* 01420 We use the remove_const to extract the Type from the mutable class, 01421 because the const class is not instantiated yet. 01422 */ 01423 IdType subId ( const typename remove_const< GridImp >::type::Traits::template Codim< 0 >::Entity &e, 01424 int i, unsigned int codim ) const 01425 { 01426 return grid.getRealImplementation(e).subPersistentIndex(codim, i); 01427 } 01428 01429 private: 01430 const GridImp& grid; 01431 }; 01432 01433 01434 template<int dim, int dimworld, class ctype> 01435 struct SGridFamily 01436 { 01437 typedef GridTraits<dim,dimworld,Dune::SGrid<dim,dimworld,ctype>, 01438 SGeometry,SEntity, 01439 SEntityPointer,SLevelIterator, 01440 SIntersection, // leaf intersection 01441 SIntersection, // level intersection 01442 SIntersectionIterator, // leaf intersection iter 01443 SIntersectionIterator, // level intersection iter 01444 SHierarchicIterator, 01445 SLevelIterator, 01446 SGridLevelIndexSet<const SGrid<dim,dimworld,ctype> >, 01447 SGridLeafIndexSet<const SGrid<dim,dimworld,ctype> >, 01448 SGridGlobalIdSet<const SGrid<dim,dimworld,ctype> >, 01449 bigunsignedint<dim*sgrid_dim_bits+sgrid_level_bits+sgrid_codim_bits>, 01450 SGridGlobalIdSet<const SGrid<dim,dimworld,ctype> >, 01451 bigunsignedint<dim*sgrid_dim_bits+sgrid_level_bits+sgrid_codim_bits>, 01452 CollectiveCommunication<Dune::SGrid<dim,dimworld,ctype> >, 01453 DefaultLevelGridViewTraits, DefaultLeafGridViewTraits, 01454 SEntitySeed> 01455 Traits; 01456 }; 01457 01458 01459 //************************************************************************ 01513 template<int dim, int dimworld, typename _ctype = sgrid_ctype> 01514 class SGrid : public GridDefaultImplementation <dim,dimworld,_ctype,SGridFamily<dim,dimworld,_ctype> > 01515 { 01516 public: 01517 typedef SGridFamily<dim,dimworld,_ctype> GridFamily; 01518 typedef bigunsignedint<dim*sgrid_dim_bits+sgrid_level_bits+sgrid_codim_bits> PersistentIndexType; 01519 01520 // need for friend declarations in entity 01521 typedef SGridLevelIndexSet<SGrid<dim,dimworld> > LevelIndexSetType; 01522 typedef SGridLeafIndexSet<SGrid<dim,dimworld> > LeafIndexSetType; 01523 typedef SGridGlobalIdSet<SGrid<dim,dimworld> > GlobalIdSetType; 01524 01525 typedef typename SGridFamily<dim,dimworld,_ctype>::Traits Traits; 01526 01528 enum { MAXL=32 }; 01529 01531 typedef _ctype ctype; 01532 01533 // constructors 01534 01542 SGrid (const int * const N_, const ctype * const H_); 01543 01551 SGrid (const int * const N_, const ctype * const L_, const ctype * const H_); 01552 01562 SGrid (FieldVector<int,dim> N_, FieldVector<ctype,dim> L_, FieldVector<ctype,dim> H_); 01563 01565 SGrid (); 01566 01568 ~SGrid (); 01569 01572 int maxLevel() const; 01573 01575 template<int cd, PartitionIteratorType pitype> 01576 typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator lbegin (int level) const; 01577 01579 template<int cd, PartitionIteratorType pitype> 01580 typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator lend (int level) const; 01581 01583 template<int cd> 01584 typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator lbegin (int level) const 01585 { 01586 return lbegin<cd,All_Partition>(level); 01587 } 01588 01590 template<int cd> 01591 typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator lend (int level) const 01592 { 01593 return lend<cd,All_Partition>(level); 01594 } 01595 01597 template<int cd, PartitionIteratorType pitype> 01598 typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator leafbegin () const; 01599 01601 template<int cd, PartitionIteratorType pitype> 01602 typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator leafend () const; 01603 01605 template<int cd> 01606 typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator leafbegin () const 01607 { 01608 return leafbegin<cd,All_Partition>(); 01609 } 01610 01612 template<int cd> 01613 typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator leafend () const 01614 { 01615 return leafend<cd,All_Partition>(); 01616 } 01617 01618 // \brief obtain EntityPointer from EntitySeed. */ 01619 template <typename Seed> 01620 typename Traits::template Codim<Seed::codimension>::EntityPointer 01621 entityPointer(const Seed& seed) const 01622 { 01623 enum { codim = Seed::codimension }; 01624 return SEntityPointer<codim,const SGrid<dim,dimworld> >(this,seed.level(),seed.index()); 01625 } 01626 01640 template<class T, template<class> class P, int codim> 01641 void communicate (T& t, InterfaceType iftype, CommunicationDirection dir, int level) 01642 { 01643 // SGrid is sequential and has no periodic boundaries, so do nothing ... 01644 return; 01645 } 01646 01648 int size (int level, int codim) const; 01649 01651 int size (int codim) const 01652 { 01653 return size(maxLevel(),codim); 01654 } 01655 01657 int size (int level, GeometryType type) const 01658 { 01659 return (type.isCube()) ? size(level,dim-type.dim()) : 0; 01660 } 01661 01663 int size (GeometryType type) const 01664 { 01665 return size(maxLevel(),type); 01666 } 01667 01669 size_t numBoundarySegments () const 01670 { 01671 return boundarysize; 01672 } 01673 01675 int global_size (int codim) const; 01676 01678 int overlapSize (int level, int codim) 01679 { 01680 return 0; 01681 } 01682 01684 int ghostSize (int level, int codim) 01685 { 01686 return 0; 01687 } 01688 01689 // these are all members specific to sgrid 01690 01692 void globalRefine (int refCount); 01693 01695 const array<int, dim>& dims(int level) const { 01696 return N[level]; 01697 } 01698 01700 const FieldVector<ctype, dimworld>& lowerLeft() const { 01701 return low; 01702 } 01703 01705 FieldVector<ctype, dimworld> upperRight() const { 01706 return H; 01707 } 01708 01710 bool adapt () 01711 { 01712 globalRefine(1); 01713 return true; 01714 } 01715 01716 // The new index sets from DDM 11.07.2005 01717 const typename Traits::GlobalIdSet& globalIdSet() const 01718 { 01719 return *theglobalidset; 01720 } 01721 01722 const typename Traits::LocalIdSet& localIdSet() const 01723 { 01724 return *theglobalidset; 01725 } 01726 01727 const typename Traits::LevelIndexSet& levelIndexSet(int level) const 01728 { 01729 assert(level>=0 && level<=maxLevel()); 01730 return *(indexsets[level]); 01731 } 01732 01733 const typename Traits::LeafIndexSet& leafIndexSet() const 01734 { 01735 return *theleafindexset; 01736 } 01737 01742 template<class DataHandle> 01743 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int level) const 01744 { 01745 } 01746 01747 template<class DataHandle> 01748 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const 01749 { 01750 } 01751 01752 const CollectiveCommunication<SGrid>& comm () const 01753 { 01754 return ccobj; 01755 } 01756 01758 int overlapSize (int level, int codim) const 01759 { 01760 return 0; 01761 } 01762 01764 int overlapSize (int codim) const 01765 { 01766 return 0; 01767 } 01768 01770 int ghostSize (int level, int codim) const 01771 { 01772 return 0; 01773 } 01774 01776 int ghostSize (int codim) const 01777 { 01778 return 0; 01779 } 01780 01781 /* 01782 @} 01783 */ 01784 01785 private: 01786 /* 01787 Make associated classes friends to grant access to the real entity 01788 */ 01789 friend class Dune::SGridLevelIndexSet<Dune::SGrid<dim,dimworld> >; 01790 friend class Dune::SGridLeafIndexSet<Dune::SGrid<dim,dimworld> >; 01791 friend class Dune::SGridGlobalIdSet<Dune::SGrid<dim,dimworld> >; 01792 friend class Dune::SIntersectionIterator<Dune::SGrid<dim,dimworld> >; 01793 friend class Dune::SHierarchicIterator<Dune::SGrid<dim,dimworld> >; 01794 friend class Dune::SEntity<0,dim,Dune::SGrid<dim,dimworld> >; 01795 01796 friend class Dune::SGridLevelIndexSet<const Dune::SGrid<dim,dimworld> >; 01797 friend class Dune::SGridLeafIndexSet<const Dune::SGrid<dim,dimworld> >; 01798 friend class Dune::SGridGlobalIdSet<const Dune::SGrid<dim,dimworld> >; 01799 friend class Dune::SIntersectionIterator<const Dune::SGrid<dim,dimworld> >; 01800 friend class Dune::SHierarchicIterator<const Dune::SGrid<dim,dimworld> >; 01801 friend class Dune::SEntity<0,dim,const Dune::SGrid<dim,dimworld> >; 01802 01803 template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_> 01804 friend class Dune::SEntityBase; 01805 01806 template<int codim_, class GridImp_> 01807 friend class Dune::SEntityPointer; 01808 01809 template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_> 01810 friend class Entity; 01811 01813 FieldVector<ctype, dimworld> pos (int level, array<int,dim>& z) const; 01814 01816 int calc_codim (int level, const array<int,dim>& z) const; 01817 01819 int n (int level, const array<int,dim>& z) const; 01820 01822 array<int,dim> z (int level, int i, int codim) const; 01823 01825 array<int,dim> subz (const array<int,dim> & z, int i, int codim) const; 01826 01828 array<int,dim> compress (int level, const array<int,dim>& z) const; 01829 01831 array<int,dim> expand (int level, const array<int,dim>& r, int b) const; 01832 01836 int partition (int level, const array<int,dim>& z) const; 01837 01839 bool exists (int level, const array<int,dim>& zred) const; 01840 01841 // compute boundary segment index for a given zentity and a face 01842 int boundarySegmentIndex (int l, int face, const array<int,dim> & zentity) const 01843 { 01844 array<int,dim-1> zface; 01845 int dir = face/2; 01846 int side = face%2; 01847 // compute z inside the global face 01848 for (int i=0; i<dir; i++) zface[i] = zentity[i]/(1<<l); 01849 for (int i=dir+1; i<dim; i++) zface[i-1] = zentity[i]/(1<<l); 01850 zface = boundarymapper[dir].expand(zface, 0); 01851 // compute index in the face 01852 int index = boundarymapper[dir].n(zface); 01853 // compute offset 01854 for (int i=0; i<dir; i++) 01855 index += 2*boundarymapper[i].elements(0); 01856 index += side*boundarymapper[dir].elements(0); 01857 return index; 01858 } 01859 01860 // compute persistent index for a given zentity 01861 PersistentIndexType persistentIndex (int l, int codim, const array<int,dim> & zentity) const 01862 { 01863 if (codim!=dim) 01864 { 01865 // encode codim, this would actually not be necessary 01866 // because z is unique in codim 01867 PersistentIndexType id(codim); 01868 01869 // encode level 01870 id = id << sgrid_level_bits; 01871 id = id+PersistentIndexType(l); 01872 01873 // encode coordinates 01874 for (int i=dim-1; i>=0; i--) 01875 { 01876 id = id << sgrid_dim_bits; 01877 id = id+PersistentIndexType(zentity[i]); 01878 } 01879 01880 return id; 01881 } 01882 else 01883 { 01884 // determine min number of trailing zeroes 01885 // consider that z is on the doubled grid ! 01886 int trailing = 1000; 01887 for (int i=0; i<dim; i++) 01888 { 01889 // count trailing zeros 01890 int zeros = 0; 01891 for (int j=0; j<l; j++) 01892 if (zentity[i]&(1<<(j+1))) 01893 break; 01894 else 01895 zeros++; 01896 trailing = std::min(trailing,zeros); 01897 } 01898 01899 // determine the level of this vertex 01900 int level = l-trailing; 01901 01902 // encode codim 01903 PersistentIndexType id(dim); 01904 01905 // encode level 01906 id = id << sgrid_level_bits; 01907 id = id+PersistentIndexType(level); 01908 01909 // encode coordinates 01910 for (int i=dim-1; i>=0; i--) 01911 { 01912 id = id << sgrid_dim_bits; 01913 id = id+PersistentIndexType(zentity[i]>>trailing); 01914 } 01915 01916 return id; 01917 } 01918 } 01919 01920 // disable copy and assign 01921 SGrid(const SGrid &) {}; 01922 SGrid & operator = (const SGrid &) { return *this; }; 01923 // generate SGrid 01924 void makeSGrid (const array<int,dim>& N_, const FieldVector<ctype, dim>& L_, const FieldVector<ctype, dim>& H_); 01925 01926 /* 01927 internal data 01928 */ 01929 CollectiveCommunication<SGrid> ccobj; 01930 01931 ReservedVector<SGridLevelIndexSet<const SGrid<dim,dimworld> >*, MAXL> indexsets; 01932 SGridLeafIndexSet<const SGrid<dim,dimworld> > *theleafindexset; 01933 SGridGlobalIdSet<const SGrid<dim,dimworld> > *theglobalidset; 01934 01935 int L; // number of levels in hierarchic mesh 0<=level<L 01936 FieldVector<ctype, dim> low; // lower left corner of the grid 01937 FieldVector<ctype, dim> H; // length of cube per direction 01938 array<int,dim> *N; // number of elements per direction 01939 FieldVector<ctype, dim> *h; // mesh size per direction 01940 mutable CubeMapper<dim> *mapper;// a mapper for each level 01941 01942 // boundary segement index set 01943 array<CubeMapper<dim-1>, dim> boundarymapper; // a mapper for each coarse grid face 01944 int boundarysize; 01945 01946 // faster implementation of subIndex 01947 mutable array <int,dim> zrefStatic; // for subIndex of SEntity 01948 mutable array <int,dim> zentityStatic; // for subIndex of SEntity 01949 }; 01950 01951 namespace Capabilities 01952 { 01953 01965 template<int dim, int dimw> 01966 struct hasSingleGeometryType< SGrid<dim,dimw> > 01967 { 01968 static const bool v = true; 01969 static const unsigned int topologyId = GenericGeometry :: CubeTopology< dim > :: type :: id ; 01970 }; 01971 01975 template<int dim, int dimw> 01976 struct isCartesian< SGrid<dim,dimw> > 01977 { 01978 static const bool v = true; 01979 }; 01980 01984 template<int dim, int dimw, int cdim> 01985 struct hasEntity< SGrid<dim,dimw>, cdim> 01986 { 01987 static const bool v = true; 01988 }; 01989 01993 template<int dim, int dimw> 01994 struct isLevelwiseConforming< SGrid<dim,dimw> > 01995 { 01996 static const bool v = true; 01997 }; 01998 02002 template<int dim, int dimw> 02003 struct isLeafwiseConforming< SGrid<dim,dimw> > 02004 { 02005 static const bool v = true; 02006 }; 02007 02008 } // end namespace Capabilities 02009 02010 } // end namespace Dune 02011 02012 #include"sgrid/sgrid.cc" 02013 02014 #endif