dune-grid  2.1.1
io/file/dgfparser/dgfparser.hh
Go to the documentation of this file.
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 > &parameters ( 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 > &params ( 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