dune-grid  2.2.1
dgfalu.hh
Go to the documentation of this file.
1 #ifndef DUNE_DGFPARSERALU_HH
2 #define DUNE_DGFPARSERALU_HH
3 
4 // only include if ALUGrid is used
5 #if HAVE_ALUGRID
6 
7 #include <dune/grid/alugrid.hh>
11 
13 
14 namespace Dune
15 {
16 
17  // forward declaration
18  // -------------------
19 
20  template< class GridImp, template< class > class IntersectionImp >
21  class Intersection;
22 
23 // DGFGridInfo (specialization for ALUGrid)
24 // ----------------------------------------
25 
27 template<>
28 struct DGFGridInfo< ALUCubeGrid< 3, 3 > >
29 {
30  static int refineStepsForHalf () { return 1; }
31  static double refineWeight () { return 0.125; }
32 };
33 
34 template<>
35 struct DGFGridInfo< ALUSimplexGrid< 3, 3 > >
36 {
37  static int refineStepsForHalf () { return 1; }
38  static double refineWeight () { return 0.125; }
39 };
40 
41 template<int dimw>
42 struct DGFGridInfo< ALUSimplexGrid< 2, dimw > >
43 {
44  static int refineStepsForHalf () { return 1; }
45  static double refineWeight () { return 0.25; }
46 };
47 
48 template<int dimw>
49 struct DGFGridInfo< ALUCubeGrid< 2, dimw > >
50 {
51  static int refineStepsForHalf () { return 1; }
52  static double refineWeight () { return 0.25; }
53 };
54 
55 template<int dimw>
56 struct DGFGridInfo< Dune::ALUConformGrid< 2, dimw > >
57 {
58  static int refineStepsForHalf () { return 2; }
59  static double refineWeight () { return 0.5; }
60 };
61 
62 template<int dimg, int dimw, ALUGridElementType eltype, ALUGridRefinementType refinementtype, class Comm >
63 struct DGFGridInfo< Dune::ALUGrid< dimg, dimw, eltype, refinementtype, Comm > >
64 {
65  static int refineStepsForHalf () { return ( refinementtype == conforming ) ? dimg : 1; }
66  static double refineWeight () { return ( refinementtype == conforming ) ? 0.5 : 1.0/(std::pow( 2.0, double(dimg))); }
67 };
70  // DGFGridFactory for AluGrid
71  // --------------------------
72 
73  // template< int dim, int dimworld > // for a first version
74  template < class G >
75  struct DGFBaseFactory
76  {
77  typedef G Grid;
78  const static int dimension = Grid::dimension;
79  typedef MPIHelper::MPICommunicator MPICommunicatorType;
80  typedef typename Grid::template Codim<0>::Entity Element;
81  typedef typename Grid::template Codim<dimension>::Entity Vertex;
82  typedef Dune::GridFactory<Grid> GridFactory;
83 
84  DGFBaseFactory ()
85  : factory_( ),
86  dgf_( 0, 1 )
87  {
88  }
89 
90  explicit DGFBaseFactory ( MPICommunicatorType comm )
91  : factory_(),
92  dgf_( rank(comm), size(comm) )
93  {
94  }
95 
96  Grid *grid () const
97  {
98  return grid_;
99  }
100 
101  template< class Intersection >
102  bool wasInserted ( const Intersection &intersection ) const
103  {
104  return factory_.wasInserted( intersection );
105  }
106 
107  template < class GG, template < class > class II >
108  int boundaryId ( const Intersection< GG, II > & intersection ) const
109  {
110  typedef Dune::Intersection< GG, II > Intersection;
111  typename Intersection::EntityPointer inside = intersection.inside();
112  const typename Intersection::Entity & entity = *inside;
113  const int face = intersection.indexInInside();
114 
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 )
120  {
121  const int k = refElem.subEntity( face, 1, i, dimension );
122  bound[ i ] = factory_.insertionIndex( *entity.template subEntity< dimension >( k ) );
123  }
124 
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;
129  else
130  return (intersection.boundary() ? 1 : 0);
131  }
132 
133  template < class GG, template < class > class II >
134  const typename DGFBoundaryParameter::type &
135  boundaryParameter ( const Intersection< GG, II > & intersection ) const
136  {
137  typedef Dune::Intersection< GG, II > Intersection;
138  typename Intersection::EntityPointer inside = intersection.inside();
139  const typename Intersection::Entity & entity = *inside;
140  const int face = intersection.indexInInside();
141 
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 )
147  {
148  const int k = refElem.subEntity( face, 1, i, dimension );
149  bound[ i ] = factory_.insertionIndex( *entity.template subEntity< dimension >( k ) );
150  }
151 
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;
156  else
158  }
159 
160  template< int codim >
161  int numParameters () const
162  {
163  if( codim == 0 )
164  return dgf_.nofelparams;
165  else if( codim == dimension )
166  return dgf_.nofvtxparams;
167  else
168  return 0;
169  }
170 
171  // return true if boundary parameters found
172  bool haveBoundaryParameters () const
173  {
174  return dgf_.haveBndParameters;
175  }
176 
177  std::vector< double > &parameter ( const Element &element )
178  {
179  if( numParameters< 0 >() <= 0 )
180  {
181  DUNE_THROW( InvalidStateException,
182  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
183  }
184  return dgf_.elParams[ factory_.insertionIndex( element ) ];
185  }
186 
187  std::vector< double > &parameter ( const Vertex &vertex )
188  {
189  if( numParameters< dimension >() <= 0 )
190  {
191  DUNE_THROW( InvalidStateException,
192  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
193  }
194  return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
195  }
196 
197  protected:
198  bool generateALUGrid( const ALUGridElementType eltype,
199  std::istream &file,
200  MPICommunicatorType communicator,
201  const std::string &filename );
202 
203  bool generateALU2dGrid( const ALUGridElementType eltype,
204  std::istream &file,
205  MPICommunicatorType communicator,
206  const std::string &filename );
207 
208  static Grid* callDirectly( const char* gridname,
209  const int rank,
210  const char *filename,
211  MPICommunicatorType communicator )
212  {
213  #if ALU3DGRID_PARALLEL
214  // in parallel runs add rank to filename
215  std :: stringstream tmps;
216  tmps << filename << "." << rank;
217  const std :: string &tmp = tmps.str();
218 
219  // if file exits then use it
220  if( fileExists( tmp.c_str() ) )
221  return new Grid( tmp.c_str(), communicator );
222  #endif
223  // for rank 0 we also check the normal file name
224  if( rank == 0 )
225  {
226  if( fileExists( filename ) )
227  return new Grid( filename , communicator );
228 
229  // only throw this exception on rank 0 because
230  // for the other ranks we can still create empty grids
231  DUNE_THROW( GridError, "Unable to create " << gridname << " from '"
232  << filename << "'." );
233  }
234  else
235  dwarn << "WARNING: P[" << rank << "]: Creating empty grid!" << std::endl;
236 
237  // return empty grid on all other processes
238  return new Grid( communicator );
239  }
240  static bool fileExists ( const char *fileName )
241  {
242  std :: ifstream testfile( fileName );
243  if( !testfile )
244  return false;
245  testfile.close();
246  return true;
247  }
248  static int rank( MPICommunicatorType MPICOMM )
249  {
250  int rank = 0;
251 #if HAVE_MPI
252  MPI_Comm_rank( MPICOMM, &rank );
253 #endif
254  return rank;
255  }
256  static int size( MPICommunicatorType MPICOMM )
257  {
258  int size = 1;
259 #if HAVE_MPI
260  MPI_Comm_size( MPICOMM, &size );
261 #endif
262  return size;
263  }
264  Grid *grid_;
265  GridFactory factory_;
266  DuneGridFormatParser dgf_;
267  };
268 
269  template <>
270  struct DGFGridFactory< ALUSimplexGrid<3,3> > :
271  public DGFBaseFactory< ALUSimplexGrid<3,3> >
272  {
273  explicit DGFGridFactory ( std::istream &input,
274  MPICommunicatorType comm = MPIHelper::getCommunicator() )
275  : DGFBaseFactory< ALUSimplexGrid<3,3> >( comm )
276  {
277  input.clear();
278  input.seekg( 0 );
279  if( !input )
280  DUNE_THROW( DGFException, "Error resetting input stream." );
281  generate( input, comm );
282  }
283 
284  explicit DGFGridFactory ( const std::string &filename,
285  MPICommunicatorType comm = MPIHelper::getCommunicator())
286  : DGFBaseFactory< ALUSimplexGrid<3,3> >( comm )
287  {
288  std::ifstream input( filename.c_str() );
289 
290  bool fileFound = input.is_open() ;
291  if( fileFound )
292  fileFound = generate( input, comm, filename );
293 
294  if( ! fileFound )
295  grid_ = callDirectly( "ALUSimplexGrid< 3 , 3 >", this->rank( comm ), filename.c_str(), comm );
296  }
297 
298  protected:
299  bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
300  };
301 
302  template <>
303  struct DGFGridFactory< ALUCubeGrid<3,3> > :
304  public DGFBaseFactory< ALUCubeGrid<3,3> >
305  {
306  explicit DGFGridFactory ( std::istream &input,
307  MPICommunicatorType comm = MPIHelper::getCommunicator() )
308  : DGFBaseFactory< ALUCubeGrid<3,3> >( comm )
309  {
310  input.clear();
311  input.seekg( 0 );
312  if( !input )
313  DUNE_THROW( DGFException, "Error resetting input stream." );
314  generate( input, comm );
315  }
316 
317  explicit DGFGridFactory ( const std::string &filename,
318  MPICommunicatorType comm = MPIHelper::getCommunicator())
319  : DGFBaseFactory< ALUCubeGrid<3,3> >( comm )
320  {
321  std::ifstream input( filename.c_str() );
322  bool fileFound = input.is_open() ;
323  if( fileFound )
324  fileFound = generate( input, comm, filename );
325 
326  if( ! fileFound )
327  grid_ = callDirectly( "ALUCubeGrid< 3 , 3 >", this->rank( comm ), filename.c_str(), comm );
328  }
329 
330  protected:
331  bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
332  };
333 
334  template < ALUGridElementType eltype, ALUGridRefinementType refinementtype, class Comm >
335  struct DGFGridFactory< ALUGrid<3,3, eltype, refinementtype, Comm > > :
336  public DGFBaseFactory< ALUGrid<3,3, eltype, refinementtype, Comm > >
337  {
338  typedef ALUGrid<3,3, eltype, refinementtype, Comm > DGFGridType;
339  typedef DGFBaseFactory< DGFGridType > BaseType;
340  typedef typename BaseType :: MPICommunicatorType MPICommunicatorType;
341  protected:
342  using BaseType :: grid_;
343  using BaseType :: callDirectly;
344  public:
345  explicit DGFGridFactory ( std::istream &input,
346  MPICommunicatorType comm = MPIHelper::getCommunicator() )
347  : BaseType( comm )
348  {
349  input.clear();
350  input.seekg( 0 );
351  if( !input )
352  DUNE_THROW( DGFException, "Error resetting input stream." );
353  generate( input, comm );
354  }
355 
356  explicit DGFGridFactory ( const std::string &filename,
357  MPICommunicatorType comm = MPIHelper::getCommunicator())
358  : BaseType( comm )
359  {
360  std::ifstream input( filename.c_str() );
361  bool fileFound = input.is_open() ;
362  if( fileFound )
363  fileFound = generate( input, comm, filename );
364 
365  if( ! fileFound )
366  grid_ = callDirectly( "ALUGrid< 3, 3, eltype, ref, comm >", this->rank( comm ), filename.c_str(), comm );
367  }
368 
369  protected:
370  bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
371  };
372 
373  template <int dimw>
374  struct DGFGridFactory< ALUConformGrid<2,dimw> > :
375  public DGFBaseFactory< ALUConformGrid<2,dimw> >
376  {
377  typedef DGFBaseFactory< ALUConformGrid<2,dimw> > Base;
378  typedef typename Base:: MPICommunicatorType MPICommunicatorType;
379  typedef typename Base::Grid Grid;
380  const static int dimension = Grid::dimension;
381  typedef typename Base::GridFactory GridFactory;
382 
383  explicit DGFGridFactory ( std::istream &input,
384  MPICommunicatorType comm = MPIHelper::getCommunicator() )
385  : Base()
386  {
387  input.clear();
388  input.seekg( 0 );
389  if( !input )
390  DUNE_THROW( DGFException, "Error resetting input stream." );
391  generate( input, comm );
392  }
393 
394  explicit DGFGridFactory ( const std::string &filename,
395  MPICommunicatorType comm = MPIHelper::getCommunicator())
396  : DGFBaseFactory< ALUConformGrid<2,dimw> >()
397  {
398  std::ifstream input( filename.c_str() );
399  if( !input )
400  DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
401  if( !generate( input, comm, filename ) )
402  {
403  if( Base::fileExists( filename.c_str() ) )
404  Base::grid_ = new Grid( filename );
405  else
406  DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
407  }
408  }
409  protected:
410  bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
411  using Base::grid_;
412  using Base::factory_;
413  using Base::dgf_;
414  };
415 
416  template <int dimw>
417  struct DGFGridFactory< ALUSimplexGrid<2,dimw> > :
418  public DGFBaseFactory< ALUSimplexGrid<2,dimw> >
419  {
420  typedef DGFBaseFactory< ALUSimplexGrid<2,dimw> > Base;
421  typedef typename Base::MPICommunicatorType MPICommunicatorType;
422  typedef typename Base::Grid Grid;
423  const static int dimension = Grid::dimension;
424  typedef typename Base::GridFactory GridFactory;
425 
426  explicit DGFGridFactory ( std::istream &input,
427  MPICommunicatorType comm = MPIHelper::getCommunicator() )
428  : Base()
429  {
430  input.clear();
431  input.seekg( 0 );
432  if( !input )
433  DUNE_THROW(DGFException, "Error resetting input stream." );
434  generate( input, comm );
435  }
436 
437  explicit DGFGridFactory ( const std::string &filename,
438  MPICommunicatorType comm = MPIHelper::getCommunicator())
439  : DGFBaseFactory< ALUSimplexGrid<2,dimw> >()
440  {
441  std::ifstream input( filename.c_str() );
442  if( !input )
443  DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
444  if( !generate( input, comm, filename ) )
445  {
446  if( Base::fileExists( filename.c_str() ) )
447  Base::grid_ = new Grid( filename );
448  else
449  DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
450  }
451  }
452 
453  protected:
454  bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
455  using Base::grid_;
456  using Base::factory_;
457  using Base::dgf_;
458  };
459 
460  template <int dimw>
461  struct DGFGridFactory< ALUCubeGrid<2,dimw> > :
462  public DGFBaseFactory< ALUCubeGrid<2,dimw> >
463  {
464  typedef DGFBaseFactory< ALUCubeGrid<2,dimw> > Base;
465  typedef typename Base::MPICommunicatorType MPICommunicatorType;
466  typedef typename Base::Grid Grid;
467  const static int dimension = Grid::dimension;
468  typedef typename Base::GridFactory GridFactory;
469 
470  explicit DGFGridFactory ( std::istream &input,
471  MPICommunicatorType comm = MPIHelper::getCommunicator() )
472  : Base()
473  {
474  input.clear();
475  input.seekg( 0 );
476  if( !input )
477  DUNE_THROW(DGFException, "Error resetting input stream." );
478  generate( input, comm );
479  }
480 
481  explicit DGFGridFactory ( const std::string &filename,
482  MPICommunicatorType comm = MPIHelper::getCommunicator())
483  : DGFBaseFactory< ALUCubeGrid<2,dimw> >()
484  {
485  std::ifstream input( filename.c_str() );
486  if( !input )
487  DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
488  if( !generate( input, comm, filename ) )
489  {
490  if( Base::fileExists( filename.c_str() ) )
491  Base::grid_ = new Grid( filename );
492  else
493  DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
494  }
495  }
496 
497  protected:
498  bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
499  using Base::grid_;
500  using Base::factory_;
501  using Base::dgf_;
502  };
503 
504  template < int dimw, ALUGridElementType eltype, ALUGridRefinementType refinementtype, class Comm >
505  struct DGFGridFactory< ALUGrid<2, dimw, eltype, refinementtype, Comm > > :
506  public DGFBaseFactory< ALUGrid< 2, dimw, eltype, refinementtype, Comm > >
507  {
508  typedef ALUGrid< 2, dimw, eltype, refinementtype, Comm > DGFGridType;
509  typedef DGFBaseFactory< DGFGridType > BaseType;
510  typedef typename BaseType :: MPICommunicatorType MPICommunicatorType;
511  typedef typename BaseType::Grid Grid;
512  const static int dimension = Grid::dimension;
513  typedef typename BaseType::GridFactory GridFactory;
514 
515  explicit DGFGridFactory ( std::istream &input,
516  MPICommunicatorType comm = MPIHelper::getCommunicator() )
517  : BaseType( comm )
518  {
519  input.clear();
520  input.seekg( 0 );
521  if( !input )
522  DUNE_THROW(DGFException, "Error resetting input stream." );
523  generate( input, comm );
524  }
525 
526  explicit DGFGridFactory ( const std::string &filename,
527  MPICommunicatorType comm = MPIHelper::getCommunicator())
528  : BaseType( comm )
529  {
530  std::ifstream input( filename.c_str() );
531  if( !input )
532  DUNE_THROW( DGFException, "Macrofile '" << filename << "' not found." );
533  if( !generate( input, comm, filename ) )
534  {
535  if( BaseType::fileExists( filename.c_str() ) )
536  grid_ = new Grid( filename );
537  else
538  DUNE_THROW( GridError, "Unable to create a 2d ALUGrid from '" << filename << "'." );
539  }
540  }
541 
542  protected:
543  bool generate( std::istream &file, MPICommunicatorType comm, const std::string &filename = "" );
544  using BaseType::grid_;
545  using BaseType::factory_;
546  using BaseType::dgf_;
547  };
548 
549 }
550 
551 #include "dgfalu.cc"
552 
553 #endif // #if HAVE_ALUGRID
554 
555 #endif // #ifndef DUNE_DGFPARSERALU_HH