dune-grid  2.1.1
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/dgfprojectionblock.hh>
00010 
00011 namespace Dune
00012 {
00013 
00014 // DGFGridInfo (specialization for ALUGrid)
00015 // ----------------------------------------
00016 
00018 template<>
00019 struct DGFGridInfo< ALUCubeGrid< 3, 3 > >
00020 {
00021   static int refineStepsForHalf () { return 1; }
00022   static double refineWeight () { return 0.125; }
00023 };
00024 
00025 template<>
00026 struct DGFGridInfo< ALUSimplexGrid< 3, 3 > >
00027 {
00028   static int refineStepsForHalf () { return 1; }
00029   static double refineWeight () { return 0.125; }
00030 };
00031 
00032 template<int dimw>
00033 struct DGFGridInfo< ALUSimplexGrid< 2, dimw > >
00034 {
00035   static int refineStepsForHalf () { return 1; }
00036   static double refineWeight () { return 0.25; }
00037 };
00038 
00039 template<int dimw>
00040 struct DGFGridInfo< ALUCubeGrid< 2, dimw > >
00041 {
00042   static int refineStepsForHalf () { return 1; }
00043   static double refineWeight () { return 0.25; }
00044 };
00045 
00046 template<int dimw>
00047 struct DGFGridInfo< Dune::ALUConformGrid< 2, dimw > >
00048 {
00049   static int refineStepsForHalf () { return 2; }
00050   static double refineWeight () { return 0.5; }
00051 };
00054   // DGFGridFactory for AluSimplexGrid
00055   // ------------------------------
00056 
00057   // template< int dim, int dimworld > // for a first version
00058   template < class G >
00059   struct DGFBaseFactory
00060   {
00061     typedef G  Grid;
00062     const static int dimension = Grid::dimension;
00063     typedef MPIHelper::MPICommunicator MPICommunicatorType;
00064     typedef typename Grid::template Codim<0>::Entity Element;
00065     typedef typename Grid::template Codim<dimension>::Entity Vertex;
00066     typedef Dune::GridFactory<Grid> GridFactory;
00067 
00068     DGFBaseFactory ()
00069       : factory_( ),
00070         dgf_( 0, 1 )
00071     {
00072     }
00073 
00074     explicit DGFBaseFactory ( MPICommunicatorType comm ) 
00075       : factory_(comm),
00076         dgf_( rank(comm), size(comm) )
00077     {
00078     }
00079 
00080     Grid *grid () const
00081     {
00082       return grid_;
00083     }
00084 
00085     template< class Intersection >
00086     bool wasInserted ( const Intersection &intersection ) const
00087     {
00088       return factory_.wasInserted( intersection );
00089     }
00090 
00091     template< class Intersection >
00092     int boundaryId ( const Intersection &intersection ) const
00093     {
00094       return intersection.boundaryId();
00095     }
00096 
00097     template< int codim >
00098     int numParameters () const
00099     {
00100       if( codim == 0 )
00101         return dgf_.nofelparams;
00102       else if( codim == dimension )
00103         return dgf_.nofvtxparams;
00104       else
00105         return 0;
00106     }
00107 
00108     std::vector< double > &parameter ( const Element &element )
00109     {
00110       if( numParameters< 0 >() <= 0 )
00111       {
00112         DUNE_THROW( InvalidStateException,
00113                     "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
00114       }
00115       return dgf_.elParams[ factory_.insertionIndex( element ) ];
00116     }
00117 
00118     std::vector< double > &parameter ( const Vertex &vertex )
00119     {
00120       if( numParameters< dimension >() <= 0 )
00121       {
00122         DUNE_THROW( InvalidStateException,
00123                     "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
00124       }
00125       return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
00126     }
00127 
00128   protected:
00129     static Grid* callDirectly( const char* gridname,
00130                                const int rank, 
00131                                const char *filename,
00132                                MPICommunicatorType communicator )
00133     {
00134   #if ALU3DGRID_PARALLEL
00135       // in parallel runs add rank to filename 
00136       std :: stringstream tmps;
00137       tmps << filename << "." << rank;
00138       const std :: string &tmp = tmps.str();
00139 
00140       // if file exits then use it 
00141       if( fileExists( tmp.c_str() ) )
00142         return new Grid( tmp.c_str(), communicator );
00143   #endif
00144       // for rank 0 we also check the normal file name 
00145       if( rank == 0 ) 
00146       {
00147         if( fileExists( filename ) )
00148           return new Grid( filename , communicator );
00149 
00150         // only throw this exception on rank 0 because 
00151         // for the other ranks we can still create empty grids 
00152         DUNE_THROW( GridError, "Unable to create " << gridname << " from '"
00153                     << filename << "'." );
00154       }
00155       else 
00156         dwarn << "WARNING:  P[" << rank << "]: Creating empty grid!" << std::endl;
00157 
00158       // return empty grid on all other processes 
00159       return new Grid( communicator );
00160     }
00161     static bool fileExists ( const char *fileName )
00162     {
00163       std :: ifstream testfile( fileName );
00164       if( !testfile )
00165         return false;
00166       testfile.close();
00167       return true;
00168     }
00169     static int rank( MPICommunicatorType MPICOMM )
00170     {
00171       int rank = 0;
00172 #if HAVE_MPI
00173       MPI_Comm_rank( MPICOMM, &rank );
00174 #endif
00175       return rank;
00176     }
00177     static int size( MPICommunicatorType MPICOMM )
00178     {
00179       int size = 1;
00180 #if HAVE_MPI
00181       MPI_Comm_size( MPICOMM, &size );
00182 #endif
00183       return size;
00184     }
00185     Grid *grid_;
00186     GridFactory factory_;
00187     DuneGridFormatParser dgf_;
00188   };
00189 
00190   template <>
00191   struct DGFGridFactory< ALUSimplexGrid<3,3> > : 
00192     public DGFBaseFactory< ALUSimplexGrid<3,3> >
00193   {
00194     explicit DGFGridFactory ( std::istream &input,
00195                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00196     : DGFBaseFactory< ALUSimplexGrid<3,3> >( comm )
00197     {
00198       input.clear();
00199       input.seekg( 0 );
00200       if( !input )
00201         DUNE_THROW( DGFException, "Error resetting input stream." );
00202       generate( input, comm );
00203     }
00204 
00205     explicit DGFGridFactory ( const std::string &filename, 
00206                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00207     : DGFBaseFactory< ALUSimplexGrid<3,3> >( comm ) 
00208     {
00209       std::ifstream input( filename.c_str() );
00210 
00211       bool fileFound = input.is_open() ;
00212       if( fileFound ) 
00213         fileFound = generate( input, comm, filename );
00214 
00215       if( ! fileFound ) 
00216         grid_ = callDirectly( "ALUSimplexGrid< 3 , 3 >", rank( comm ), filename.c_str(), comm );
00217     }
00218 
00219   protected:
00220     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00221   };
00222 
00223   template <>
00224   struct DGFGridFactory< ALUCubeGrid<3,3> > : 
00225     public DGFBaseFactory< ALUCubeGrid<3,3> >
00226   {
00227     explicit DGFGridFactory ( std::istream &input,
00228                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00229     : DGFBaseFactory< ALUCubeGrid<3,3> >( comm )
00230     {
00231       input.clear();
00232       input.seekg( 0 );
00233       if( !input )
00234         DUNE_THROW( DGFException, "Error resetting input stream." );
00235       generate( input, comm );
00236     }
00237 
00238     explicit DGFGridFactory ( const std::string &filename, 
00239                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00240     : DGFBaseFactory< ALUCubeGrid<3,3> >( comm ) 
00241     {
00242       std::ifstream input( filename.c_str() );
00243       bool fileFound = input.is_open() ;
00244       if( fileFound ) 
00245         fileFound = generate( input, comm, filename );
00246 
00247       if( ! fileFound ) 
00248         grid_ = callDirectly( "ALUCubeGrid< 3 , 3 >", rank( comm ), filename.c_str(), comm );
00249     }
00250 
00251   protected:
00252     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00253   };
00254 
00255   template <int dimw>
00256   struct DGFGridFactory< ALUConformGrid<2,dimw> > : 
00257     public DGFBaseFactory< ALUConformGrid<2,dimw> >
00258   {
00259     typedef DGFBaseFactory< ALUConformGrid<2,dimw> > Base;
00260     typedef typename Base:: MPICommunicatorType  MPICommunicatorType;
00261     typedef typename Base::Grid Grid;
00262     const static int dimension = Grid::dimension;
00263     typedef typename Base::GridFactory GridFactory;
00264 
00265     explicit DGFGridFactory ( std::istream &input,
00266                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00267     : Base() 
00268     {
00269       input.clear();
00270       input.seekg( 0 );
00271       if( !input )
00272         DUNE_THROW( DGFException, "Error resetting input stream." );
00273       generate( input, comm );
00274     }
00275 
00276     explicit DGFGridFactory ( const std::string &filename, 
00277                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00278     : DGFBaseFactory< ALUConformGrid<2,dimw> >()
00279     {
00280       std::ifstream input( filename.c_str() );
00281       if( !input )
00282         DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
00283       if( !generate( input, comm, filename ) )
00284       {
00285         if( Base::fileExists( filename.c_str() ) )
00286           Base::grid_ = new Grid( filename );
00287         else
00288           DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
00289       }
00290     }
00291   protected:
00292     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00293     using Base::grid_;
00294     using Base::factory_;
00295     using Base::dgf_;
00296   };
00297 
00298   template <int dimw>
00299   struct DGFGridFactory< ALUSimplexGrid<2,dimw> > : 
00300     public DGFBaseFactory< ALUSimplexGrid<2,dimw> >
00301   {
00302     typedef DGFBaseFactory< ALUSimplexGrid<2,dimw> > Base;
00303     typedef typename Base::MPICommunicatorType  MPICommunicatorType;
00304     typedef typename Base::Grid Grid;
00305     const static int dimension = Grid::dimension;
00306     typedef typename Base::GridFactory GridFactory;
00307 
00308     explicit DGFGridFactory ( std::istream &input,
00309                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00310     : Base() 
00311     {
00312       input.clear();
00313       input.seekg( 0 );
00314       if( !input )
00315         DUNE_THROW(DGFException, "Error resetting input stream." );
00316       generate( input, comm );
00317     }
00318 
00319     explicit DGFGridFactory ( const std::string &filename, 
00320                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00321     : DGFBaseFactory< ALUSimplexGrid<2,dimw> >()
00322     {
00323       std::ifstream input( filename.c_str() );
00324       if( !input )
00325         DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
00326       if( !generate( input, comm, filename ) )
00327       {
00328         if( Base::fileExists( filename.c_str() ) )
00329           Base::grid_ = new Grid( filename );
00330         else
00331           DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
00332       }
00333     }
00334 
00335   protected:
00336     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00337     using Base::grid_;
00338     using Base::factory_;
00339     using Base::dgf_;
00340   };
00341 
00342   template <int dimw>
00343   struct DGFGridFactory< ALUCubeGrid<2,dimw> > : 
00344     public DGFBaseFactory< ALUCubeGrid<2,dimw> >
00345   {
00346     typedef DGFBaseFactory< ALUCubeGrid<2,dimw> > Base;
00347     typedef typename Base::MPICommunicatorType  MPICommunicatorType;
00348     typedef typename Base::Grid Grid;
00349     const static int dimension = Grid::dimension;
00350     typedef typename Base::GridFactory GridFactory;
00351 
00352     explicit DGFGridFactory ( std::istream &input,
00353                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00354     : Base() 
00355     {
00356       input.clear();
00357       input.seekg( 0 );
00358       if( !input )
00359         DUNE_THROW(DGFException, "Error resetting input stream." );
00360       generate( input, comm );
00361     }
00362 
00363     explicit DGFGridFactory ( const std::string &filename, 
00364                               MPICommunicatorType comm = MPIHelper::getCommunicator()) 
00365     : DGFBaseFactory< ALUCubeGrid<2,dimw> >()
00366     {
00367       std::ifstream input( filename.c_str() );
00368       if( !input )
00369         DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
00370       if( !generate( input, comm, filename ) )
00371       {
00372         if( Base::fileExists( filename.c_str() ) )
00373           Base::grid_ = new Grid( filename );
00374         else
00375           DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
00376       }
00377     }
00378 
00379   protected:
00380     bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
00381     using Base::grid_;
00382     using Base::factory_;
00383     using Base::dgf_;
00384   };
00385 
00386 }
00387 
00388 #include "dgfalu.cc"
00389 
00390 #endif // #if defined ENABLE_ALUGRID
00391 
00392 #endif // #ifndef DUNE_DGFPARSERALU_HH