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