dune-grid  2.2.0
albertagrid/dgfparser.hh
Go to the documentation of this file.
00001 #ifndef DUNE_ALBERTA_DGFPARSER_HH
00002 #define DUNE_ALBERTA_DGFPARSER_HH
00003 
00004 #include <vector>
00005 
00006 #include <dune/grid/albertagrid.hh>
00007 #include <dune/grid/albertagrid/gridfactory.hh>
00008 
00009 #include <dune/grid/io/file/dgfparser/dgfparser.hh>
00010 #include <dune/grid/io/file/dgfparser/blocks/projection.hh>
00011 
00012 #include <dune/grid/common/intersection.hh>
00013 #include <dune/grid/io/file/dgfparser/parser.hh>
00014 
00015 #if HAVE_ALBERTA 
00016 
00017 namespace Dune
00018 {
00019 
00020   // forward declaration
00021   // -------------------
00022 
00023   template < class GridImp, template < class > class IntersectionImp >
00024   class Intersection;
00025 
00026   // DGFGridFactory for AlbertaGrid
00027   // ------------------------------
00028 
00029   template< int dim, int dimworld >
00030   struct DGFGridFactory< AlbertaGrid< dim, dimworld > >
00031   {
00032     typedef AlbertaGrid<dim,dimworld>  Grid;
00033     const static int dimension = Grid::dimension;
00034     typedef MPIHelper::MPICommunicator MPICommunicatorType;
00035     typedef typename Grid::template Codim<0>::Entity Element;
00036     typedef typename Grid::template Codim<dimension>::Entity Vertex;
00037     typedef Dune::GridFactory<Grid> GridFactory;
00038 
00039     explicit DGFGridFactory ( std::istream &input,
00040                               MPICommunicatorType comm = MPIHelper::getCommunicator() );
00041     explicit DGFGridFactory ( const std::string &filename, 
00042                               MPICommunicatorType comm = MPIHelper::getCommunicator() );
00043 
00044     Grid *grid () const
00045     {
00046       return grid_;
00047     }
00048 
00049     template< class Intersection >
00050     bool wasInserted ( const Intersection &intersection ) const
00051     {
00052       return factory_.wasInserted( intersection );
00053     }
00054 
00055     template< class Intersection >
00056     int boundaryId ( const Intersection &intersection ) const
00057     {
00058       return Grid::getRealImplementation( intersection ).boundaryId();
00059     }
00060 
00061     // return true if boundary paramters found
00062     bool haveBoundaryParameters () const
00063     {
00064       return dgf_.haveBndParameters;
00065     }
00066 
00067     template < class GG, template < class > class II >
00068     const DGFBoundaryParameter::type &
00069       boundaryParameter ( const Intersection< GG, II > & intersection ) const
00070     {
00071       typedef Dune::Intersection< GG, II > Intersection;
00072       typename Intersection::EntityPointer inside = intersection.inside();
00073       const typename Intersection::Entity & entity = *inside;
00074       const int face = intersection.indexInInside();
00075 
00076       const GenericReferenceElement< double, dimension > & refElem =
00077         GenericReferenceElements< double, dimension >::general( entity.type() );
00078       int corners = refElem.size( face, 1, dimension );
00079       std :: vector< unsigned int > bound( corners );
00080       for( int i=0; i < corners; ++i )
00081       {
00082         const int k =  refElem.subEntity( face, 1, i, dimension );
00083         bound[ i ] = factory_.insertionIndex( *entity.template subEntity< dimension >( k ) );
00084       }
00085 
00086       DuneGridFormatParser::facemap_t::key_type key( bound, false );
00087       const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
00088       if( pos != dgf_.facemap.end() )
00089         return dgf_.facemap.find( key )->second.second;
00090       else
00091         return DGFBoundaryParameter::defaultValue(); 
00092     }
00093 
00094     template< int codim >
00095     int numParameters () const
00096     {
00097       if( codim == 0 )
00098         return dgf_.nofelparams;
00099       else if( codim == dimension )
00100         return dgf_.nofvtxparams;
00101       else
00102         return 0;
00103     }
00104 
00105     std::vector< double > &parameter ( const Element &element )
00106     {
00107       if( numParameters< 0 >() <= 0 )
00108       {
00109         DUNE_THROW( InvalidStateException,
00110                     "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
00111       }
00112       return dgf_.elParams[ factory_.insertionIndex( element ) ];
00113     }
00114 
00115     std::vector< double > &parameter ( const Vertex &vertex )
00116     {
00117       if( numParameters< dimension >() <= 0 )
00118       {
00119         DUNE_THROW( InvalidStateException,
00120                     "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
00121       }
00122       return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
00123     }
00124 
00125   private:
00126     bool generate( std::istream &input );
00127 
00128     Grid *grid_;
00129     GridFactory factory_;
00130     DuneGridFormatParser dgf_;
00131   };
00132 
00133 
00134 
00135   // DGFGridInfo for AlbertaGrid
00136   // ---------------------------
00137 
00138   template< int dim, int dimworld >
00139   struct DGFGridInfo< AlbertaGrid< dim, dimworld > >
00140   {
00141     static int refineStepsForHalf ()
00142     {
00143       return dim;
00144     }
00145 
00146     static double refineWeight ()
00147     {
00148       return 0.5;
00149     }
00150   };
00151 
00152 
00153 
00154   // Implementation of DGFGridFactory for AlbertaGrid
00155   // ------------------------------------------------
00156 
00157   template< int dim, int dimworld >
00158   inline DGFGridFactory< AlbertaGrid< dim, dimworld > >
00159     ::DGFGridFactory ( std::istream &input, MPICommunicatorType comm )
00160   : dgf_( 0, 1 )
00161   {
00162     input.clear();
00163     input.seekg( 0 );
00164     if( !input )
00165       DUNE_THROW(DGFException, "Error resetting input stream." );
00166     generate( input );
00167   }
00168 
00169 
00170   template< int dim, int dimworld >
00171   inline DGFGridFactory< AlbertaGrid< dim, dimworld > >
00172     ::DGFGridFactory ( const std::string &filename, MPICommunicatorType comm )
00173   : dgf_( 0, 1 )
00174   {
00175     std::ifstream input( filename.c_str() );
00176     if( !input )
00177       DUNE_THROW( DGFException, "Macrofile " << filename << " not found." );
00178     if( !generate( input ) )
00179       grid_ = new AlbertaGrid< dim, dimworld >( filename.c_str() );
00180     input.close();
00181   }
00182 
00183 }
00184 
00185 #endif // #if HAVE_ALBERTA
00186 
00187 #endif // #ifndef DUNE_ALBERTA_DGFPARSER_HH