dune-grid
2.2.0
|
00001 #ifndef DUNE_SIZECACHE_HH 00002 #define DUNE_SIZECACHE_HH 00003 00004 #include <cassert> 00005 #include <vector> 00006 #include <set> 00007 00008 #include <dune/common/forloop.hh> 00009 #include <dune/common/exceptions.hh> 00010 00011 #include <dune/geometry/type.hh> 00012 #include <dune/geometry/referenceelements.hh> 00013 00014 #include <dune/grid/common/gridenums.hh> 00015 #include <dune/grid/common/capabilities.hh> 00016 00023 namespace Dune { 00024 00026 template <class GridImp> 00027 class SizeCache 00028 { 00029 typedef SizeCache<GridImp> ThisType; 00031 enum { dim = GridImp::dimension }; 00032 00034 enum { nCodim = GridImp::dimension+1 }; 00035 00036 // type of grid 00037 typedef GridImp GridType; 00038 00039 // coordinate type 00040 typedef typename GridType :: ctype ctype ; 00041 00042 // stores all sizes of the levels 00043 mutable std::vector< int > levelSizes_[nCodim]; 00044 00045 // stores all sizes of the levels 00046 mutable std::vector< std::vector< int > > levelTypeSizes_[nCodim]; 00047 00048 // stores all sizes of leafs 00049 mutable int leafSizes_[nCodim]; 00050 00051 // stores all sizes of leafs 00052 mutable std::vector< int > leafTypeSizes_[nCodim]; 00053 00054 // the grid 00055 const GridType & grid_; 00056 00057 // count elements of set by iterating the grid 00058 template < int codim, bool gridHasCodim > 00059 struct CountLevelEntitiesBase 00060 { 00061 template < class SzCacheType > 00062 static void apply(const SzCacheType & sc, int level, int cd) 00063 { 00064 if( cd == codim ) 00065 { 00066 sc.template countLevelEntities<All_Partition,codim> (level); 00067 } 00068 } 00069 }; 00070 00071 template < int codim > 00072 struct CountLevelEntitiesBase< codim, false > 00073 { 00074 template < class SzCacheType > 00075 static void apply(const SzCacheType & sc, int level, int cd) 00076 { 00077 if( cd == codim ) 00078 { 00079 sc.template countLevelEntitiesNoCodim<All_Partition,codim> (level); 00080 } 00081 } 00082 }; 00083 00084 template < int codim > 00085 struct CountLevelEntities 00086 : public CountLevelEntitiesBase< codim, Capabilities :: hasEntity< GridType, codim > :: v > 00087 { 00088 }; 00089 00090 // count elements of set by iterating the grid 00091 template < int codim, bool gridHasCodim > 00092 struct CountLeafEntitiesBase 00093 { 00094 template <class SzCacheType> 00095 static void apply(const SzCacheType & sc, int cd) 00096 { 00097 if( cd == codim ) 00098 { 00099 sc.template countLeafEntities<All_Partition,codim> (); 00100 } 00101 } 00102 }; 00103 00104 // count elements of set by iterating the grid 00105 template < int codim > 00106 struct CountLeafEntitiesBase< codim, false > 00107 { 00108 template <class SzCacheType> 00109 static void apply(const SzCacheType & sc, int cd) 00110 { 00111 if( cd == codim ) 00112 { 00113 sc.template countLeafEntitiesNoCodim<All_Partition,codim> (); 00114 } 00115 } 00116 }; 00117 00118 template < int codim > 00119 struct CountLeafEntities 00120 : public CountLeafEntitiesBase< codim, Capabilities :: hasEntity< GridType, codim > :: v > 00121 { 00122 }; 00123 00124 int gtIndex( const GeometryType& type ) const 00125 { 00126 return type.id() >> 1 ; 00127 } 00128 00129 int sizeCodim( const int codim ) const 00130 { 00131 const int mydim = GridType :: dimension - codim; 00132 return ((1 << mydim) + 1) / 2; 00133 } 00134 00135 // private copy constructor 00136 SizeCache (const SizeCache & ); 00137 public: 00139 SizeCache (const GridType & grid) : grid_( grid ) 00140 { 00141 reset(); 00142 } 00143 00145 void reset() 00146 { 00147 for(int codim=0; codim<nCodim; ++codim) 00148 { 00149 leafSizes_[ codim ] = -1; 00150 leafTypeSizes_[ codim ].resize( sizeCodim( codim ), -1 ); 00151 } 00152 00153 const int numMxl = grid_.maxLevel()+1; 00154 for(int codim=0; codim<nCodim; ++codim) 00155 { 00156 std::vector<int> & vec = levelSizes_[codim]; 00157 vec.resize(numMxl); 00158 levelTypeSizes_[codim].resize( numMxl ); 00159 for(int level = 0; level<numMxl; ++level) 00160 { 00161 vec[level] = -1; 00162 levelTypeSizes_[codim][level].resize( sizeCodim( codim ), -1 ); 00163 } 00164 } 00165 } 00166 00167 //******************************************************************** 00168 // level sizes 00169 //******************************************************************** 00171 int size (int level, int codim) const 00172 { 00173 assert( codim >= 0 ); 00174 assert( codim < nCodim ); 00175 assert( level >= 0 ); 00176 if( level >= (int) levelSizes_[codim].size() ) return 0; 00177 00178 if( levelSizes_[codim][level] < 0) 00179 ForLoop< CountLevelEntities, 0, dim > :: apply( *this, level, codim ); 00180 00181 // CountLevelEntities<ThisType,All_Partition,dim>::count(*this,level,codim); 00182 00183 assert( levelSizes_[codim][level] >= 0 ); 00184 return levelSizes_[codim][level]; 00185 } 00186 00188 int size (int level, GeometryType type) const 00189 { 00190 const int codim = GridType ::dimension - type.dim(); 00191 if( levelSizes_[codim][level] < 0) 00192 ForLoop< CountLevelEntities, 0, dim > :: apply( *this, level, codim ); 00193 00194 assert( levelTypeSizes_[codim][level][gtIndex( type )] >= 0 ); 00195 return levelTypeSizes_[codim][level][gtIndex( type )]; 00196 } 00197 00198 //******************************************************************** 00199 // leaf sizes 00200 //******************************************************************** 00202 int size (int codim) const 00203 { 00204 assert( codim >= 0 ); 00205 assert( codim < nCodim ); 00206 if( leafSizes_[codim] < 0 ) 00207 ForLoop< CountLeafEntities, 0, dim > :: apply( *this, codim ); 00208 00209 assert( leafSizes_[codim] >= 0 ); 00210 return leafSizes_[codim]; 00211 }; 00212 00214 int size ( const GeometryType type ) const 00215 { 00216 const int codim = GridType :: dimension - type.dim(); 00217 if( leafSizes_[codim] < 0 ) 00218 ForLoop< CountLeafEntities, 0, dim > :: apply( *this, codim ); 00219 00220 assert( leafTypeSizes_[codim][ gtIndex( type )] >= 0 ); 00221 return leafTypeSizes_[codim][ gtIndex( type )]; 00222 } 00223 00224 private: 00225 template <PartitionIteratorType pitype, int codim> 00226 void countLevelEntities(int level) const 00227 { 00228 typedef typename GridType :: LevelGridView GridView ; 00229 typedef typename GridView :: template Codim< codim > :: template Partition<pitype> :: Iterator Iterator ; 00230 GridView gridView = grid_.levelView( level ); 00231 Iterator it = gridView.template begin<codim,pitype> (); 00232 Iterator end = gridView.template end<codim,pitype> (); 00233 levelSizes_[codim][level] = countElements(it,end, levelTypeSizes_[codim][level]); 00234 } 00235 00236 template <PartitionIteratorType pitype, int codim> 00237 void countLeafEntities() const 00238 { 00239 // count All_Partition entities 00240 typedef typename GridType :: LeafGridView GridView ; 00241 typedef typename GridView :: template Codim< codim > :: template Partition<pitype> :: Iterator Iterator ; 00242 GridView gridView = grid_.leafView(); 00243 Iterator it = gridView.template begin<codim,pitype> (); 00244 Iterator end = gridView.template end<codim,pitype> (); 00245 leafSizes_[codim] = countElements(it,end, leafTypeSizes_[codim] ); 00246 } 00247 00248 // counts entities with given type for given iterator 00249 template <class IteratorType> 00250 int countElements(IteratorType & it, const IteratorType & end, std::vector<int>& typeSizes) const 00251 { 00252 int overall = 0; 00253 const size_t types = typeSizes.size(); 00254 for(size_t i=0; i<types; ++i) typeSizes[i] = 0; 00255 for( ; it != end; ++it ) 00256 { 00257 const GeometryType type = it->type(); 00258 ++typeSizes[ gtIndex( type ) ]; 00259 ++overall; 00260 } 00261 00262 int sumtypes = 0; 00263 for(size_t i=0; i<types; ++i) sumtypes += typeSizes[i]; 00264 00265 assert( overall == sumtypes ); 00266 return overall; 00267 } 00268 00269 template <PartitionIteratorType pitype, int codim> 00270 void countLevelEntitiesNoCodim(int level) const 00271 { 00272 typedef typename GridType :: LevelGridView GridView ; 00273 typedef typename GridView :: template Codim< 0 > :: template Partition<pitype> :: Iterator Iterator ; 00274 GridView gridView = grid_.levelView( level ); 00275 Iterator it = gridView.template begin< 0, pitype> (); 00276 Iterator end = gridView.template end< 0, pitype> (); 00277 levelSizes_[codim][level] = countElementsNoCodim< codim >(it,end, levelTypeSizes_[codim][level]); 00278 } 00279 00280 template <PartitionIteratorType pitype, int codim> 00281 void countLeafEntitiesNoCodim() const 00282 { 00283 // count All_Partition entities 00284 typedef typename GridType :: LeafGridView GridView ; 00285 typedef typename GridView :: template Codim< 0 > :: template Partition<pitype> :: Iterator Iterator ; 00286 GridView gridView = grid_.leafView(); 00287 Iterator it = gridView.template begin< 0, pitype > (); 00288 Iterator end = gridView.template end< 0, pitype > (); 00289 leafSizes_[codim] = countElementsNoCodim< codim >(it,end, leafTypeSizes_[codim] ); 00290 } 00291 00292 // counts entities with given type for given iterator 00293 template < int codim, class IteratorType > 00294 int countElementsNoCodim(IteratorType & it, const IteratorType & end, std::vector<int>& typeSizes) const 00295 { 00296 typedef typename GridType :: LocalIdSet LocalIdSet ; 00297 typedef typename LocalIdSet :: IdType IdType ; 00298 00299 typedef GenericReferenceElement< ctype, dim > ReferenceElementType; 00300 typedef GenericReferenceElements< ctype, dim > ReferenceElementContainerType; 00301 00302 typedef std::set< IdType > CodimIdSetType ; 00303 00304 typedef typename IteratorType :: Entity ElementType ; 00305 00306 // get id set 00307 const LocalIdSet& idSet = grid_.localIdSet(); 00308 00309 const size_t types = typeSizes.size(); 00310 for(size_t i=0; i<types; ++i) typeSizes[ i ] = 0; 00311 00312 std::vector< CodimIdSetType > typeCount( types ); 00313 00314 // count all elements of codimension codim 00315 for( ; it != end; ++it ) 00316 { 00317 // get entity 00318 const ElementType& element = *it ; 00319 // get reference element 00320 const ReferenceElementType& refElem = 00321 ReferenceElementContainerType :: general( element.type() ); 00322 00323 // count all sub entities of codimension codim 00324 const int count = element.template count< codim > (); 00325 for( int i=0; i< count; ++ i ) 00326 { 00327 // get geometry type 00328 const GeometryType geomType = refElem.type( i, codim ); 00329 // get id of sub entity 00330 const IdType id = idSet.subId( element, i, codim ); 00331 // insert id into set 00332 typeCount[ gtIndex( geomType ) ].insert( id ); 00333 } 00334 } 00335 00336 // accumulate numbers 00337 int overall = 0; 00338 for(size_t i=0; i<types; ++i) 00339 { 00340 typeSizes[ i ] = typeCount[ i ].size(); 00341 overall += typeSizes[ i ]; 00342 } 00343 00344 return overall; 00345 } 00346 }; 00347 00349 template <class GridImp> 00350 struct SingleTypeSizeCache : public SizeCache< GridImp > 00351 { 00352 DUNE_DEPRECATED SingleTypeSizeCache( const GridImp& grid, bool isSimplex , bool isCube, bool notWorry = false ) : SizeCache< GridImp >( grid ) {} 00353 }; 00354 00355 } // end namespace Dune 00356 #endif