dune-grid
2.1.1
|
00001 #ifndef DUNE_MACROGRIDPARSER_HH 00002 #define DUNE_MACROGRIDPARSER_HH 00003 00004 #include <iostream> 00005 #include <fstream> 00006 00007 #include <sstream> 00008 #include <string> 00009 #include <cstring> 00010 #include <vector> 00011 #include <memory> 00012 #include <map> 00013 #include <assert.h> 00014 #include <cmath> 00015 00016 //- Dune includes 00017 #include <dune/common/mpihelper.hh> 00018 #include <dune/common/stdstreams.hh> 00019 #include <dune/grid/common/gridenums.hh> 00020 #include <dune/grid/common/datahandleif.hh> 00021 //- local includes 00022 00023 #include "dgfexception.hh" 00024 #include "entitykey.hh" 00025 #include "dgfparserblocks.hh" 00026 00027 namespace Dune 00028 { 00029 00030 class DGFPrintInfo; 00031 00034 class DuneGridFormatParser 00035 { 00036 template< class GridType > 00037 friend class DGFGridFactory; 00038 template< class GridType > 00039 friend class DGFBaseFactory; 00040 public: 00041 typedef enum {Simplex,Cube,General} element_t; 00042 typedef enum {counterclockwise=1,clockwise=-1} orientation_t; 00043 00045 DuneGridFormatParser ( int rank, int size ); 00046 00047 static bool isDuneGridFormat ( std::istream & ); 00048 00053 bool readDuneGrid( std::istream &, int dimG, int dimW ); 00054 00055 bool readDuneGrid( std::istream &input, int dimG = -1 ) DUNE_DEPRECATED 00056 { 00057 return readDuneGrid( input, dimG, -1 ); 00058 } 00059 00061 void writeTetgenPoly ( const std :: string &, std ::string &, std :: string & ); 00062 void writeTetgenPoly ( std::ostream &out, const bool writeSegments = true ); 00063 00064 protected: 00065 // dimension of world and problem: set through the readDuneGrid() method 00066 int dimw, dimgrid; 00067 // vector of vertex coordinates 00068 std::vector < std::vector <double> > vtx; 00069 int nofvtx; 00070 int vtxoffset; 00071 double minVertexDistance; // min. L^1 distance of distinct points 00072 // vector of elements 00073 std :: vector< std :: vector< unsigned int > > elements; 00074 int nofelements; 00075 // vector of boundary segments + identifier 00076 std::vector < std::vector <int> > bound; 00077 int nofbound; 00078 // map to generate and find boundary segments 00079 typedef std :: map< DGFEntityKey< unsigned int >, int > facemap_t; 00080 facemap_t facemap; 00081 // set by generator depending on element type wanted 00082 element_t element; 00083 // set by the readDuneGrid method depending 00084 // on what type the elements were generated 00085 bool simplexgrid; 00086 // true if grid is generated using the intervall Block 00087 bool cube2simplex; 00088 // parameters on elements 00089 int nofvtxparams,nofelparams; 00090 std::vector<std::vector<double> > vtxParams,elParams; 00091 // write information about generation process 00092 DGFPrintInfo* info; 00093 int rank_; 00094 int size_; 00095 00096 void generateBoundaries ( std::istream &, bool ); 00097 00098 // call to tetgen/triangle 00099 void generateSimplexGrid ( std::istream & ); 00100 void readTetgenTriangle ( const std :: string & ); 00101 00102 // helper methods 00103 void removeCopies (); 00104 void setOrientation ( int use1, int use2, 00105 orientation_t orientation=counterclockwise ); 00106 void setRefinement ( int use1, int use2, int is1=-1, int is2=-1 ); 00107 double testTriang ( int snr ); 00108 00109 std :: vector< double > &getElParam ( int i, std::vector< double > &coord ); 00110 std :: vector< double > &getVtxParam ( int i, std::vector< double > &coord ); 00111 00112 static std::string temporaryFileName (); 00113 }; 00114 00115 00116 class MacroGrid : protected DuneGridFormatParser 00117 { 00118 template< class GridType > 00119 friend class DGFGridFactory; 00120 00121 public: 00122 typedef MPIHelper::MPICommunicator MPICommunicatorType; 00123 00124 protected: 00126 MacroGrid(const char* filename, MPICommunicatorType MPICOMM = MPIHelper::getCommunicator()) 00127 : DuneGridFormatParser( rank(MPICOMM), size(MPICOMM) ) 00128 , filename_(filename) 00129 , MPICOMM_(MPICOMM) {} 00130 00132 MacroGrid(MPICommunicatorType MPICOMM = MPIHelper::getCommunicator()) 00133 : DuneGridFormatParser( rank(MPICOMM), size(MPICOMM) ) 00134 , filename_(0) 00135 , MPICOMM_(MPICOMM) {} 00136 00138 template <class GridType> 00139 inline GridType * createGrid () 00140 { 00141 return Impl<GridType>::generate(*this,filename_,MPICOMM_); 00142 } 00143 private: 00144 static int rank( MPICommunicatorType MPICOMM ) 00145 { 00146 int rank = 0; 00147 #if HAVE_MPI 00148 MPI_Comm_rank( MPICOMM, &rank ); 00149 #endif 00150 return rank; 00151 } 00152 static int size( MPICommunicatorType MPICOMM ) 00153 { 00154 int size = 1; 00155 #if HAVE_MPI 00156 MPI_Comm_size( MPICOMM, &size ); 00157 #endif 00158 return size; 00159 } 00171 template< class GridType > 00172 struct Impl 00173 { 00174 static GridType* generate(MacroGrid& mg, 00175 const char* filename, MPICommunicatorType MPICOMM = MPIHelper::getCommunicator() ) 00176 { 00177 // make assertion depend on the template argument but always evaluate to false 00178 dune_static_assert( GridType::dimension<0,"dgf grid factory missing - did you forget to add the corresponding dgf header or dgfgridtype.hh ?"); 00179 } 00180 }; 00181 00182 const char* filename_; 00183 MPICommunicatorType MPICOMM_; 00184 }; 00185 00186 template <class G> 00187 struct DGFGridFactory 00188 { 00189 typedef G Grid; 00190 const static int dimension = Grid::dimension; 00191 typedef MPIHelper::MPICommunicator MPICommunicatorType; 00192 00193 private: 00194 typedef typename Grid::template Codim< 0 >::Entity Element; 00195 typedef typename Grid::template Codim< dimension >::Entity Vertex; 00196 00197 public: 00198 explicit DGFGridFactory ( std::istream &input, 00199 MPICommunicatorType comm = MPIHelper::getCommunicator() ) 00200 : macroGrid_( comm ) 00201 { 00202 DUNE_THROW( DGFException, "DGF factories using old MacroGrid implementation" 00203 "don't support creation from std::istream." ); 00204 } 00205 00206 explicit DGFGridFactory ( const std::string &filename, 00207 MPICommunicatorType comm = MPIHelper::getCommunicator() ) 00208 : macroGrid_( filename.c_str(), comm ) 00209 { 00210 grid_ = macroGrid_.template createGrid< Grid >(); 00211 00212 if( macroGrid_.nofelparams > 0 ) 00213 { 00214 const size_t nofElements = macroGrid_.elements.size(); 00215 for( size_t i = 0; i < nofElements; ++i ) 00216 { 00217 std::vector< double > coord; 00218 00219 DomainType p(0); 00220 const size_t nofCorners = macroGrid_.elements[i].size(); 00221 for (size_t k=0;k<nofCorners;++k) 00222 for (int j=0;j<DomainType::dimension;++j) 00223 p[j]+=macroGrid_.vtx[macroGrid_.elements[i][k]][j]; 00224 p/=double(nofCorners); 00225 00226 elInsertOrder_.insert( std::make_pair( p, i ) ); 00227 } 00228 } 00229 00230 if( macroGrid_.nofvtxparams > 0 ) 00231 { 00232 const size_t nofVertices = macroGrid_.vtx.size(); 00233 for( size_t i = 0; i < nofVertices; ++i ) 00234 { 00235 std::vector< double > coord; 00236 00237 DomainType p; 00238 for( int k = 0; k < DomainType::dimension; ++k ) 00239 p[ k ] = macroGrid_.vtx[i][k]; 00240 00241 vtxInsertOrder_.insert( std::make_pair( p, i ) ); 00242 } 00243 } 00244 } 00245 Grid *grid() 00246 { 00247 return grid_; 00248 } 00249 template <class Intersection> 00250 bool wasInserted(const Intersection &intersection) const 00251 { 00252 return intersection.boundary(); 00253 } 00254 template <class Intersection> 00255 int boundaryId(const Intersection &intersection) const 00256 { 00257 return intersection.boundaryId(); 00258 } 00259 00260 template< int codim > 00261 int numParameters () const 00262 { 00263 if( codim == 0 ) 00264 return macroGrid_.nofelparams; 00265 else if( codim == dimension ) 00266 return macroGrid_.nofvtxparams; 00267 else 00268 return 0; 00269 } 00270 00271 std::vector<double>& parameter(const Element &element) 00272 { 00273 const typename Element::Geometry &geo = element.geometry(); 00274 DomainType coord( geo.corner( 0 ) ); 00275 for( int i = 1; i < geo.corners(); ++i ) 00276 coord += geo.corner( i ); 00277 coord /= double( geo.corners() ); 00278 00279 InsertOrderIterator it = elInsertOrder_.find( coord ); 00280 if( it != elInsertOrder_.end() ) 00281 return macroGrid_.elParams[ it->second ]; 00282 assert(0); 00283 return emptyParam; 00284 } 00285 00286 std::vector<double>& parameter(const Vertex &vertex) 00287 { 00288 const typename Vertex::Geometry &geo = vertex.geometry(); 00289 DomainType coord( geo.corner( 0 ) ); 00290 00291 InsertOrderIterator it = vtxInsertOrder_.find( coord ); 00292 if( it != vtxInsertOrder_.end() ) 00293 return macroGrid_.vtxParams[ it->second ]; 00294 return emptyParam; 00295 } 00296 00297 private: 00298 typedef FieldVector<typename Grid::ctype,Grid::dimensionworld> DomainType; 00299 struct Compare 00300 { 00301 bool operator() ( const DomainType &a, const DomainType &b ) const 00302 { 00303 // returns true, if a < b; c[i] < -eps; 00304 const DomainType c = a - b; 00305 const double eps = 1e-8; 00306 00307 for( int i = 0; i < DomainType::dimension; ++i ) 00308 { 00309 if( c[ i ] <= -eps ) 00310 return true; 00311 if( c[ i ] >= eps ) 00312 return false; 00313 } 00314 return false; 00315 } 00316 }; 00317 typedef std::map< DomainType, size_t, Compare > InsertOrderMap; 00318 typedef typename InsertOrderMap::const_iterator InsertOrderIterator; 00319 MacroGrid macroGrid_; 00320 Grid *grid_; 00321 InsertOrderMap elInsertOrder_; 00322 InsertOrderMap vtxInsertOrder_; 00323 std::vector<double> emptyParam; 00324 }; 00325 00338 template< class GridType > 00339 struct GridPtr 00340 { 00341 typedef MPIHelper::MPICommunicator MPICommunicatorType; 00342 static const int dimension = GridType::dimension; 00343 00345 explicit GridPtr ( const std::string &filename, 00346 MPICommunicatorType comm = MPIHelper::getCommunicator() ) 00347 : gridPtr_( 0 ), 00348 elParam_( 0 ), 00349 vtxParam_( 0 ), 00350 bndId_( 0 ), 00351 nofElParam_( 0 ), 00352 nofVtxParam_( 0 ) 00353 { 00354 DGFGridFactory< GridType > dgfFactory( filename, comm ); 00355 initialize( dgfFactory ); 00356 } 00357 00359 explicit GridPtr ( std::istream &input, 00360 MPICommunicatorType comm = MPIHelper::getCommunicator() ) 00361 : gridPtr_( 0 ), 00362 elParam_( 0 ), 00363 vtxParam_( 0 ), 00364 bndId_( 0 ), 00365 nofElParam_( 0 ), 00366 nofVtxParam_( 0 ) 00367 { 00368 DGFGridFactory< GridType > dgfFactory( input, comm ); 00369 initialize( dgfFactory ); 00370 } 00371 00373 GridPtr() 00374 : gridPtr_(0), 00375 nofElParam_(0), 00376 nofVtxParam_(0) 00377 {} 00378 00380 GridPtr( GridType *grd ) 00381 : gridPtr_(grd), 00382 nofElParam_(0), 00383 nofVtxParam_(0) 00384 {} 00385 00387 GridPtr( const GridPtr &org ) 00388 : gridPtr_(org.gridPtr_), 00389 elParam_(org.elParam_), 00390 vtxParam_(org.vtxParam_), 00391 bndId_(org.bndId_), 00392 nofElParam_(org.nofElParam_), 00393 nofVtxParam_(org.nofVtxParam_) 00394 {} 00395 00397 GridPtr &operator= ( const GridPtr &org ) 00398 { 00399 gridPtr_ = org.gridPtr_; 00400 elParam_ = org.elParam_; 00401 vtxParam_ = org.vtxParam_; 00402 bndId_ = org.bndId_; 00403 nofVtxParam_ = org.nofVtxParam_; 00404 nofElParam_ = org.nofElParam_; 00405 return *this; 00406 } 00407 00409 GridPtr & operator = (GridType * grd) 00410 { 00411 gridPtr_ = std::auto_ptr<GridType>(grd); 00412 nofVtxParam_ = 0; 00413 nofElParam_ = 0; 00414 elParam_.resize(0); 00415 vtxParam_.resize(0); 00416 bndId_.resize(0); 00417 return *this; 00418 } 00419 00421 GridType& operator*() { 00422 return *gridPtr_; 00423 } 00424 00426 GridType* operator->() { 00427 return gridPtr_.operator -> (); 00428 } 00429 00431 const GridType& operator*() const { 00432 return *gridPtr_; 00433 } 00434 00436 const GridType* operator->() const { 00437 return gridPtr_.operator -> (); 00438 } 00439 00441 GridType* release () { 00442 return gridPtr_.release(); 00443 } 00444 00446 int nofParameters(int cdim) const { 00447 switch (cdim) { 00448 case 0: return nofElParam_; break; 00449 case GridType::dimension: return nofVtxParam_; break; 00450 } 00451 return 0; 00452 } 00454 template <class Entity> 00455 const std::vector< double > ¶meters ( const Entity &entity ) const 00456 { 00457 switch( (int)Entity::codimension ) 00458 { 00459 case 0: 00460 if( nofElParam_ > 0 ) 00461 { 00462 assert( (unsigned int)gridPtr_->leafView().indexSet().index( entity ) < elParam_.size() ); 00463 return elParam_[ gridPtr_->leafView().indexSet().index( entity ) ]; 00464 } 00465 break; 00466 case GridType::dimension: 00467 if( nofVtxParam_ > 0 ) 00468 { 00469 assert( (unsigned int)gridPtr_->leafView().indexSet().index( entity ) < vtxParam_.size() ); 00470 return vtxParam_[ gridPtr_->leafView().indexSet().index( entity ) ]; 00471 } 00472 break; 00473 } 00474 return emptyParam_; 00475 } 00476 00477 void loadBalance() 00478 { 00479 if ( gridPtr_->comm().size() == 1 ) 00480 return; 00481 int params = nofElParam_ + nofVtxParam_; 00482 if ( gridPtr_->comm().max( params ) > 0 ) 00483 { 00484 DataHandle dh(*this); 00485 gridPtr_->loadBalance( dh.interface() ); 00486 gridPtr_->communicate( dh.interface(), InteriorBorder_All_Interface,ForwardCommunication); 00487 } else 00488 { 00489 gridPtr_->loadBalance(); 00490 } 00491 } 00492 00493 protected: 00494 void initialize ( DGFGridFactory< GridType > &dgfFactory ) 00495 { 00496 gridPtr_ = std::auto_ptr< GridType >( dgfFactory.grid() ); 00497 00498 typedef typename GridType::LeafGridView GridView; 00499 GridView gridView = gridPtr_->leafView(); 00500 const typename GridView::IndexSet &indexSet = gridView.indexSet(); 00501 00502 nofElParam_ = dgfFactory.template numParameters< 0 >(); 00503 nofVtxParam_ = dgfFactory.template numParameters< dimension >(); 00504 if ( nofElParam_ > 0 ) 00505 elParam_.resize( indexSet.size(0) ); 00506 if ( nofVtxParam_ > 0 ) 00507 vtxParam_.resize( indexSet.size(dimension) ); 00508 bndId_.resize( indexSet.size(1) ); 00509 00510 const PartitionIteratorType partType = Interior_Partition; 00511 typedef typename GridView::template Codim< 0 >::template Partition< partType >::Iterator Iterator; 00512 const Iterator enditer = gridView.template end< 0, partType >(); 00513 for( Iterator iter = gridView.template begin< 0, partType >(); iter != enditer; ++iter ) 00514 { 00515 const typename Iterator::Entity &el = *iter; 00516 if ( nofElParam_ > 0 ) { 00517 std::swap( elParam_[ indexSet.index(el) ], dgfFactory.parameter(el) ); 00518 assert( elParam_[ indexSet.index(el) ].size() == (size_t)nofElParam_ ); 00519 } 00520 if ( nofVtxParam_ > 0 ) 00521 { 00522 for ( int v = 0; v < el.template count<dimension>(); ++v) 00523 { 00524 typename GridView::IndexSet::IndexType index = indexSet.subIndex(el,v,dimension); 00525 if ( vtxParam_[ index ].empty() ) 00526 std::swap( vtxParam_[ index ], dgfFactory.parameter(*el.template subEntity<dimension>(v) ) ); 00527 assert( vtxParam_[ index ].size() == (size_t)nofVtxParam_ ); 00528 } 00529 } 00530 if ( el.hasBoundaryIntersections() ) 00531 { 00532 typedef typename GridView::IntersectionIterator IntersectionIterator; 00533 const IntersectionIterator iend = gridView.iend(el); 00534 for( IntersectionIterator iiter = gridView.ibegin(el); iiter != iend; ++iiter ) 00535 { 00536 const typename IntersectionIterator::Intersection &inter = *iiter; 00537 // dirty hack: check for "none" to make corner point grid work 00538 if ( inter.boundary() && !inter.type().isNone() ) 00539 { 00540 bndId_[ indexSet.subIndex(el,inter.indexInInside(),1) ] 00541 = dgfFactory.boundaryId( inter ); 00542 } 00543 } 00544 } 00545 } 00546 } 00547 00548 template <class Entity> 00549 std::vector< double > ¶ms ( const Entity &entity ) 00550 { 00551 switch( (int)Entity::codimension ) 00552 { 00553 case 0: 00554 if( nofElParam_ > 0 ) { 00555 if ( gridPtr_->leafView().indexSet().index( entity ) >= elParam_.size() ) 00556 elParam_.resize( gridPtr_->leafView().indexSet().index( entity ) ); 00557 return elParam_[ gridPtr_->leafView().indexSet().index( entity ) ]; 00558 } 00559 break; 00560 case GridType::dimension: 00561 if( nofVtxParam_ > 0 ) { 00562 if ( gridPtr_->leafView().indexSet().index( entity ) >= vtxParam_.size() ) 00563 vtxParam_.resize( gridPtr_->leafView().indexSet().index( entity ) ); 00564 return vtxParam_[ gridPtr_->leafView().indexSet().index( entity ) ]; 00565 } 00566 break; 00567 } 00568 return emptyParam_; 00569 } 00570 void setNofParams( int cdim, int nofP ) 00571 { 00572 switch (cdim) { 00573 case 0: nofElParam_ = nofP; break; 00574 case GridType::dimension: nofVtxParam_ = nofP; break; 00575 } 00576 } 00577 struct DataHandle : public CommDataHandleIF<DataHandle,double> 00578 { 00579 DataHandle( GridPtr& gridPtr) : 00580 gridPtr_(gridPtr), 00581 idSet_(gridPtr->localIdSet()) 00582 { 00583 typedef typename GridType::LeafGridView GridView; 00584 GridView gridView = gridPtr_->leafView(); 00585 const typename GridView::IndexSet &indexSet = gridView.indexSet(); 00586 00587 const PartitionIteratorType partType = Interior_Partition; 00588 typedef typename GridView::template Codim< 0 >::template Partition< partType >::Iterator Iterator; 00589 const Iterator enditer = gridView.template end< 0, partType >(); 00590 for( Iterator iter = gridView.template begin< 0, partType >(); iter != enditer; ++iter ) 00591 { 00592 const typename Iterator::Entity &el = *iter; 00593 if ( gridPtr_.nofElParam_ > 0 ) 00594 std::swap( gridPtr_.elParam_[ indexSet.index(el) ], elData_[ idSet_.id(el) ] ); 00595 if ( gridPtr_.nofVtxParam_ > 0 ) 00596 { 00597 for ( int v = 0; v < el.template count<dimension>(); ++v) 00598 { 00599 typename GridView::IndexSet::IndexType index = indexSet.subIndex(el,v,dimension); 00600 if ( ! gridPtr_.vtxParam_[ index ].empty() ) 00601 std::swap( gridPtr_.vtxParam_[ index ], vtxData_[ idSet_.subId(el,v,dimension) ] ); 00602 } 00603 } 00604 } 00605 } 00606 ~DataHandle() 00607 { 00608 typedef typename GridType::LeafGridView GridView; 00609 GridView gridView = gridPtr_->leafView(); 00610 const typename GridView::IndexSet &indexSet = gridView.indexSet(); 00611 00612 if ( gridPtr_.nofElParam_ > 0 ) 00613 gridPtr_.elParam_.resize( indexSet.size(0) ); 00614 if ( gridPtr_.nofVtxParam_ > 0 ) 00615 gridPtr_.vtxParam_.resize( indexSet.size(dimension) ); 00616 00617 const PartitionIteratorType partType = All_Partition; 00618 typedef typename GridView::template Codim< 0 >::template Partition< partType >::Iterator Iterator; 00619 const Iterator enditer = gridView.template end< 0, partType >(); 00620 for( Iterator iter = gridView.template begin< 0, partType >(); iter != enditer; ++iter ) 00621 { 00622 const typename Iterator::Entity &el = *iter; 00623 if ( gridPtr_.nofElParam_ > 0 ) 00624 { 00625 std::swap( gridPtr_.elParam_[ indexSet.index(el) ], elData_[ idSet_.id(el) ] ); 00626 assert( gridPtr_.elParam_[ indexSet.index(el) ].size() == (unsigned int)gridPtr_.nofElParam_ ); 00627 } 00628 if ( gridPtr_.nofVtxParam_ > 0 ) 00629 { 00630 for ( int v = 0; v < el.template count<dimension>(); ++v) 00631 { 00632 typename GridView::IndexSet::IndexType index = indexSet.subIndex(el,v,dimension); 00633 if ( gridPtr_.vtxParam_[ index ].empty() ) 00634 std::swap( gridPtr_.vtxParam_[ index ], vtxData_[ idSet_.subId(el,v,dimension) ] ); 00635 assert( gridPtr_.vtxParam_[ index ].size() == (unsigned int)gridPtr_.nofVtxParam_ ); 00636 } 00637 } 00638 } 00639 } 00640 CommDataHandleIF<DataHandle,double> &interface() 00641 { 00642 return *this; 00643 } 00644 bool contains (int dim, int codim) const 00645 { 00646 return (codim==dim || codim==0); 00647 } 00648 bool fixedsize (int dim, int codim) const 00649 { 00650 return false; 00651 } 00652 template<class EntityType> 00653 size_t size (const EntityType& e) const 00654 { 00655 return gridPtr_.nofParameters(e.codimension); 00656 } 00657 template<class MessageBufferImp, class EntityType> 00658 void gather (MessageBufferImp& buff, const EntityType& e) const 00659 { 00660 const std::vector<double> &v = (e.codimension==0)? elData_[idSet_.id(e)]:vtxData_[idSet_.id(e)]; 00661 const size_t s = v.size(); 00662 for (size_t i=0;i<s;++i) 00663 buff.write( v[i] ); 00664 assert( s == (size_t)gridPtr_.nofParameters(e.codimension) ); 00665 } 00666 template<class MessageBufferImp, class EntityType> 00667 void scatter (MessageBufferImp& buff, const EntityType& e, size_t n) 00668 { 00669 std::vector<double> &v = (e.codimension==0)? elData_[idSet_.id(e)]:vtxData_[idSet_.id(e)]; 00670 v.resize( n ); 00671 gridPtr_.setNofParams( e.codimension, n ); 00672 for (size_t i=0;i<n;++i) 00673 buff.read( v[i] ); 00674 } 00675 private: 00676 typedef typename GridType::LocalIdSet IdSet; 00677 GridPtr &gridPtr_; 00678 const IdSet &idSet_; 00679 mutable std::map< typename IdSet::IdType, std::vector<double> > elData_, vtxData_; 00680 }; 00681 00682 // grid auto pointer 00683 mutable std::auto_ptr<GridType> gridPtr_; 00684 // element and vertex parameters 00685 std::vector< std::vector< double > > elParam_; 00686 std::vector< std::vector< double > > vtxParam_; 00687 std::vector< int > bndId_; 00688 int nofElParam_,nofVtxParam_; 00689 std::vector < double > emptyParam_; 00690 }; // end of class GridPtr 00691 00692 } 00693 01370 /* 01371 Dune::Alberta with \c dimworld=3: \n 01372 if Tetgen is used to construct a 01373 tetrahedral grid for Dune::Alberta then the bisection routine does 01374 not necessarily terminate. This problem does not occur 01375 if the grid is constructed using the \b Interval block. 01376 */ 01377 01378 namespace Dune { 01379 01382 template <class GridType> 01383 struct DGFGridInfo 01384 { 01386 static int refineStepsForHalf(); 01389 static double refineWeight(); 01390 }; 01391 01392 } // end namespace Dune 01393 #endif