1 #ifndef DUNE_DGFPARSERALU_HH
2 #define DUNE_DGFPARSERALU_HH
20 template<
class Gr
idImp,
template<
class >
class IntersectionImp >
28 struct DGFGridInfo< ALUCubeGrid< 3, 3 > >
35 struct DGFGridInfo< ALUSimplexGrid< 3, 3 > >
42 struct DGFGridInfo< ALUSimplexGrid< 2, dimw > >
49 struct DGFGridInfo< ALUCubeGrid< 2, dimw > >
56 struct DGFGridInfo< Dune::ALUConformGrid< 2, dimw > >
62 template<
int dimg,
int dimw, ALUGr
idElementType eltype, ALUGr
idRefinementType refinementtype,
class Comm >
63 struct DGFGridInfo< Dune::ALUGrid< dimg, dimw, eltype, refinementtype, Comm > >
66 static double refineWeight () {
return ( refinementtype ==
conforming ) ? 0.5 : 1.0/(std::pow( 2.0,
double(dimg))); }
79 typedef MPIHelper::MPICommunicator MPICommunicatorType;
80 typedef typename Grid::template Codim<0>::Entity
Element;
81 typedef typename Grid::template Codim<dimension>::Entity Vertex;
90 explicit DGFBaseFactory ( MPICommunicatorType comm )
92 dgf_( rank(comm), size(comm) )
101 template<
class Intersection >
102 bool wasInserted (
const Intersection &intersection )
const
104 return factory_.wasInserted( intersection );
107 template <
class GG,
template <
class >
class II >
108 int boundaryId (
const Intersection< GG, II > & intersection )
const
113 const int face = intersection.indexInInside();
115 const GenericReferenceElement< double, dimension > & refElem =
116 GenericReferenceElements< double, dimension >::general( entity.type() );
117 int corners = refElem.size( face, 1, dimension );
118 std :: vector< unsigned int > bound( corners );
119 for(
int i=0; i < corners; ++i )
121 const int k = refElem.subEntity( face, 1, i, dimension );
122 bound[ i ] = factory_.insertionIndex( *entity.template subEntity< dimension >( k ) );
125 DuneGridFormatParser::facemap_t::key_type key( bound,
false );
126 const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
127 if( pos != dgf_.facemap.end() )
128 return dgf_.facemap.find( key )->second.first;
130 return (intersection.boundary() ? 1 : 0);
133 template <
class GG,
template <
class >
class II >
135 boundaryParameter (
const Intersection< GG, II > & intersection )
const
140 const int face = intersection.indexInInside();
142 const GenericReferenceElement< double, dimension > & refElem =
143 GenericReferenceElements< double, dimension >::general( entity.type() );
144 int corners = refElem.size( face, 1, dimension );
145 std :: vector< unsigned int > bound( corners );
146 for(
int i=0; i < corners; ++i )
148 const int k = refElem.subEntity( face, 1, i, dimension );
149 bound[ i ] = factory_.insertionIndex( *entity.template subEntity< dimension >( k ) );
152 DuneGridFormatParser::facemap_t::key_type key( bound,
false );
153 const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
154 if( pos != dgf_.facemap.end() )
155 return dgf_.facemap.find( key )->second.second;
160 template<
int codim >
161 int numParameters ()
const
164 return dgf_.nofelparams;
165 else if( codim == dimension )
166 return dgf_.nofvtxparams;
172 bool haveBoundaryParameters ()
const
174 return dgf_.haveBndParameters;
177 std::vector< double > ¶meter (
const Element &element )
179 if( numParameters< 0 >() <= 0 )
181 DUNE_THROW( InvalidStateException,
182 "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
184 return dgf_.elParams[ factory_.insertionIndex( element ) ];
187 std::vector< double > ¶meter (
const Vertex &
vertex )
189 if( numParameters< dimension >() <= 0 )
191 DUNE_THROW( InvalidStateException,
192 "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
194 return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
200 MPICommunicatorType communicator,
201 const std::string &filename );
205 MPICommunicatorType communicator,
206 const std::string &filename );
208 static Grid* callDirectly(
const char* gridname,
210 const char *filename,
211 MPICommunicatorType communicator )
213 #if ALU3DGRID_PARALLEL
215 std :: stringstream tmps;
216 tmps << filename <<
"." << rank;
217 const std :: string &tmp = tmps.str();
220 if( fileExists( tmp.c_str() ) )
221 return new Grid( tmp.c_str(), communicator );
226 if( fileExists( filename ) )
227 return new Grid( filename , communicator );
231 DUNE_THROW( GridError,
"Unable to create " << gridname <<
" from '"
232 << filename <<
"'." );
235 dwarn <<
"WARNING: P[" << rank <<
"]: Creating empty grid!" << std::endl;
238 return new Grid( communicator );
240 static bool fileExists (
const char *fileName )
242 std :: ifstream testfile( fileName );
248 static int rank( MPICommunicatorType MPICOMM )
252 MPI_Comm_rank( MPICOMM, &rank );
256 static int size( MPICommunicatorType MPICOMM )
260 MPI_Comm_size( MPICOMM, &size );
265 GridFactory factory_;
266 DuneGridFormatParser dgf_;
270 struct DGFGridFactory< ALUSimplexGrid<3,3> > :
271 public DGFBaseFactory< ALUSimplexGrid<3,3> >
275 : DGFBaseFactory< ALUSimplexGrid<3,3> >( comm )
280 DUNE_THROW( DGFException,
"Error resetting input stream." );
281 generate( input, comm );
286 : DGFBaseFactory< ALUSimplexGrid<3,3> >( comm )
288 std::ifstream input( filename.c_str() );
290 bool fileFound = input.is_open() ;
292 fileFound = generate( input, comm, filename );
295 grid_ = callDirectly(
"ALUSimplexGrid< 3 , 3 >", this->rank( comm ), filename.c_str(), comm );
299 bool generate( std::istream &file,
MPICommunicatorType comm,
const std::string &filename =
"" );
303 struct DGFGridFactory< ALUCubeGrid<3,3> > :
304 public DGFBaseFactory< ALUCubeGrid<3,3> >
308 : DGFBaseFactory< ALUCubeGrid<3,3> >( comm )
313 DUNE_THROW( DGFException,
"Error resetting input stream." );
314 generate( input, comm );
319 : DGFBaseFactory< ALUCubeGrid<3,3> >( comm )
321 std::ifstream input( filename.c_str() );
322 bool fileFound = input.is_open() ;
324 fileFound = generate( input, comm, filename );
327 grid_ = callDirectly(
"ALUCubeGrid< 3 , 3 >", this->rank( comm ), filename.c_str(), comm );
331 bool generate( std::istream &file,
MPICommunicatorType comm,
const std::string &filename =
"" );
334 template < ALUGr
idElementType eltype, ALUGr
idRefinementType refinementtype,
class Comm >
335 struct DGFGridFactory< ALUGrid<3,3, eltype, refinementtype, Comm > > :
336 public DGFBaseFactory< ALUGrid<3,3, eltype, refinementtype, Comm > >
338 typedef ALUGrid<3,3, eltype, refinementtype, Comm > DGFGridType;
339 typedef DGFBaseFactory< DGFGridType > BaseType;
342 using BaseType :: grid_;
343 using BaseType :: callDirectly;
352 DUNE_THROW( DGFException,
"Error resetting input stream." );
353 generate( input, comm );
360 std::ifstream input( filename.c_str() );
361 bool fileFound = input.is_open() ;
363 fileFound = generate( input, comm, filename );
366 grid_ = callDirectly(
"ALUGrid< 3, 3, eltype, ref, comm >", this->rank( comm ), filename.c_str(), comm );
370 bool generate( std::istream &file,
MPICommunicatorType comm,
const std::string &filename =
"" );
374 struct DGFGridFactory< ALUConformGrid<2,dimw> > :
375 public DGFBaseFactory< ALUConformGrid<2,dimw> >
377 typedef DGFBaseFactory< ALUConformGrid<2,dimw> > Base;
379 typedef typename Base::Grid
Grid;
381 typedef typename Base::GridFactory GridFactory;
390 DUNE_THROW( DGFException,
"Error resetting input stream." );
391 generate( input, comm );
396 : DGFBaseFactory< ALUConformGrid<2,dimw> >()
398 std::ifstream input( filename.c_str() );
400 DUNE_THROW( DGFException,
"Macrofile '" << filename <<
"' not found." );
401 if( !generate( input, comm, filename ) )
403 if( Base::fileExists( filename.c_str() ) )
404 Base::grid_ =
new Grid( filename );
406 DUNE_THROW( GridError,
"Unable to create a 2d ALUGrid from '" << filename <<
"'." );
410 bool generate( std::istream &file,
MPICommunicatorType comm,
const std::string &filename =
"" );
412 using Base::factory_;
417 struct DGFGridFactory< ALUSimplexGrid<2,dimw> > :
418 public DGFBaseFactory< ALUSimplexGrid<2,dimw> >
420 typedef DGFBaseFactory< ALUSimplexGrid<2,dimw> > Base;
422 typedef typename Base::Grid
Grid;
424 typedef typename Base::GridFactory GridFactory;
433 DUNE_THROW(DGFException,
"Error resetting input stream." );
434 generate( input, comm );
439 : DGFBaseFactory< ALUSimplexGrid<2,dimw> >()
441 std::ifstream input( filename.c_str() );
443 DUNE_THROW( DGFException,
"Macrofile '" << filename <<
"' not found." );
444 if( !generate( input, comm, filename ) )
446 if( Base::fileExists( filename.c_str() ) )
447 Base::grid_ =
new Grid( filename );
449 DUNE_THROW( GridError,
"Unable to create a 2d ALUGrid from '" << filename <<
"'." );
454 bool generate( std::istream &file,
MPICommunicatorType comm,
const std::string &filename =
"" );
456 using Base::factory_;
461 struct DGFGridFactory< ALUCubeGrid<2,dimw> > :
462 public DGFBaseFactory< ALUCubeGrid<2,dimw> >
464 typedef DGFBaseFactory< ALUCubeGrid<2,dimw> > Base;
466 typedef typename Base::Grid
Grid;
468 typedef typename Base::GridFactory GridFactory;
477 DUNE_THROW(DGFException,
"Error resetting input stream." );
478 generate( input, comm );
483 : DGFBaseFactory< ALUCubeGrid<2,dimw> >()
485 std::ifstream input( filename.c_str() );
487 DUNE_THROW( DGFException,
"Macrofile '" << filename <<
"' not found." );
488 if( !generate( input, comm, filename ) )
490 if( Base::fileExists( filename.c_str() ) )
491 Base::grid_ =
new Grid( filename );
493 DUNE_THROW( GridError,
"Unable to create a 2d ALUGrid from '" << filename <<
"'." );
498 bool generate( std::istream &file,
MPICommunicatorType comm,
const std::string &filename =
"" );
500 using Base::factory_;
504 template <
int dimw, ALUGr
idElementType eltype, ALUGr
idRefinementType refinementtype,
class Comm >
505 struct DGFGridFactory< ALUGrid<2, dimw, eltype, refinementtype, Comm > > :
506 public DGFBaseFactory< ALUGrid< 2, dimw, eltype, refinementtype, Comm > >
508 typedef ALUGrid< 2, dimw, eltype, refinementtype, Comm > DGFGridType;
509 typedef DGFBaseFactory< DGFGridType > BaseType;
511 typedef typename BaseType::Grid
Grid;
513 typedef typename BaseType::GridFactory GridFactory;
522 DUNE_THROW(DGFException,
"Error resetting input stream." );
523 generate( input, comm );
530 std::ifstream input( filename.c_str() );
532 DUNE_THROW( DGFException,
"Macrofile '" << filename <<
"' not found." );
533 if( !generate( input, comm, filename ) )
535 if( BaseType::fileExists( filename.c_str() ) )
536 grid_ =
new Grid( filename );
538 DUNE_THROW( GridError,
"Unable to create a 2d ALUGrid from '" << filename <<
"'." );
543 bool generate( std::istream &file,
MPICommunicatorType comm,
const std::string &filename =
"" );
544 using BaseType::grid_;
545 using BaseType::factory_;
546 using BaseType::dgf_;
553 #endif // #if HAVE_ALUGRID
555 #endif // #ifndef DUNE_DGFPARSERALU_HH