dune-grid  2.1.1
sgrid.hh
Go to the documentation of this file.
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/finitestack.hh>
00011 #include <dune/grid/common/capabilities.hh>
00012 #include <dune/common/bigunsignedint.hh>
00013 #include <dune/common/collectivecommunication.hh>
00014 #include <dune/grid/common/grid.hh>
00015 #include <dune/grid/sgrid/numbering.hh>
00016 #include <dune/grid/common/indexidset.hh>
00017 
00018 #include <dune/grid/genericgeometry/topologytypes.hh>
00019 
00025 namespace Dune {
00026 
00027 //************************************************************************
00031   typedef double sgrid_ctype; 
00032 
00033   // globally define the persistent index type
00034   const int sgrid_dim_bits = 24;   // bits for encoding each dimension
00035   const int sgrid_level_bits = 6;  // bits for encoding level number
00036   const int sgrid_codim_bits = 4;  // bits for encoding codimension
00037 
00038 //************************************************************************
00039 // forward declaration of templates
00040 
00041 template<int dim, int dimworld, class GridImp> class SGeometry;
00042 template<int codim, int dim, class GridImp> class SEntity;
00043 template<int codim, class GridImp> class SEntityPointer;
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 //************************************************************************
00081 template<int mydim, int cdim, class GridImp> 
00082 class SGeometry
00083 : public GeometryDefaultImplementation<mydim,cdim,GridImp,SGeometry>
00084 {
00085 public:
00087         typedef typename GridImp::ctype ctype;
00088 
00090         GeometryType type () const
00091         {
00092           static const GeometryType cubeType(GeometryType::cube,mydim);
00093           return cubeType;
00094         }
00095 
00097         bool affine() const { return true ; }
00098 
00100         int corners () const
00101         {
00102           return 1<<mydim;
00103         }
00104 
00106         FieldVector< ctype, cdim > corner ( const int i ) const
00107         {
00108           return c[i];
00109         }
00110 
00112         FieldVector<ctype, cdim > center ( ) const
00113         {
00114           return centroid;
00115         }
00116 
00118         FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const;
00119 
00121         FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& global) const;
00122 
00142         ctype integrationElement (const FieldVector<ctype, mydim>& local) const
00143         {
00144           return volume();
00145         }
00146 
00148         ctype volume() const;
00149         
00150         const FieldMatrix<ctype, mydim, cdim > &jacobianTransposed ( const FieldVector< ctype, mydim > &local ) const;
00151         const FieldMatrix<ctype,cdim,mydim>& jacobianInverseTransposed (const FieldVector<ctype, mydim>& local) const;
00152 
00154         void print (std::ostream& ss, int indent) const;
00155 
00160         void make (FieldMatrix<ctype,mydim+1,cdim>& __As);
00161 
00163         SGeometry () : builtinverse(false) {};
00164 
00165 private:
00166         FieldVector<ctype, cdim> s;          
00167         FieldVector<ctype, cdim> centroid;   
00168         FieldMatrix<ctype,mydim,cdim> A;     
00169         array<FieldVector<ctype, cdim>, 1<<mydim> c; 
00170         mutable FieldMatrix<ctype,cdim,mydim> Jinv;       
00171         mutable bool builtinverse;
00172 };
00173 
00175 template<int cdim, class GridImp> 
00176 class SGeometry<0,cdim,GridImp>
00177 : public GeometryDefaultImplementation<0,cdim,GridImp,SGeometry>
00178 {
00179 public:
00181         typedef typename GridImp::ctype ctype;
00182 
00184         GeometryType type () const
00185         {
00186           static const GeometryType cubeType(GeometryType::cube,0);
00187           return cubeType;
00188         }
00189 
00191         bool affine() const { return true ; }
00192 
00194         int corners () const
00195         {
00196           return 1;
00197         }
00198 
00200         FieldVector<ctype, cdim > corner ( const int i ) const
00201         {
00202           return s;
00203         }
00204 
00206         FieldVector<ctype, cdim > center ( ) const
00207         {
00208           return s;
00209         }
00210 
00212         void print (std::ostream& ss, int indent) const;
00213 
00215         void make (FieldMatrix<ctype,1,cdim>& __As);
00216 
00218         FieldVector<ctype, cdim> global (const FieldVector<ctype, 0>& local) const { return corner(0); }
00219 
00221         FieldVector<ctype, 0> local (const FieldVector<ctype, cdim>& global) const { return FieldVector<ctype,0> (0.0); } 
00222 
00244 
00245         SGeometry () {};
00246 
00252         ctype volume() const 
00253         {
00254           return 1;
00255         }
00256 
00262         ctype integrationElement(const FieldVector<ctype, 0>& local) const {
00263           return volume();
00264         }
00265         
00266         const FieldMatrix<ctype, 0, cdim > &jacobianTransposed ( const FieldVector<ctype, 0 > &local ) const
00267         {
00268           static const FieldMatrix<ctype, 0, cdim > dummy ( ctype( 0 ) );
00269           return dummy;
00270         }
00271 
00272         const FieldMatrix<ctype,cdim,0>& jacobianInverseTransposed (const FieldVector<ctype, 0>& local) const
00273         {
00274           static const FieldMatrix<ctype,cdim,0> dummy( ctype(0) );
00275           return dummy;
00276         }
00277 
00278 protected:
00279         FieldVector<ctype, cdim> s;         
00280 };
00281 
00282 template <int mydim, int cdim, class GridImp>
00283 inline std::ostream& operator<< (std::ostream& s, SGeometry<mydim,cdim,GridImp>& e)
00284 {
00285         e.print(s,0);
00286         return s;
00287 }
00288 
00289 //************************************************************************
00294 template<int codim, int dim, class GridImp, template<int,int,class> class EntityImp>
00295 class SEntityBase :
00296   public EntityDefaultImplementation<codim,dim,GridImp,EntityImp>
00297 {
00298   friend class SEntityPointer<codim,GridImp>;
00299   friend class SIntersectionIterator<GridImp>;
00300   enum { dimworld = GridImp::dimensionworld };
00301 public:
00302   typedef typename GridImp::ctype ctype;
00303   typedef typename GridImp::template Codim<codim>::Geometry Geometry;
00304   typedef typename GridImp::PersistentIndexType PersistentIndexType;
00305 
00307   int level () const
00308   {
00309     return l;
00310   }
00311 
00313   int globalIndex() const;
00314   
00316   GeometryType type () const
00317   {
00318     static const GeometryType cubeType(GeometryType::cube,dim-codim);
00319     return cubeType;
00320   }
00321 
00323   const Geometry& geometry () const
00324   {
00325       if (!builtgeometry) makegeometry();
00326       
00327       // return result
00328       return geo;
00329   }
00330   
00331   PartitionType partitionType () const { return InteriorEntity; }
00332 
00334   SEntityBase (GridImp* _grid, int _l, int _index) :
00335     grid(_grid),
00336     l(_l),
00337     index(_index),
00338     z(grid->z(l,index,codim)),
00339     geo(SGeometry<dim-codim,dimworld,GridImp>()),
00340     builtgeometry(false) {}
00341 
00343   SEntityBase () :
00344     geo(SGeometry<dim-codim,dimworld,GridImp>()),
00345     builtgeometry(false) // mark geometry as not built
00346   {}
00347 
00349   SEntityBase ( const SEntityBase& other ) : 
00350     grid(other.grid), 
00351     l(other.l),
00352     index(other.index), 
00353     z(other.z), 
00354     geo(SGeometry<dim-codim,dimworld,GridImp>()), // do not copy geometry 
00355     builtgeometry(false) // mark geometry as not built
00356   {}
00357     
00359   void make (GridImp* _grid, int _l, int _id);
00360 
00362   void make (int _l, int _id);
00363 
00365   void makegeometry () const;
00366   
00368   PersistentIndexType persistentIndex () const 
00369   {
00370     return grid->persistentIndex(l, codim, z);
00371   }
00372 
00374   int compressedIndex () const
00375   {
00376     return index;
00377   }
00378 
00380   int compressedLeafIndex () const
00381   {
00382     // codim != dim -> there are no copies of entities
00383     // maxlevel -> ids are fine
00384     if (codim<dim || l==grid->maxLevel())
00385       return compressedIndex();
00386     
00387     // this is a vertex which is not on the finest level
00388     // move coordinates up to maxlevel (multiply by 2 for each level
00389     array<int,dim> coord;
00390     for (int k=0; k<dim; k++)
00391       coord[k] = z[k]*(1<<(grid->maxLevel()-l));
00392     
00393     // compute number with respect to maxLevel
00394     return grid->n(grid->maxLevel(),coord);
00395   }
00396 
00397 protected:
00398   // this is how we implement our elements
00399   GridImp* grid;         
00400   int l;                 
00401   int index;             
00402   array<int,dim> z;      
00403   mutable Geometry geo;  
00404   mutable bool builtgeometry;   
00405 };
00406 
00407 
00413 template<int codim, int dim, class GridImp> 
00414 class SEntity : public SEntityBase<codim,dim,GridImp,SEntity>
00415 {
00416   typedef Dune::SEntityBase<codim,dim,GridImp,Dune::SEntity> SEntityBase;
00417 public:
00419   SEntity (GridImp* _grid, int _l, int _id) :
00420     SEntityBase(_grid,_l,_id) {};
00421 };
00422 
00449 template<int dim, class GridImp>
00450 class SEntity<0,dim,GridImp> : public SEntityBase<0,dim,GridImp,SEntity>
00451 {
00452   enum { dimworld = GridImp::dimensionworld };
00453   typedef Dune::SEntityBase<0,dim,GridImp,Dune::SEntity> SEntityBase;
00454   using SEntityBase::grid;
00455   using SEntityBase::l;
00456   using SEntityBase::index;
00457   using SEntityBase::z;
00458 public:
00459   typedef typename GridImp::ctype ctype;
00460   typedef typename GridImp::template Codim<0>::Geometry Geometry;
00461   typedef typename GridImp::template Codim<0>::LocalGeometry LocalGeometry;
00462   template <int cd>
00463   struct Codim
00464   {
00465     typedef typename GridImp::template Codim<cd>::EntityPointer EntityPointer;
00466   };
00467   typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer;
00468   typedef typename GridImp::LeafIntersectionIterator IntersectionIterator;
00469   typedef typename GridImp::HierarchicIterator HierarchicIterator;
00470   typedef typename GridImp::PersistentIndexType PersistentIndexType;
00471 
00473   friend class SHierarchicIterator<GridImp>;
00474   
00479   template<int cc> int count () const; 
00480 
00485   template<int cc> typename Codim<cc>::EntityPointer subEntity (int i) const;
00486 
00488   int subCompressedIndex (int codim, int i) const
00489   {
00490         if (codim==0) return this->compressedIndex();
00491     // compute subIndex
00492     return (this->grid)->n(this->l, this->grid->subz(this->z,i,codim));
00493   }
00494 
00498   int subCompressedLeafIndex (int codim, int i) const
00499   {
00500         if (codim==0) return this->compressedLeafIndex();
00501 
00502     assert(this->l == this->grid->maxLevel());
00503     // compute subIndex
00504     return (this->grid)->n(this->l, this->grid->subz(this->z,i,codim));
00505   }
00506 
00508   PersistentIndexType subPersistentIndex (int codim, int i) const
00509   {
00510         if (codim==0) return this->persistentIndex();
00511     // compute subId
00512     return this->grid->persistentIndex(this->l, codim, this->grid->subz(this->z,i,codim));
00513   }
00514   
00522   IntersectionIterator ibegin () const;
00523   IntersectionIterator ileafbegin () const;
00524   IntersectionIterator ilevelbegin () const;
00526   IntersectionIterator iend () const;
00527   IntersectionIterator ileafend () const;
00528   IntersectionIterator ilevelend () const;
00529 
00535   EntityPointer father () const;
00536 
00538   bool hasFather () const
00539   {
00540     return (this->level()>0);
00541   }
00542 
00544   bool isLeaf () const
00545     {
00546       return ( this->grid->maxLevel() == this->level() );
00547     }
00548 
00560   const LocalGeometry& geometryInFather () const;
00561 
00568   HierarchicIterator hbegin (int maxLevel) const;
00569 
00571   HierarchicIterator hend (int maxLevel) const;
00572 
00573   // members specific to SEntity
00575   SEntity (GridImp* _grid, int _l, int _index) : 
00576     SEntityBase(_grid,_l,_index),
00577     built_father(false),
00578     in_father_local(SGeometry<dim,dim,GridImp>())
00579     {}
00580 
00581   SEntity (const SEntity& other ) :
00582     SEntityBase(other.grid, other.l, other.index ),
00583     built_father(false),
00584     in_father_local(SGeometry<dim,dim,GridImp>())
00585     {}
00586 
00588   void make (GridImp* _grid, int _l, int _id)
00589     {
00590       SEntityBase::make(_grid,_l,_id);
00591       built_father = false;
00592     }
00593 
00595   void make (int _l, int _id)
00596     {
00597       SEntityBase::make(_l,_id);
00598       built_father = false;
00599     }
00600   
00601 private:
00602 
00603   SEntity();
00604   
00605   mutable bool built_father;
00606   mutable int father_index;
00607   mutable LocalGeometry in_father_local;
00608   void make_father() const;
00609 };
00610 
00611 
00612 //************************************************************************
00620 struct SHierarchicStackElem {
00621   int l;
00622   int index;
00623   SHierarchicStackElem () : l(-1), index(-1) {}
00624   SHierarchicStackElem (int _l, int _index) {l=_l; index=_index;}
00625   bool operator== (const SHierarchicStackElem& s) const {return !operator!=(s);}
00626   bool operator!= (const SHierarchicStackElem& s) const {return l!=s.l || index!=s.index;}
00627 };
00628 
00629 template<class GridImp>
00630 class SHierarchicIterator :
00631   public Dune::SEntityPointer <0,GridImp>
00632 {
00633   friend class SHierarchicIterator<const GridImp>;
00634   enum { dim = GridImp::dimension };
00635   enum { dimworld = GridImp::dimensionworld };
00636   typedef Dune::SEntityPointer<0,GridImp> SEntityPointer;
00637   using SEntityPointer::realEntity;
00638   using SEntityPointer::grid;
00639   using SEntityPointer::l;
00640   using SEntityPointer::index;
00641 public:
00642   typedef typename GridImp::template Codim<0>::Entity Entity;
00643   typedef typename GridImp::ctype ctype;
00644 
00646   void increment();
00647 
00654   SHierarchicIterator (GridImp* _grid,
00655                        const Dune::SEntity<0,GridImp::dimension,GridImp>& _e,
00656                        int _maxLevel, bool makeend) :
00657     SEntityPointer(_grid,_e.level(),_e.compressedIndex())
00658     {
00659       // without sons, we are done
00660       // (the end iterator is equal to the calling iterator)
00661       if (makeend) return;
00662       
00663       // remember element where begin has been called
00664       orig_l = this->entity().level();
00665       orig_index = _grid->getRealImplementation(this->entity()).compressedIndex();
00666       
00667       // push original element on stack
00668       SHierarchicStackElem originalElement(orig_l, orig_index);
00669       stack.push(originalElement);
00670       
00671       // compute maxLevel
00672       maxLevel = std::min(_maxLevel,this->grid->maxLevel());
00673       
00674       // ok, push all the sons as well
00675       push_sons(orig_l,orig_index);
00676       
00677       // and pop the first son
00678       increment();
00679     }
00680 
00681 private:
00682   int maxLevel;                
00683   int orig_l, orig_index;         
00684 
00685   FiniteStack<SHierarchicStackElem,GridImp::MAXL> stack;      
00686   
00687   void push_sons (int level, int fatherid); 
00688 };
00689 
00690 //************************************************************************
00697 template<class GridImp>
00698 class SIntersectionIterator
00699 {
00700   enum { dim=GridImp::dimension };
00701   enum { dimworld=GridImp::dimensionworld };
00702 public:
00703   typedef typename GridImp::template Codim<0>::Entity Entity;
00704   typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer;
00705   typedef typename GridImp::template Codim<1>::Geometry Geometry;
00706   typedef typename GridImp::template Codim<1>::LocalGeometry LocalGeometry;
00707   typedef Dune::SIntersection<GridImp> IntersectionImp;
00708   typedef Dune::Intersection<const GridImp, Dune::SIntersection> Intersection;
00710   enum { dimension=dim };
00712   enum { dimensionworld=dimworld };
00714   typedef typename GridImp::ctype ctype;
00715 
00717   bool equals(const SIntersectionIterator<GridImp>& i) const;
00719   void increment();
00720 
00722   const Intersection & dereference() const
00723   {
00724     return intersection;
00725   }
00726 
00729   EntityPointer inside() const;
00730 
00733   EntityPointer outside() const;
00734   
00736   bool boundary () const;
00737 
00739   bool conforming () const;
00740 
00741   int boundaryId () const {
00742     if (boundary()) return count + 1;
00743     return 0;
00744   };
00745 
00746   int boundarySegmentIndex () const {
00747     if (boundary())
00748       return grid->boundarySegmentIndex(self.level(), count, zred);
00749     return -1;
00750   };
00751   
00753   bool neighbor () const;
00754 
00756   FieldVector<ctype, GridImp::dimensionworld> outerNormal (const FieldVector<ctype, GridImp::dimension-1>& local) const
00757     {
00758       return unitOuterNormal(local);
00759     }
00761   FieldVector<ctype, GridImp::dimensionworld> unitOuterNormal (const FieldVector<ctype, GridImp::dimension-1>& local) const
00762     {
00763       return centerUnitOuterNormal();
00764     }
00766   FieldVector<ctype, GridImp::dimensionworld> centerUnitOuterNormal () const
00767     {
00768       // while we are at it, compute normal direction
00769       FieldVector<ctype, dimworld> normal(0.0);
00770       if (count%2)
00771         normal[count/2] =  1.0; // odd
00772       else
00773         normal[count/2] = -1.0; // even
00774       
00775       return normal; 
00776     }
00778   FieldVector<ctype, GridImp::dimensionworld> integrationOuterNormal (const FieldVector<ctype, GridImp::dimension-1>& local) const
00779     {
00780         FieldVector<ctype, dimworld> n = unitOuterNormal(local);
00781         n *= geometry().integrationElement(local);
00782         return n;
00783     }
00784 
00788   const LocalGeometry &geometryInInside () const;
00792   const LocalGeometry &geometryInOutside () const;
00796   const Geometry &geometry () const;
00797 
00799   GeometryType type () const
00800   {
00801     static const GeometryType cubeType(GeometryType::cube,dim-1);
00802     return cubeType;
00803   }
00804 
00806   int indexInInside () const;
00808   int indexInOutside () const;
00809 
00811   SIntersectionIterator (GridImp* _grid, const SEntity<0,dim,GridImp >* _self, int _count) :
00812     self(*_self), ne(self), grid(_grid),
00813     partition(_grid->partition(grid->getRealImplementation(ne).l,_self->z)),
00814     zred(_grid->compress(grid->getRealImplementation(ne).l,_self->z)),
00815     is_self_local(SGeometry<dim-1, dim, GridImp>()),
00816     is_global(SGeometry<dim-1, dimworld, GridImp>()),
00817     is_nb_local(SGeometry<dim-1, dim, GridImp>()),
00818     intersection(IntersectionImp(*this))
00819   {
00820     // make neighbor
00821     make(_count);
00822   }
00823 
00824   SIntersectionIterator (const SIntersectionIterator & other) :
00825     self(other.self), ne(other.ne), grid(other.grid),
00826     partition(other.partition), zred(other.zred),
00827     count(other.count), valid_count(other.valid_count),
00828     valid_nb(other.valid_nb), is_on_boundary(other.is_on_boundary),
00829     built_intersections(false),
00830     is_self_local(SGeometry<dim-1, dim, GridImp>()),
00831     is_global(SGeometry<dim-1, dimworld, GridImp>()),
00832     is_nb_local(SGeometry<dim-1, dim, GridImp>()),
00833     intersection(IntersectionImp(*this))
00834   {
00835   }
00836   
00838   SIntersectionIterator& operator = (const SIntersectionIterator& it)
00839     {
00840       assert(grid == it.grid);
00841       
00842       /* Assert same Iterator Context */
00843       if (partition != it.partition)
00844         DUNE_THROW(GridError, "assignment of SIntersectionIterator "
00845                    << "with different partition");
00846 
00847       /* Assign current position and update ne */
00848       self = it.self;
00849       make(it.count);
00850 
00851       return *this;
00852     }
00853   
00854 private:
00855   void make (int _count) const;           
00856   void makeintersections () const;        
00857   EntityPointer self;                     
00858   mutable EntityPointer ne;               
00859   const GridImp * grid;                   
00860   const int partition;                    
00861   const array<int,dim> zred;              
00862   mutable int count;                      
00863   mutable bool valid_count;               
00864   mutable bool valid_nb;                  
00865   mutable bool is_on_boundary;            
00866   mutable bool built_intersections;       
00867   mutable LocalGeometry is_self_local;    
00868   mutable Geometry is_global;             
00869   mutable LocalGeometry is_nb_local;      
00870   Intersection intersection;
00871 };
00872 
00873 template<class GridImp>
00874 class SIntersection
00875 {
00876   enum { dim=GridImp::dimension };
00877   enum { dimworld=GridImp::dimensionworld };
00878 public:
00879   typedef typename GridImp::template Codim<0>::Entity Entity;
00880   typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer;
00881   typedef typename GridImp::template Codim<1>::Geometry Geometry;
00882   typedef typename Geometry::LocalCoordinate LocalCoordinate;
00883   typedef typename Geometry::GlobalCoordinate GlobalCoordinate;
00884   typedef typename GridImp::template Codim<1>::LocalGeometry LocalGeometry;
00885   typedef Dune::Intersection<const GridImp, Dune::SIntersectionIterator> Intersection;
00887   enum { dimension=dim };
00889   enum { dimensionworld=dimworld };
00891   typedef typename GridImp::ctype ctype;
00892 
00893   bool boundary () const
00894   {
00895     return is.boundary();
00896   }
00897 
00899   int boundaryId () const
00900   {
00901     return is.boundaryId();
00902   }
00903 
00905   size_t boundarySegmentIndex () const
00906   {
00907     return is.boundarySegmentIndex();
00908   }
00909 
00911   bool neighbor () const 
00912   {
00913     return is.neighbor();
00914   }
00915 
00917   EntityPointer inside() const
00918   {
00919     return is.inside();
00920   }
00921 
00923   EntityPointer outside() const
00924   {
00925     return is.outside();
00926   }
00927   
00929   bool conforming () const 
00930   {
00931     return is.conforming();
00932   }
00933 
00935   const LocalGeometry &geometryInInside () const
00936   {
00937     return is.geometryInInside();
00938   }
00939 
00941   const LocalGeometry &geometryInOutside () const
00942   {
00943     return is.geometryInOutside();
00944   }
00945 
00947   const Geometry &geometry () const
00948   {
00949     return is.geometry();
00950   }
00951 
00953   GeometryType type () const
00954   {
00955     return is.type();
00956   }
00957 
00959   int indexInInside () const
00960   {
00961     return is.indexInInside();
00962   }
00963 
00965   int indexInOutside () const
00966   {
00967     return is.indexInOutside();
00968   }
00969 
00971   GlobalCoordinate outerNormal (const LocalCoordinate& local) const
00972   {
00973     return is.outerNormal(local);
00974   }
00975 
00977   GlobalCoordinate integrationOuterNormal (const LocalCoordinate& local) const
00978   {
00979     return is.integrationOuterNormal(local);
00980   }
00981 
00983   GlobalCoordinate unitOuterNormal (const LocalCoordinate& local) const
00984   {
00985     return is.unitOuterNormal(local);
00986   }
00987 
00989   GlobalCoordinate centerUnitOuterNormal () const
00990   {
00991     return is.centerUnitOuterNormal();
00992   }
00993 
00995   SIntersection (const SIntersectionIterator<GridImp> & is_) : is(is_) {}
00996   
00997 private:
00998 #ifndef DOXYGEN // doxygen can't handle this recursive usage
00999   const SIntersectionIterator<GridImp> & is;
01000 #endif
01001 };
01002 
01003 //************************************************************************
01004 
01008 template <class T>
01009 class AutoPtrStack : public std::stack<T*>
01010 {
01011 public:
01012     ~AutoPtrStack()
01013     {
01014         while(! this->empty())
01015         {
01016             T* e = this->top();
01017             delete e;
01018             this->pop();
01019         }
01020     }
01021 };
01022 
01025 template<int codim, class GridImp>
01026 class SEntityPointer
01027 {
01028   enum { dim = GridImp::dimension };
01029   friend class SIntersectionIterator<GridImp>;
01030 public:
01031   typedef SEntityPointer<codim,GridImp> EntityPointerImp;
01032   typedef typename GridImp::template Codim<codim>::Entity Entity;
01034   enum { codimension = codim };
01035   
01037   bool equals(const SEntityPointer<codim,GridImp>& i) const;
01039   Entity& dereference() const;
01041   int level () const;
01042 
01044   SEntityPointer (GridImp * _grid, int _l, int _index) :
01045     grid(_grid), l(_l), index(_index), 
01046     e(0) 
01047   {
01048   }
01049 
01051   SEntityPointer (const SEntity<codim,dim,GridImp> & _e) :
01052     grid(_e.grid), l(_e.l), index(_e.index),
01053     e(0)
01054   {
01055   }
01056 
01058   SEntityPointer (const SEntityPointer<codim,GridImp>& other) :
01059     grid(other.grid), l(other.l), index(other.index),
01060     e( 0 )
01061   {
01062   }
01063 
01065   ~SEntityPointer()
01066   {
01067     compactify();
01068 #ifndef NDEBUG 
01069     index = -1;
01070 #endif
01071   }
01072 
01074   SEntityPointer& operator = (const SEntityPointer& other)
01075   {
01076     grid = other.grid;
01077     l = other.l;
01078     index = other.index;
01079 
01080     compactify();
01081     return *this;
01082   }
01083 
01085   inline void compactify()
01086   {
01087     if( e )
01088     {
01089       enStack().push( e );
01090       e = 0;
01091     }
01092   }
01093   
01094 protected:
01095   inline SEntity<codim,dim,GridImp>& realEntity() const
01096   {
01097     return grid->getRealImplementation(entity());
01098   }
01099   
01100   inline Entity& entity() const
01101   {
01102     if( ! e )
01103     {
01104       e = getEntity( grid, l, index );
01105     }
01106     return *e;
01107   }
01108 
01109   typedef AutoPtrStack< Entity > EntityStackType;
01110   static inline EntityStackType& enStack()
01111   {
01112     static EntityStackType eStack;
01113     return eStack;
01114   }
01115 
01116   inline Entity* getEntity(GridImp* _grid, int _l, int _id ) const
01117   {
01118     // get stack reference 
01119     EntityStackType& enSt = enStack();
01120 
01121     if( enSt.empty() )
01122     {
01123       return (new Entity(SEntity<codim,dim,GridImp>(_grid, _l, _id)));
01124     }
01125     else
01126     {
01127       Entity* e = enSt.top();
01128       enSt.pop();
01129       grid->getRealImplementation(*e).make(_grid, _l,_id);
01130       return e;
01131     }
01132   }
01133 
01134   GridImp* grid;                 
01135   int l;                         
01136   mutable int index;             
01137   mutable Entity* e;             
01138 };
01139 //************************************************************************
01140 
01141 
01144 template<int codim, PartitionIteratorType pitype, class GridImp>
01145 class SLevelIterator :
01146   public Dune::SEntityPointer <codim,GridImp>
01147 {
01148   friend class SLevelIterator<codim, pitype,const GridImp>;
01149   enum { dim = GridImp::dimension };
01150   typedef Dune::SEntityPointer<codim,GridImp> SEntityPointer;
01151   using SEntityPointer::realEntity;
01152   using SEntityPointer::l;
01153   using SEntityPointer::index;
01154 public:
01155   typedef typename GridImp::template Codim<codim>::Entity Entity;
01156 
01158   void increment();
01159 
01161   SLevelIterator (GridImp * _grid, int _l, int _id) :
01162     SEntityPointer(_grid,_l,_id) {}
01163 };
01164 
01165 
01166 //========================================================================
01171 //========================================================================
01172 
01173 template<class GridImp>
01174 class SGridLevelIndexSet : public IndexSet<GridImp,SGridLevelIndexSet<GridImp> >
01175 {
01176   typedef SGridLevelIndexSet< GridImp > This;
01177   typedef IndexSet< GridImp, This > Base;
01178 
01179   enum { dim = GridImp::dimension };
01180 
01181 public:
01182 
01184   SGridLevelIndexSet ( const GridImp &g, int l )
01185   : grid( g ),
01186     level( l )
01187   {
01188     // TODO move list of geometrytypes to grid, can be computed static (singleton)
01189     // contains a single element type;
01190         for (int codim=0; codim<=GridImp::dimension; codim++)
01191           mytypes[codim].push_back(GeometryType(GeometryType::cube,GridImp::dimension-codim));
01192   }
01193 
01195   template<int cd>
01196   int index (const typename GridImp::Traits::template Codim<cd>::Entity& e) const 
01197   {
01198       return grid.getRealImplementation(e).compressedIndex(); 
01199   }
01200 
01201   template< int cc >
01202   int subIndex ( const typename GridImp::Traits::template Codim< cc >::Entity &e,
01203                  int i, unsigned int codim ) const
01204   {
01205       if( cc == 0 )
01206           return grid.getRealImplementation(e).subCompressedIndex(codim, i);
01207       else
01208           DUNE_THROW( NotImplemented, "subIndex for higher codimension entity not implemented for SGrid." );
01209   }
01210 
01211   // return true if the given entity is contained in \f$E\f$.
01212   template< class EntityType >
01213   bool contains ( const EntityType &e ) const
01214   {
01215     return (e.level() == level);
01216   }
01217 
01219   int size (GeometryType type) const
01220   {
01221     return grid.size( level, type );
01222   }
01223   
01225   int size (int codim) const
01226   {
01227     return grid.size( level, codim );
01228   }
01229 
01231   const std::vector<GeometryType>& geomTypes (int codim) const
01232   {
01233         return mytypes[codim];
01234   }
01235 
01236 private:
01237   const GridImp& grid;
01238   int level;
01239   std::vector<GeometryType> mytypes[GridImp::dimension+1];
01240 };
01241 
01242 // Leaf Index Set
01243 
01244 template<class GridImp>
01245 class SGridLeafIndexSet : public IndexSet<GridImp,SGridLeafIndexSet<GridImp> >
01246 {
01247   typedef SGridLeafIndexSet< GridImp > This;
01248   typedef IndexSet< GridImp, This > Base;
01249 
01250   enum { dim = GridImp::dimension };
01251 
01252 public:
01253 
01255   explicit SGridLeafIndexSet ( const GridImp &g )
01256   : grid( g )
01257   {
01258     // TODO move list of geometrytypes to grid, can be computed static (singleton)
01259     // contains a single element type;
01260         for (int codim=0; codim<=dim; codim++)
01261           mytypes[codim].push_back(GeometryType(GeometryType::cube,dim-codim));
01262   }
01263 
01265   /*
01266     We use the remove_const to extract the Type from the mutable class,
01267     because the const class is not instantiated yet.
01268   */
01269   template<int cd>
01270   int index (const typename remove_const<GridImp>::type::Traits::template Codim<cd>::Entity& e) const 
01271   {
01272         return grid.getRealImplementation(e).compressedLeafIndex(); 
01273   }
01274 
01275   template< int cc >
01276   int subIndex ( const typename GridImp::Traits::template Codim< cc >::Entity &e,
01277                  int i, unsigned int codim ) const
01278   {
01279       if( cc == 0 )
01280           return grid.getRealImplementation(e).subCompressedIndex(codim, i);
01281       else
01282           DUNE_THROW( NotImplemented, "subIndex for higher codimension entity not implemented for SGrid." );
01283   }
01284 
01286   int size (GeometryType type) const
01287   {
01288     return grid.size( grid.maxLevel(), type );
01289   }
01290 
01292   int size (int codim) const
01293   {
01294     return grid.size( grid.maxLevel(), codim );
01295   }
01296 
01297   // return true if the given entity is contained in \f$E\f$.
01298   template< class EntityType >
01299   bool contains ( const EntityType &e ) const
01300   {
01301     return (e.level() == grid.maxLevel());
01302   }
01303 
01305   const std::vector<GeometryType>& geomTypes (int codim) const
01306   {
01307         return mytypes[codim];
01308   }
01309 
01310 private:
01311   const GridImp& grid;
01312   std::vector<GeometryType> mytypes[dim+1];
01313 };
01314 
01315 
01316 //========================================================================
01321 //========================================================================
01322 
01323 template<class GridImp>
01324 class SGridGlobalIdSet :
01325   public IdSet<GridImp,SGridGlobalIdSet<GridImp>, typename remove_const<GridImp>::type::PersistentIndexType>
01326   /*
01327     We used the remove_const to extract the Type from the mutable class,
01328     because the const class is not instantiated yet.
01329   */
01330 {
01331   typedef SGridGlobalIdSet< GridImp > This;
01332 
01333 public:
01334   
01337   using IdSet<GridImp,SGridGlobalIdSet<GridImp>, typename remove_const<GridImp>::type::PersistentIndexType>::subId;
01338 
01340   /*
01341     We use the remove_const to extract the Type from the mutable class,
01342     because the const class is not instantiated yet.
01343   */
01344   typedef typename remove_const<GridImp>::type::PersistentIndexType IdType;
01345 
01347   explicit SGridGlobalIdSet ( const GridImp &g )
01348   : grid( g )
01349   {}
01350 
01352   /*
01353     We use the remove_const to extract the Type from the mutable class,
01354     because the const class is not instantiated yet.
01355   */
01356   template<int cd>
01357   IdType id (const typename remove_const<GridImp>::type::Traits::template Codim<cd>::Entity& e) const 
01358   {
01359           return grid.getRealImplementation(e).persistentIndex();
01360   }
01361 
01363   /*
01364     We use the remove_const to extract the Type from the mutable class,
01365     because the const class is not instantiated yet.
01366   */
01367   IdType subId ( const typename remove_const< GridImp >::type::Traits::template Codim< 0 >::Entity &e,
01368                  int i, unsigned int codim ) const
01369   {
01370       return grid.getRealImplementation(e).subPersistentIndex(codim, i);
01371   }
01372 
01373 private:
01374   const GridImp& grid;
01375 };
01376 
01377 
01378 template<int dim, int dimworld, class ctype>
01379 struct SGridFamily
01380 {
01381   typedef GridTraits<dim,dimworld,Dune::SGrid<dim,dimworld,ctype>,
01382                      SGeometry,SEntity,
01383                      SEntityPointer,SLevelIterator,
01384                      SIntersection, // leaf intersection
01385                      SIntersection, // level intersection
01386                      SIntersectionIterator, // leaf  intersection iter 
01387                      SIntersectionIterator, // level intersection iter
01388                      SHierarchicIterator,
01389                      SLevelIterator,
01390                                          SGridLevelIndexSet<const SGrid<dim,dimworld,ctype> >,
01391                                          SGridLeafIndexSet<const SGrid<dim,dimworld,ctype> >,
01392                                          SGridGlobalIdSet<const SGrid<dim,dimworld,ctype> >,
01393                                          bigunsignedint<dim*sgrid_dim_bits+sgrid_level_bits+sgrid_codim_bits>,
01394                                          SGridGlobalIdSet<const SGrid<dim,dimworld,ctype> >,
01395                                          bigunsignedint<dim*sgrid_dim_bits+sgrid_level_bits+sgrid_codim_bits>, 
01396                                          CollectiveCommunication<Dune::SGrid<dim,dimworld,ctype> > > 
01397   Traits;
01398 };
01399 
01400 
01401 //************************************************************************
01455 template<int dim, int dimworld, typename _ctype = sgrid_ctype>
01456 class SGrid : public GridDefaultImplementation <dim,dimworld,_ctype,SGridFamily<dim,dimworld,_ctype> >
01457 {
01458 public:
01459   typedef SGridFamily<dim,dimworld,_ctype> GridFamily;
01460   typedef bigunsignedint<dim*sgrid_dim_bits+sgrid_level_bits+sgrid_codim_bits> PersistentIndexType;
01461 
01462   // need for friend declarations in entity
01463   typedef SGridLevelIndexSet<SGrid<dim,dimworld> > LevelIndexSetType;
01464   typedef SGridLeafIndexSet<SGrid<dim,dimworld> > LeafIndexSetType;
01465   typedef SGridGlobalIdSet<SGrid<dim,dimworld> > GlobalIdSetType;
01466   
01467   typedef typename SGridFamily<dim,dimworld,_ctype>::Traits Traits;
01468 
01470   enum { MAXL=32 };
01471 
01473   typedef _ctype ctype;
01474 
01475   // constructors
01476 
01484   SGrid (const int * const N_, const ctype * const H_);
01485   
01493   SGrid (const int * const N_, const ctype * const L_, const ctype * const H_);
01494 
01504   SGrid (FieldVector<int,dim> N_, FieldVector<ctype,dim> L_, FieldVector<ctype,dim> H_);
01505 
01507   SGrid ();
01508 
01510   ~SGrid ();
01511 
01514   int maxLevel() const;
01515 
01517   template<int cd, PartitionIteratorType pitype>
01518   typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator lbegin (int level) const;
01519 
01521   template<int cd, PartitionIteratorType pitype>
01522   typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator lend (int level) const;
01523 
01525   template<int cd>
01526   typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator lbegin (int level) const
01527     {
01528       return lbegin<cd,All_Partition>(level);
01529     }
01530 
01532   template<int cd>
01533   typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator lend (int level) const
01534     {
01535       return lend<cd,All_Partition>(level);
01536     }
01537 
01539   template<int cd, PartitionIteratorType pitype>
01540   typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator leafbegin () const;
01541 
01543   template<int cd, PartitionIteratorType pitype>
01544   typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator leafend () const;
01545 
01547   template<int cd>
01548   typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator leafbegin () const
01549     {
01550       return leafbegin<cd,All_Partition>();
01551     }
01552 
01554   template<int cd>
01555   typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator leafend () const
01556     {
01557       return leafend<cd,All_Partition>();
01558     }
01559 
01573   template<class T, template<class> class P, int codim>
01574   void communicate (T& t, InterfaceType iftype, CommunicationDirection dir, int level)
01575   {
01576           // SGrid is sequential and has no periodic boundaries, so do nothing ...
01577           return;
01578   }
01579 
01581   int size (int level, int codim) const;
01582 
01584   int size (int codim) const
01585   {
01586         return size(maxLevel(),codim);
01587   }
01588 
01590   int size (int level, GeometryType type) const
01591   {
01592       return (type.isCube()) ? size(level,dim-type.dim()) : 0;
01593   }
01594 
01596   int size (GeometryType type) const
01597   {
01598         return size(maxLevel(),type);
01599   }
01600 
01602   size_t numBoundarySegments () const
01603   {
01604     return boundarysize;
01605   }
01606 
01608   int global_size (int codim) const;
01609 
01611   int overlapSize (int level, int codim)
01612   {
01613         return 0;
01614   }
01615 
01617   int ghostSize (int level, int codim)
01618   {
01619         return 0;
01620   }
01621 
01622   // these are all members specific to sgrid
01623 
01625   void globalRefine (int refCount);
01626 
01628     const array<int, dim>& dims(int level) const {
01629         return N[level];
01630     }
01631 
01633     const FieldVector<ctype, dimworld>& lowerLeft() const {
01634         return low;
01635     }
01636 
01638     FieldVector<ctype, dimworld> upperRight() const {
01639         return H;
01640     }
01641 
01643   bool adapt ()    
01644   {
01645         globalRefine(1);
01646         return true; 
01647   }
01648   
01649   // The new index sets from DDM 11.07.2005
01650   const typename Traits::GlobalIdSet& globalIdSet() const
01651   {
01652         return *theglobalidset;
01653   }
01654   
01655   const typename Traits::LocalIdSet& localIdSet() const
01656   {
01657         return *theglobalidset;
01658   }
01659 
01660   const typename Traits::LevelIndexSet& levelIndexSet(int level) const
01661   {
01662         assert(level>=0 && level<=maxLevel());
01663         return *(indexsets[level]);
01664   }
01665 
01666   const typename Traits::LeafIndexSet& leafIndexSet() const
01667   {
01668         return *theleafindexset;
01669   }
01670 
01675   template<class DataHandle>
01676   void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int level) const
01677   {
01678   }
01679 
01680   template<class DataHandle>
01681   void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const
01682   {
01683   }
01684 
01685   const CollectiveCommunication<SGrid>& comm () const
01686   {
01687         return ccobj;
01688   }
01689 
01691   int overlapSize (int level, int codim) const
01692   {
01693         return 0;
01694   }
01695 
01697   int overlapSize (int codim) const
01698   {
01699         return 0;
01700   }
01701 
01703   int ghostSize (int level, int codim) const
01704   {
01705         return 0;
01706   }
01707 
01709   int ghostSize (int codim) const
01710   {
01711         return 0;
01712   }
01713 
01714   /*
01715     @}
01716    */
01717 
01718 private:
01719   /*
01720     Make associated classes friends to grant access to the real entity
01721   */
01722   friend class Dune::SGridLevelIndexSet<Dune::SGrid<dim,dimworld> >;
01723   friend class Dune::SGridLeafIndexSet<Dune::SGrid<dim,dimworld> >;
01724   friend class Dune::SGridGlobalIdSet<Dune::SGrid<dim,dimworld> >;
01725   friend class Dune::SIntersectionIterator<Dune::SGrid<dim,dimworld> >;
01726   friend class Dune::SHierarchicIterator<Dune::SGrid<dim,dimworld> >;
01727   friend class Dune::SEntity<0,dim,Dune::SGrid<dim,dimworld> >;
01728   
01729   friend class Dune::SGridLevelIndexSet<const Dune::SGrid<dim,dimworld> >;
01730   friend class Dune::SGridLeafIndexSet<const Dune::SGrid<dim,dimworld> >;
01731   friend class Dune::SGridGlobalIdSet<const Dune::SGrid<dim,dimworld> >;
01732   friend class Dune::SIntersectionIterator<const Dune::SGrid<dim,dimworld> >;
01733   friend class Dune::SHierarchicIterator<const Dune::SGrid<dim,dimworld> >;
01734   friend class Dune::SEntity<0,dim,const Dune::SGrid<dim,dimworld> >;
01735 
01736   template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
01737   friend class Dune::SEntityBase;
01738   
01739   template<int codim_, class GridImp_>
01740   friend class Dune::SEntityPointer;
01741   
01742   template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
01743   friend class Entity;
01744 
01746   FieldVector<ctype, dimworld> pos (int level, array<int,dim>& z) const;
01747  
01749   int calc_codim (int level, const array<int,dim>& z) const;
01750 
01752   int n (int level, const array<int,dim>& z) const;
01753 
01755   array<int,dim> z (int level, int i, int codim) const;
01756         
01758   array<int,dim> subz (const array<int,dim> & z, int i, int codim) const;
01759   
01761   array<int,dim> compress (int level, const array<int,dim>& z) const; 
01762 
01764   array<int,dim> expand (int level, const array<int,dim>& r, int b) const; 
01765 
01769   int partition (int level, const array<int,dim>& z) const; 
01770 
01772   bool exists (int level, const array<int,dim>& zred) const;
01773 
01774   // compute boundary segment index for a given zentity and a face
01775   int boundarySegmentIndex (int l, int face, const array<int,dim> & zentity) const
01776   {
01777     array<int,dim-1> zface;
01778     int dir = face/2;
01779     int side = face%2;
01780     // compute z inside the global face
01781     for (int i=0; i<dir; i++) zface[i] = zentity[i]/(1<<l);
01782     for (int i=dir+1; i<dim; i++) zface[i-1] = zentity[i]/(1<<l);
01783     zface = boundarymapper[dir].expand(zface, 0);
01784     // compute index in the face
01785     int index = boundarymapper[dir].n(zface);
01786     // compute offset
01787     for (int i=0; i<dir; i++)
01788       index += 2*boundarymapper[i].elements(0);
01789     index += side*boundarymapper[dir].elements(0);
01790     return index;
01791   }
01792   
01793   // compute persistent index for a given zentity
01794   PersistentIndexType persistentIndex (int l, int codim, const array<int,dim> & zentity) const
01795   {
01796           if (codim!=dim)
01797                 {
01798                   // encode codim, this would actually not be necessary
01799           // because z is unique in codim
01800                   PersistentIndexType id(codim);
01801 
01802                   // encode level
01803                   id = id << sgrid_level_bits;
01804                   id = id+PersistentIndexType(l);
01805         
01806                   // encode coordinates
01807                   for (int i=dim-1; i>=0; i--)
01808                         {
01809                           id = id << sgrid_dim_bits;
01810                           id = id+PersistentIndexType(zentity[i]);
01811                         }
01812                   
01813                   return id;
01814                 }
01815           else
01816                 {
01817                   // determine min number of trailing zeroes
01818                   // consider that z is on the doubled grid !
01819                   int trailing = 1000;
01820                   for (int i=0; i<dim; i++)
01821                         {
01822                           // count trailing zeros
01823                           int zeros = 0;
01824                           for (int j=0; j<l; j++)
01825                                 if (zentity[i]&(1<<(j+1)))
01826                                   break;
01827                                 else
01828                                   zeros++;
01829                           trailing = std::min(trailing,zeros);
01830                         }
01831 
01832                   // determine the level of this vertex
01833                   int level = l-trailing;
01834 
01835                   // encode codim
01836                   PersistentIndexType id(dim);
01837 
01838                   // encode level
01839                   id = id << sgrid_level_bits;
01840                   id = id+PersistentIndexType(level);
01841         
01842                   // encode coordinates
01843                   for (int i=dim-1; i>=0; i--)
01844                         {
01845                           id = id << sgrid_dim_bits;
01846                           id = id+PersistentIndexType(zentity[i]>>trailing);
01847                         }
01848         
01849                   return id;
01850                 }
01851   }
01852     
01853   // disable copy and assign
01854   SGrid(const SGrid &) {};
01855   SGrid & operator = (const SGrid &) { return *this; };
01856   // generate SGrid 
01857   void makeSGrid (const array<int,dim>& N_, const FieldVector<ctype, dim>& L_, const FieldVector<ctype, dim>& H_);
01858 
01859   /*
01860     internal data
01861    */
01862   CollectiveCommunication<SGrid> ccobj;
01863 
01864   ReservedVector<SGridLevelIndexSet<const SGrid<dim,dimworld> >*, MAXL> indexsets;
01865   SGridLeafIndexSet<const SGrid<dim,dimworld> > *theleafindexset;
01866   SGridGlobalIdSet<const SGrid<dim,dimworld> > *theglobalidset;
01867 
01868   int L;                          // number of levels in hierarchic mesh 0<=level<L
01869   FieldVector<ctype, dim> low;    // lower left corner of the grid
01870   FieldVector<ctype, dim> H;      // length of cube per direction
01871   array<int,dim> *N;              // number of elements per direction
01872   FieldVector<ctype, dim> *h;     // mesh size per direction
01873   mutable CubeMapper<dim> *mapper;// a mapper for each level
01874 
01875   // boundary segement index set
01876   array<CubeMapper<dim-1>, dim> boundarymapper;  // a mapper for each coarse grid face
01877   int boundarysize;
01878 
01879   // faster implementation of subIndex 
01880   mutable array <int,dim> zrefStatic;     // for subIndex of SEntity
01881   mutable array <int,dim> zentityStatic;  // for subIndex of SEntity
01882 };
01883 
01884 namespace Capabilities
01885 {
01886 
01898   template<int dim, int dimw>
01899   struct hasSingleGeometryType< SGrid<dim,dimw> >
01900   {
01901     static const bool v = true;
01902     static const unsigned int topologyId = GenericGeometry :: CubeTopology< dim > :: type :: id ;
01903   };
01904   
01908   template<int dim, int dimw>
01909   struct isCartesian< SGrid<dim,dimw> >
01910   {
01911     static const bool v = true;
01912   };
01913   
01917   template<int dim, int dimw, int cdim>
01918   struct hasEntity< SGrid<dim,dimw>, cdim>
01919   {
01920     static const bool v = true;
01921   };
01922   
01926   template<int dim, int dimw>
01927   struct isLevelwiseConforming< SGrid<dim,dimw> >
01928   {
01929     static const bool v = true;
01930   };
01931 
01935   template<int dim, int dimw>
01936   struct isLeafwiseConforming< SGrid<dim,dimw> >
01937   {
01938     static const bool v = true;
01939   };
01940 
01941 } // end namespace Capabilities
01942 
01943 } // end namespace Dune
01944 
01945 #include"sgrid/sgrid.cc"
01946 
01947 #endif