dune-grid  2.2.0
dgfalu.hh
Go to the documentation of this file.
00001 #ifndef DUNE_DGFPARSERALU_HH
00002 #define DUNE_DGFPARSERALU_HH
00003 
00004 // only include if ALUGrid is used
00005 #if defined ENABLE_ALUGRID
00006 
00007 #include <dune/grid/alugrid.hh>
00008 #include <dune/grid/io/file/dgfparser/dgfparser.hh>
00009 #include <dune/grid/io/file/dgfparser/parser.hh>
00010 #include <dune/grid/io/file/dgfparser/blocks/projection.hh>
00011 
00012 #include <dune/grid/common/intersection.hh>
00013 
00014 namespace Dune
00015 {
00016 
00017   // forward declaration
00018   // -------------------
00019 
00020   template< class GridImp, template< class > class IntersectionImp >
00021   class Intersection;
00022 
00023 // DGFGridInfo (specialization for ALUGrid)
00024 // ----------------------------------------
00025 
00027 template<>
00028 struct DGFGridInfo< ALUCubeGrid< 3, 3 > >
00029 {
00030   static int refineStepsForHalf () { return 1; }
00031   static double refineWeight () { return 0.125; }
00032 };
00033 
00034 template<>
00035 struct DGFGridInfo< ALUSimplexGrid< 3, 3 > >
00036 {
00037   static int refineStepsForHalf () { return 1; }
00038   static double refineWeight () { return 0.125; }
00039 };
00040 
00041 template<int dimw>
00042 struct DGFGridInfo< ALUSimplexGrid< 2, dimw > >
00043 {
00044   static int refineStepsForHalf () { return 1; }
00045   static double refineWeight () { return 0.25; }
00046 };
00047 
00048 template<int dimw>
00049 struct DGFGridInfo< ALUCubeGrid< 2, dimw > >
00050 {
00051   static int refineStepsForHalf () { return 1; }
00052   static double refineWeight () { return 0.25; }
00053 };
00054 
00055 template<int dimw>
00056 struct DGFGridInfo< Dune::ALUConformGrid< 2, dimw > >
00057 {
00058   static int refineStepsForHalf () { return 2; }
00059   static double refineWeight () { return 0.5; }
00060 };
00061 
00062 template<int dimg, int dimw, ALUGridElementType eltype, ALUGridRefinementType refinementtype, class Comm >
00063 struct DGFGridInfo< Dune::ALUGrid< dimg, dimw, eltype, refinementtype, Comm > >
00064 {
00065   static int refineStepsForHalf () { return ( refinementtype == conforming ) ? dimg : 1; }
00066   static double refineWeight () { return ( refinementtype == conforming ) ? 0.5 : 1.0/(std::pow( 2.0, double(dimg))); }
00067 };
00070   // DGFGridFactory for AluGrid
00071   // --------------------------
00072 
00073   // template< int dim, int dimworld > // for a first version
00074   template < class G >
00075   struct DGFBaseFactory
00076   {
00077     typedef G  Grid;
00078     const static int dimension = Grid::dimension;
00079     typedef MPIHelper::MPICommunicator MPICommunicatorType;
00080     typedef typename Grid::template Codim<0>::Entity Element;
00081     typedef typename Grid::template Codim<dimension>::Entity Vertex;
00082     typedef Dune::GridFactory<Grid> GridFactory;
00083 
00084     DGFBaseFactory ()
00085       : factory_( ),
00086         dgf_( 0, 1 )
00087     {
00088     }
00089 
00090     explicit DGFBaseFactory ( MPICommunicatorType comm ) 
00091       : factory_(),
00092         dgf_( rank(comm), size(comm) )
00093     {
00094     }
00095 
00096     Grid *grid () const
00097     {
00098       return grid_;
00099     }
00100 
00101     template< class Intersection >
00102     bool wasInserted ( const Intersection &intersection ) const
00103     {
00104       return factory_.wasInserted( intersection );
00105     }
00106 
00107     template < class GG, template < class > class II >
00108     int boundaryId ( const Intersection< GG, II > & intersection ) const
00109     {
00110       typedef Dune::Intersection< GG, II > Intersection;
00111       typename Intersection::EntityPointer inside = intersection.inside();
00112       const typename Intersection::Entity & entity = *inside;
00113       const int face = intersection.indexInInside();
00114 
00115       const GenericReferenceElement< double, dimension > & refElem = 
00116         GenericReferenceElements< double, dimension >::general( entity.type() );
00117       int corners = refElem.size( face, 1, dimension );
00118       std :: vector< unsigned int > bound( corners );
00119       for( int i=0; i < corners; ++i )
00120       {
00121         const int k =  refElem.subEntity( face, 1, i, dimension );
00122         bound[ i ] = factory_.insertionIndex( *entity.template subEntity< dimension >( k ) );
00123       }
00124 
00125       DuneGridFormatParser::facemap_t::key_type key( bound, false );
00126       const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
00127       if( pos != dgf_.facemap.end() )
00128         return dgf_.facemap.find( key )->second.first;
00129       else
00130         return (intersection.boundary() ? 1 : 0);
00131     }
00132 
00133     template < class GG, template < class > class II >
00134     const typename DGFBoundaryParameter::type & 
00135       boundaryParameter ( const Intersection< GG, II > & intersection ) const
00136     {
00137       typedef Dune::Intersection< GG, II > Intersection;
00138       typename Intersection::EntityPointer inside = intersection.inside();
00139       const typename Intersection::Entity & entity = *inside;
00140       const int face = intersection.indexInInside();
00141 
00142       const GenericReferenceElement< double, dimension > & refElem = 
00143         GenericReferenceElements< double, dimension >::general( entity.type() );
00144       int corners = refElem.size( face, 1, dimension );
00145       std :: vector< unsigned int > bound( corners );
00146       for( int i=0; i < corners; ++i )
00147       {
00148         const int k =  refElem.subEntity( face, 1, i, dimension );
00149         bound[ i ] = factory_.insertionIndex( *entity.template subEntity< dimension >( k ) );
00150       }
00151 
00152       DuneGridFormatParser::facemap_t::key_type key( bound, false );
00153       const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
00154       if( pos != dgf_.facemap.end() )
00155         return dgf_.facemap.find( key )->second.second;
00156       else
00157         return DGFBoundaryParameter::defaultValue();
00158     }
00159 
00160     template< int codim >
00161     int numParameters () const
00162     {
00163       if( codim == 0 )
00164         return dgf_.nofelparams;
00165       else if( codim == dimension )
00166         return dgf_.nofvtxparams;
00167       else
00168         return 0;
00169     }
00170 
00171     // return true if boundary parameters found
00172     bool haveBoundaryParameters () const
00173     {
00174       return dgf_.haveBndParameters;
00175     }
00176 
00177     std::vector< double > &parameter ( const Element &element )
00178     {
00179       if( numParameters< 0 >() <= 0 )
00180       {
00181         DUNE_THROW( InvalidStateException,
00182                     "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
00183       }
00184       return dgf_.elParams[ factory_.insertionIndex( element ) ];
00185     }
00186 
00187     std::vector< double > &parameter ( const Vertex &vertex )
00188     {
00189       if( numParameters< dimension >() <= 0 )
00190       {
00191         DUNE_THROW( InvalidStateException,
00192                     "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
00193       }
00194       return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
00195     }
00196 
00197   protected:
00198     bool generateALUGrid( const ALUGridElementType eltype,
00199                           std::istream &file, 
00200                           MPICommunicatorType communicator,
00201                           const std::string &filename );
00202 
00203     bool generateALU2dGrid( const ALUGridElementType eltype,
00204                             std::istream &file, 
00205                             MPICommunicatorType communicator,
00206                             const std::string &filename );
00207 
00208     static Grid* callDirectly( const char* gridname,
00209                                const int rank, 
00210                                const char *filename,
00211                                MPICommunicatorType communicator )
00212     {
00213   #if ALU3DGRID_PARALLEL
00214       // in parallel runs add rank to filename 
00215       std :: stringstream tmps;
00216       tmps << filename << "." << rank;
00217       const std :: string &tmp = tmps.str();
00218 
00219       // if file exits then use it 
00220       if( fileExists( tmp.c_str() ) )
00221         return new Grid( tmp.c_str(), communicator );
00222   #endif
00223       // for rank 0 we also check the normal file name 
00224       if( rank == 0 ) 
00225       {
00226         if( fileExists( filename ) )
00227           return new Grid( filename , communicator );
00228 
00229         // only throw this exception on rank 0 because 
00230         // for the other ranks we can still create empty grids 
00231         DUNE_THROW( GridError, "Unable to create " << gridname << " from '"
00232                     << filename << "'." );
00233       }
00234       else 
00235         dwarn << "WARNING:  P[" << rank << "]: Creating empty grid!" << std::endl;
00236 
00237       // return empty grid on all other processes 
00238       return new Grid( communicator );
00239     }
00240     static bool fileExists ( const char *fileName )
00241     {
00242       std :: ifstream testfile( fileName );
00243       if( !testfile )
00244         return false;
00245       testfile.close();
00246       return true;
00247     }
00248     static int rank( MPICommunicatorType MPICOMM )
00249     {
00250       int rank = 0;
00251 #if HAVE_MPI
00252       MPI_Comm_rank( MPICOMM, &rank );
00253 #endif
00254       return rank;
00255     }
00256     static int size( MPICommunicatorType MPICOMM )
00257     {
00258       int size = 1;
00259 #if HAVE_MPI
00260       MPI_Comm_size( MPICOMM, &size );
00261 #endif
00262       return size;
00263     }
00264     Grid *grid_;
00265     GridFactory factory_;
00266     DuneGridFormatParser dgf_;
00267   };
00268 
00269   template <>
00270   struct DGFGridFactory< ALUSimplexGrid<3,3> > : 
00271     public DGFBaseFactory< ALUSimplexGrid<3,3> >
00272   {
00273     explicit DGFGridFactory ( std::istream &input,
00274                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00275     : DGFBaseFactory< ALUSimplexGrid<3,3> >( comm )
00276     {
00277       input.clear();
00278       input.seekg( 0 );
00279       if( !input )
00280         DUNE_THROW( DGFException, "Error resetting input stream." );
00281       generate( input, comm );
00282     }
00283 
00284     explicit DGFGridFactory ( const std::string &filename, 
00285                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00286     : DGFBaseFactory< ALUSimplexGrid<3,3> >( comm ) 
00287     {
00288       std::ifstream input( filename.c_str() );
00289 
00290       bool fileFound = input.is_open() ;
00291       if( fileFound ) 
00292         fileFound = generate( input, comm, filename );
00293 
00294       if( ! fileFound ) 
00295         grid_ = callDirectly( "ALUSimplexGrid< 3 , 3 >", this->rank( comm ), filename.c_str(), comm );
00296     }
00297 
00298   protected:
00299     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00300   };
00301 
00302   template <>
00303   struct DGFGridFactory< ALUCubeGrid<3,3> > : 
00304     public DGFBaseFactory< ALUCubeGrid<3,3> >
00305   {
00306     explicit DGFGridFactory ( std::istream &input,
00307                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00308     : DGFBaseFactory< ALUCubeGrid<3,3> >( comm )
00309     {
00310       input.clear();
00311       input.seekg( 0 );
00312       if( !input )
00313         DUNE_THROW( DGFException, "Error resetting input stream." );
00314       generate( input, comm );
00315     }
00316 
00317     explicit DGFGridFactory ( const std::string &filename, 
00318                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00319     : DGFBaseFactory< ALUCubeGrid<3,3> >( comm ) 
00320     {
00321       std::ifstream input( filename.c_str() );
00322       bool fileFound = input.is_open() ;
00323       if( fileFound ) 
00324         fileFound = generate( input, comm, filename );
00325 
00326       if( ! fileFound ) 
00327         grid_ = callDirectly( "ALUCubeGrid< 3 , 3 >", this->rank( comm ), filename.c_str(), comm );
00328     }
00329 
00330   protected:
00331     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00332   };
00333 
00334   template < ALUGridElementType eltype, ALUGridRefinementType refinementtype, class Comm >
00335   struct DGFGridFactory< ALUGrid<3,3, eltype, refinementtype, Comm > > : 
00336     public DGFBaseFactory< ALUGrid<3,3, eltype, refinementtype, Comm > >
00337   {
00338     typedef ALUGrid<3,3, eltype, refinementtype, Comm > DGFGridType;
00339     typedef DGFBaseFactory< DGFGridType > BaseType;
00340     typedef typename BaseType :: MPICommunicatorType MPICommunicatorType;
00341   protected:  
00342     using BaseType :: grid_;
00343     using BaseType :: callDirectly;
00344   public:  
00345     explicit DGFGridFactory ( std::istream &input,
00346                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00347     : BaseType( comm )
00348     {
00349       input.clear();
00350       input.seekg( 0 );
00351       if( !input )
00352         DUNE_THROW( DGFException, "Error resetting input stream." );
00353       generate( input, comm );
00354     }
00355 
00356     explicit DGFGridFactory ( const std::string &filename, 
00357                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00358     : BaseType( comm ) 
00359     {
00360       std::ifstream input( filename.c_str() );
00361       bool fileFound = input.is_open() ;
00362       if( fileFound ) 
00363         fileFound = generate( input, comm, filename );
00364 
00365       if( ! fileFound ) 
00366         grid_ = callDirectly( "ALUGrid< 3, 3, eltype, ref, comm >", this->rank( comm ), filename.c_str(), comm );
00367     }
00368 
00369   protected:
00370     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00371   };
00372 
00373   template <int dimw>
00374   struct DGFGridFactory< ALUConformGrid<2,dimw> > : 
00375     public DGFBaseFactory< ALUConformGrid<2,dimw> >
00376   {
00377     typedef DGFBaseFactory< ALUConformGrid<2,dimw> > Base;
00378     typedef typename Base:: MPICommunicatorType  MPICommunicatorType;
00379     typedef typename Base::Grid Grid;
00380     const static int dimension = Grid::dimension;
00381     typedef typename Base::GridFactory GridFactory;
00382 
00383     explicit DGFGridFactory ( std::istream &input,
00384                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00385     : Base() 
00386     {
00387       input.clear();
00388       input.seekg( 0 );
00389       if( !input )
00390         DUNE_THROW( DGFException, "Error resetting input stream." );
00391       generate( input, comm );
00392     }
00393 
00394     explicit DGFGridFactory ( const std::string &filename, 
00395                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00396     : DGFBaseFactory< ALUConformGrid<2,dimw> >()
00397     {
00398       std::ifstream input( filename.c_str() );
00399       if( !input )
00400         DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
00401       if( !generate( input, comm, filename ) )
00402       {
00403         if( Base::fileExists( filename.c_str() ) )
00404           Base::grid_ = new Grid( filename );
00405         else
00406           DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
00407       }
00408     }
00409   protected:
00410     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00411     using Base::grid_;
00412     using Base::factory_;
00413     using Base::dgf_;
00414   };
00415 
00416   template <int dimw>
00417   struct DGFGridFactory< ALUSimplexGrid<2,dimw> > : 
00418     public DGFBaseFactory< ALUSimplexGrid<2,dimw> >
00419   {
00420     typedef DGFBaseFactory< ALUSimplexGrid<2,dimw> > Base;
00421     typedef typename Base::MPICommunicatorType  MPICommunicatorType;
00422     typedef typename Base::Grid Grid;
00423     const static int dimension = Grid::dimension;
00424     typedef typename Base::GridFactory GridFactory;
00425 
00426     explicit DGFGridFactory ( std::istream &input,
00427                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00428     : Base() 
00429     {
00430       input.clear();
00431       input.seekg( 0 );
00432       if( !input )
00433         DUNE_THROW(DGFException, "Error resetting input stream." );
00434       generate( input, comm );
00435     }
00436 
00437     explicit DGFGridFactory ( const std::string &filename, 
00438                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00439     : DGFBaseFactory< ALUSimplexGrid<2,dimw> >()
00440     {
00441       std::ifstream input( filename.c_str() );
00442       if( !input )
00443         DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
00444       if( !generate( input, comm, filename ) )
00445       {
00446         if( Base::fileExists( filename.c_str() ) )
00447           Base::grid_ = new Grid( filename );
00448         else
00449           DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
00450       }
00451     }
00452 
00453   protected:
00454     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00455     using Base::grid_;
00456     using Base::factory_;
00457     using Base::dgf_;
00458   };
00459 
00460   template <int dimw>
00461   struct DGFGridFactory< ALUCubeGrid<2,dimw> > : 
00462     public DGFBaseFactory< ALUCubeGrid<2,dimw> >
00463   {
00464     typedef DGFBaseFactory< ALUCubeGrid<2,dimw> > Base;
00465     typedef typename Base::MPICommunicatorType  MPICommunicatorType;
00466     typedef typename Base::Grid Grid;
00467     const static int dimension = Grid::dimension;
00468     typedef typename Base::GridFactory GridFactory;
00469 
00470     explicit DGFGridFactory ( std::istream &input,
00471                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00472     : Base() 
00473     {
00474       input.clear();
00475       input.seekg( 0 );
00476       if( !input )
00477         DUNE_THROW(DGFException, "Error resetting input stream." );
00478       generate( input, comm );
00479     }
00480 
00481     explicit DGFGridFactory ( const std::string &filename, 
00482                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00483     : DGFBaseFactory< ALUCubeGrid<2,dimw> >()
00484     {
00485       std::ifstream input( filename.c_str() );
00486       if( !input )
00487         DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
00488       if( !generate( input, comm, filename ) )
00489       {
00490         if( Base::fileExists( filename.c_str() ) )
00491           Base::grid_ = new Grid( filename );
00492         else
00493           DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
00494       }
00495     }
00496 
00497   protected:
00498     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00499     using Base::grid_;
00500     using Base::factory_;
00501     using Base::dgf_;
00502   };
00503 
00504   template < int dimw, ALUGridElementType eltype, ALUGridRefinementType refinementtype, class Comm >
00505   struct DGFGridFactory< ALUGrid<2, dimw, eltype, refinementtype, Comm > > : 
00506     public DGFBaseFactory< ALUGrid< 2, dimw, eltype, refinementtype, Comm > >
00507   {
00508     typedef ALUGrid< 2, dimw, eltype, refinementtype, Comm > DGFGridType;
00509     typedef DGFBaseFactory< DGFGridType > BaseType;
00510     typedef typename BaseType :: MPICommunicatorType MPICommunicatorType;
00511     typedef typename BaseType::Grid Grid;
00512     const static int dimension = Grid::dimension;
00513     typedef typename BaseType::GridFactory GridFactory;
00514 
00515     explicit DGFGridFactory ( std::istream &input,
00516                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00517     : BaseType( comm ) 
00518     {
00519       input.clear();
00520       input.seekg( 0 );
00521       if( !input )
00522         DUNE_THROW(DGFException, "Error resetting input stream." );
00523       generate( input, comm );
00524     }
00525 
00526     explicit DGFGridFactory ( const std::string &filename, 
00527                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00528       : BaseType( comm )
00529     {
00530       std::ifstream input( filename.c_str() );
00531       if( !input )
00532         DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
00533       if( !generate( input, comm, filename ) )
00534       {
00535         if( BaseType::fileExists( filename.c_str() ) )
00536           grid_ = new Grid( filename );
00537         else
00538           DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
00539       }
00540     }
00541 
00542   protected:
00543     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00544     using BaseType::grid_;
00545     using BaseType::factory_;
00546     using BaseType::dgf_;
00547   };
00548 
00549 }
00550 
00551 #include "dgfalu.cc"
00552 
00553 #endif // #if defined ENABLE_ALUGRID
00554 
00555 #endif // #ifndef DUNE_DGFPARSERALU_HH