dune-alugrid  2.6-git
gridfactory.hh
Go to the documentation of this file.
1 #ifndef DUNE_ALU3DGRID_FACTORY_HH
2 #define DUNE_ALU3DGRID_FACTORY_HH
3 
4 #include <algorithm>
5 #include <array>
6 #include <map>
7 #include <memory>
8 #include <vector>
9 
10 #include <dune/common/shared_ptr.hh>
11 #include <dune/common/parallel/mpihelper.hh>
12 #include <dune/common/version.hh>
13 
14 #include <dune/geometry/referenceelements.hh>
15 
16 #include <dune/grid/common/gridfactory.hh>
17 #include <dune/grid/common/boundaryprojection.hh>
18 
21 
23 
24 namespace Dune
25 {
27  template< class ALUGrid >
29  : public GridFactoryInterface< ALUGrid >
30  {
32  typedef GridFactoryInterface< ALUGrid > BaseType;
33 
34  public:
35  typedef ALUGrid Grid;
36 
37  typedef typename Grid::ctype ctype;
38 
39  static const ALU3dGridElementType elementType = Grid::elementType;
40 
41  static const unsigned int dimension = Grid::dimension;
42  static const unsigned int dimensionworld = Grid::dimensionworld;
43 
44  typedef typename Grid::MPICommunicatorType MPICommunicatorType;
45 
47  typedef DuneBoundaryProjection< dimensionworld > DuneBoundaryProjectionType;
48 
49  template< int codim >
50  struct Codim
51  {
52  typedef typename Grid::template Codim< codim >::Entity Entity;
53 #if !DUNE_VERSION_NEWER(DUNE_GRID,2,5)
54  typedef typename Grid::template Codim< codim >::EntityPointer EntityPointer;
55 #endif // #if !DUNE_VERSION_NEWER(DUNE_GRID,2,5)
56  };
57 
58  typedef unsigned int VertexId;
59  typedef unsigned int GlobalIdType;
60 
62 
67 
68  private:
69  static_assert ( (elementType == tetra || elementType == hexa),
70  "ALU3dGridFactory supports only grids containing "
71  "tetrahedrons or hexahedrons exclusively." );
72 
73  typedef Dune::BoundarySegmentWrapper< dimension, dimensionworld > BoundarySegmentWrapperType;
74 
75  static const unsigned int numCorners = EntityCount< elementType >::numVertices;
76  static const unsigned int numFaces = EntityCount< elementType >::numFaces;
77  static const unsigned int numFaceCorners = EntityCount< elementType >::numVerticesPerFace;
78 
81 
82  // type of vertex coordinates put into the factory
83  typedef FieldVector< ctype, dimensionworld > VertexInputType;
85 
86  // type of vertex coordinates stored inside the factory
87  typedef FieldVector< ctype, 3 > VertexType;
88 
89  typedef std::vector< unsigned int > ElementType;
90  typedef std::array< unsigned int, numFaceCorners > FaceType;
91 
92  struct FaceLess;
93 
94  typedef std::vector< std::pair< VertexType, GlobalIdType > > VertexVector;
95  typedef std::vector< ElementType > ElementVector;
96  typedef std::pair< FaceType, int > BndPair ;
97  typedef std::map < FaceType, int > BoundaryIdMap;
98  typedef std::vector< std::pair< BndPair, BndPair > > PeriodicBoundaryVector;
99  typedef std::pair< unsigned int, int > SubEntity;
100  typedef std::map< FaceType, SubEntity, FaceLess > FaceMap;
101 
102  typedef std::map< FaceType, const DuneBoundaryProjectionType* > BoundaryProjectionMap;
103  typedef std::vector< const DuneBoundaryProjectionType* > BoundaryProjectionVector;
104 
105  typedef std::vector< Transformation > FaceTransformationVector;
106 
107  static void copy ( const std::initializer_list< unsigned int > &vertices, FaceType &faceId )
108  {
109  std::copy_n( vertices.begin(), faceId.size(), faceId.begin() );
110  }
111 
112  static FaceType makeFace ( const std::vector< unsigned int > &vertices )
113  {
114  if( vertices.size() != (dimension == 2 ? 2 : numFaceCorners) )
115  DUNE_THROW( GridError, "Wrong number of face vertices passed: " << vertices.size() << "." );
116 
117  FaceType faceId;
118  if( dimension == 2 )
119  {
120  if( elementType == tetra )
121  copy( { 0, vertices[ 1 ]+1, vertices[ 0 ]+1 }, faceId );
122  else if( elementType == hexa )
123  copy( { 2*vertices[ 0 ], 2*vertices[ 1 ], 2*vertices[ 0 ]+1, 2*vertices[ 1 ]+1 }, faceId );
124  }
125  else if( dimension == 3 )
126  std::copy_n( vertices.begin(), numFaceCorners, faceId.begin() );
127  return faceId;
128  }
129 
130  static BndPair makeBndPair ( const FaceType &face, int id )
131  {
132  BndPair bndPair;
133  for( unsigned int i = 0; i < numFaceCorners; ++i )
134  {
135  const unsigned int j = FaceTopologyMappingType::dune2aluVertex( i );
136  bndPair.first[ j ] = face[ i ];
137  }
138  bndPair.second = 1;
139  return bndPair;
140  }
141 
142  private:
143  // return grid object
144  virtual Grid* createGridObj( BoundaryProjectionVector* bndProjections, const std::string& name ) const
145  {
146  return new Grid( communicator_, globalProjection_, bndProjections , name, realGrid_ );
147  }
148 
149  protected:
151  explicit ALU3dGridFactory ( const bool verbose, const MPICommunicatorType &communicator );
152 
153  public:
155  explicit ALU3dGridFactory ( const MPICommunicatorType &communicator = Grid::defaultCommunicator(),
156  bool removeGeneratedFile = true );
157 
159  explicit ALU3dGridFactory ( const std::string &filename,
160  const MPICommunicatorType &communicator = Grid::defaultCommunicator() );
161 
163  virtual ~ALU3dGridFactory ();
164 
169  virtual void insertVertex ( const VertexInputType &pos );
170 
176  void insertVertex ( const VertexInputType &pos, const VertexId globalId );
177 
186  virtual void
187  insertElement ( const GeometryType &geometry,
188  const std::vector< VertexId > &vertices );
189 
201  virtual void
202  insertBoundary ( const GeometryType &geometry, const std::vector< VertexId > &faceVertices, int boundaryId = 1 );
203 
204 
212  virtual void insertBoundary ( int element, int face, int boundaryId = 1 );
213 
214  // for testing parallel GridFactory
215  void insertProcessBorder ( int element, int face )
216  {
218  }
219 
228  virtual void
229  insertBoundaryProjection ( const GeometryType &type,
230  const std::vector< VertexId > &vertices,
231  const DuneBoundaryProjectionType *projection );
232 
237  virtual void
238  insertBoundarySegment ( const std::vector< VertexId >& vertices ) ;
239 
240  virtual void
241  insertProcessBorder ( const std::vector< VertexId >& vertices );
242 
248  virtual void
249  insertBoundarySegment ( const std::vector< VertexId >& vertices,
250  const shared_ptr<BoundarySegment<dimension,dimensionworld> >& boundarySegment ) ;
251 
256  virtual void insertBoundaryProjection ( const DuneBoundaryProjectionType& bndProjection );
257 
267  void insertFaceTransformation ( const WorldMatrix &matrix, const WorldVector &shift );
268 
273  Grid *createGrid ();
274 
275  Grid *createGrid ( const bool addMissingBoundaries, const std::string dgfName = "" );
276 
277  Grid *createGrid ( const bool addMissingBoundaries, bool temporary, const std::string dgfName = "" );
278 
279  virtual unsigned int
280  insertionIndex ( const typename Codim< 0 >::Entity &entity ) const
281  {
282  alugrid_assert( Grid::getRealImplementation( entity ).getIndex() < int(ordering_.size()) );
283  return ordering_[ Grid::getRealImplementation( entity ).getIndex() ];
284  }
285 
286  virtual unsigned int
287  insertionIndex ( const typename Codim< dimension >::Entity &entity ) const
288  {
289  if(dimension == 2 && elementType == hexa )
290  // for quadrilaterals we simply half the number, see gridfactory.cc doInsertVertex
291  return Grid::getRealImplementation( entity ).getIndex()/2;
292  else if ( dimension == 2 && elementType == tetra )
293  // for triangles we have to substract 1, see gridfactory.cc doInsertVertex
294  return Grid::getRealImplementation( entity ).getIndex() - 1;
295  else // dimension 3
296  return Grid::getRealImplementation( entity ).getIndex();
297  }
298 
299  virtual unsigned int insertionIndex ( const typename Grid::LevelIntersection &intersection ) const
300  {
301  return boundaryInsertionIndex( intersection.inside(), intersection.indexInInside() );
302  }
303 
304  virtual unsigned int insertionIndex ( const typename Grid::LeafIntersection &intersection ) const
305  {
306  return boundaryInsertionIndex( intersection.inside(), intersection.indexInInside() );
307  }
308 
309  virtual bool wasInserted ( const typename Grid::LevelIntersection &intersection ) const
310  {
311  return intersection.boundary() && (insertionIndex(intersection) < std::numeric_limits<unsigned int>::max());
312  }
313 
314  virtual bool wasInserted ( const typename Grid::LeafIntersection &intersection ) const
315  {
316  return intersection.boundary() && (insertionIndex(intersection) < std::numeric_limits<unsigned int>::max());
317  }
318 
319  const std::vector<unsigned int>& ordering () const { return ordering_; }
320 
321  private:
322  unsigned int boundaryInsertionIndex ( const typename Codim< 0 >::Entity &entity, int face ) const
323  {
324  const auto& refElem = Dune::ReferenceElements< double, dimension >::general( entity.type() );
325  const int vxSize = refElem.size( face, 1, dimension );
326  std::vector< unsigned int > vertices( vxSize );
327  for( int i = 0; i < vxSize; ++i )
328  vertices[ i ] = insertionIndex( entity.template subEntity< dimension >( refElem.subEntity( face, 1, i, dimension ) ) );
329 
330  FaceType faceId = makeFace( vertices );
331  std::sort( faceId.begin(), faceId.end() );
332 
333  const auto pos = insertionOrder_.find( faceId );
334  return (pos != insertionOrder_.end() ? pos->second : std::numeric_limits< unsigned int >::max());
335  }
336 
337  void doInsertVertex ( const VertexInputType &pos, const GlobalIdType globalId );
338  void doInsertBoundary ( int element, int face, int boundaryId );
339 
340  GlobalIdType globalId ( const VertexId &id ) const
341  {
342  alugrid_assert ( id < vertices_.size() );
343  return vertices_[ id ].second;
344  }
345 
346  const VertexType &position ( const VertexId &id ) const
347  {
348  alugrid_assert ( id < vertices_.size() );
349  return vertices_[ id ].first;
350  }
351 
352  const VertexInputType inputPosition ( const VertexId &id ) const
353  {
354  alugrid_assert ( id < vertices_.size() );
355  VertexType vertex = vertices_[ id ].first;
356  VertexInputType iVtx(0.);
357  for(unsigned i = 0 ; i < dimensionworld ; ++i)
358  iVtx[i] = vertex[i];
359  return iVtx;
360  }
361 
362  void assertGeometryType( const GeometryType &geometry );
363  static void generateFace ( const ElementType &element, const int f, FaceType &face );
364  void generateFace ( const SubEntity &subEntity, FaceType &face ) const;
365  void correctElementOrientation ();
366  bool identifyFaces ( const Transformation &transformation, const FaceType &key1, const FaceType &key2, const int defaultId );
367  void searchPeriodicNeighbor ( FaceMap &faceMap, const typename FaceMap::iterator &pos, const int defaultId );
368  void reinsertBoundary ( const FaceMap &faceMap, const typename FaceMap::const_iterator &pos, const int id );
369  void recreateBoundaryIds ( const int defaultId = 1 );
370 
371  // sort elements according to hilbert space filling curve (if Zoltan is available)
372  void sortElements( const VertexVector& vertices, const ElementVector& elements, std::vector< unsigned int >& ordering );
373 
374  int rank_;
375 
376  VertexVector vertices_;
377  ElementVector elements_;
378  BoundaryIdMap boundaryIds_,insertionOrder_;
379  PeriodicBoundaryVector periodicBoundaries_;
380  const DuneBoundaryProjectionType* globalProjection_ ;
381  BoundaryProjectionMap boundaryProjections_;
382  FaceTransformationVector faceTransformations_;
383  unsigned int numFacesInserted_;
384  bool realGrid_;
385  const bool allowGridGeneration_;
386  bool foundGlobalIndex_ ;
387 
388  MPICommunicatorType communicator_;
389 
390  typename SpaceFillingCurveOrderingType :: CurveType curveType_;
391  std::vector< unsigned int > ordering_;
392  };
393 
394 
395 
396  template< class ALUGrid >
398  : public std::binary_function< FaceType, FaceType, bool >
399  {
400  bool operator() ( const FaceType &a, const FaceType &b ) const
401  {
402  for( unsigned int i = 0; i < numFaceCorners; ++i )
403  {
404  if( a[ i ] != b[ i ] )
405  return (a[ i ] < b[ i ]);
406  }
407  return false;
408  }
409  };
410 
411 
412  template< class ALUGrid >
413  inline void ALU3dGridFactory< ALUGrid >
414  ::assertGeometryType( const GeometryType &geometry )
415  {
416  if( elementType == tetra )
417  {
418  if( !geometry.isSimplex() )
419  DUNE_THROW( GridError, "Only simplex geometries can be inserted into "
420  "ALUGrid< 3, 3, simplex, refrule >." );
421  }
422  else
423  {
424  if( !geometry.isCube() )
425  DUNE_THROW( GridError, "Only cube geometries can be inserted into "
426  "ALUGrid< 3, 3, cube, refrule >." );
427  }
428  }
429 
433  template<int dim, int dimw, ALUGridElementType eltype, ALUGridRefinementType refinementtype , class Comm >
434  class GridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
435  : public ALU3dGridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
436  {
437  typedef GridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > > ThisType;
439 
440  public:
441  typedef typename BaseType::Grid Grid;
442 
444 
446  explicit GridFactory ( const MPICommunicatorType &communicator = Grid::defaultCommunicator() )
447  : BaseType( communicator )
448  {}
449 
451  explicit GridFactory ( const std::string &filename,
452  const MPICommunicatorType &communicator = Grid::defaultCommunicator() )
453  : BaseType( filename, communicator )
454  {}
455  };
456 
457  template< class Grid >
459 
460  // Specialization of the ReferenceGridFactory for ALUGrid
461  template<int dim, int dimw, ALUGridElementType eltype, ALUGridRefinementType refinementtype , class Comm >
462  class ReferenceGridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
463  : public ALU3dGridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
464  {
467 
468  public:
469  typedef typename BaseType::Grid Grid;
470 
472 
475  : BaseType(false, Grid::defaultCommunicator() )
476  {}
477  };
478 
479 
480 
481  // Implementation of ALU3dGridFactory
482  // ----------------------------------
483 
484  template< class ALUGrid >
485  inline
488  bool removeGeneratedFile )
489  : rank_( ALU3dGridCommunications< ALUGrid::dimension, ALUGrid::dimensionworld, elementType, MPICommunicatorType >::getRank( communicator ) ),
490  globalProjection_ ( 0 ),
491  numFacesInserted_ ( 0 ),
492  realGrid_( true ),
493  allowGridGeneration_( rank_ == 0 ),
494  foundGlobalIndex_( false ),
495  communicator_( communicator ),
496  curveType_( SpaceFillingCurveOrderingType :: DefaultCurve )
497  {}
498 
499  template< class ALUGrid >
500  inline
502  :: ALU3dGridFactory ( const std::string &filename,
503  const MPICommunicatorType &communicator )
504  : rank_( ALU3dGridCommunications< ALUGrid::dimension, ALUGrid::dimensionworld, elementType, MPICommunicatorType >::getRank( communicator ) ),
505  globalProjection_ ( 0 ),
506  numFacesInserted_ ( 0 ),
507  realGrid_( true ),
508  allowGridGeneration_( rank_ == 0 ),
509  foundGlobalIndex_( false ),
510  communicator_( communicator ),
511  curveType_( SpaceFillingCurveOrderingType :: DefaultCurve )
512  {}
513 
514  template< class ALUGrid >
515  inline
517  :: ALU3dGridFactory ( const bool realGrid,
518  const MPICommunicatorType &communicator )
519  : rank_( ALU3dGridCommunications< ALUGrid::dimension, ALUGrid::dimensionworld, elementType, MPICommunicatorType >::getRank( communicator ) ),
520  globalProjection_ ( 0 ),
521  numFacesInserted_ ( 0 ),
522  realGrid_( realGrid ),
523  allowGridGeneration_( true ),
524  foundGlobalIndex_( false ),
525  communicator_( communicator ),
526  curveType_( SpaceFillingCurveOrderingType :: DefaultCurve )
527  {}
528 
529  template< class ALUGrid >
531  insertBoundarySegment ( const std::vector< unsigned int >& vertices )
532  {
533  FaceType faceId = makeFace( vertices );
534 
535  boundaryIds_.insert( makeBndPair( faceId, 1 ) );
536 
537  std::sort( faceId.begin(), faceId.end() );
538  if( boundaryProjections_.find( faceId ) != boundaryProjections_.end() )
539  DUNE_THROW( GridError, "Only one boundary projection can be attached to a face." );
540 
541  boundaryProjections_[ faceId ] = nullptr;
542 
543  insertionOrder_.insert( std::make_pair( faceId, insertionOrder_.size() ) );
544  }
545 
546  template< class ALUGrid >
548  insertProcessBorder ( const std::vector< unsigned int >& vertices )
549  {
550  FaceType faceId = makeFace( vertices );
551 
552  boundaryIds_.insert( makeBndPair( faceId, ALU3DSPACE ProcessorBoundary_t ) );
553 
554  std::sort( faceId.begin(), faceId.end() );
555  boundaryProjections_[ faceId ] = nullptr;
556  }
557 
558  template< class ALUGrid >
560  insertBoundarySegment ( const std::vector< unsigned int >& vertices,
561  const shared_ptr<BoundarySegment<dimension,dimensionworld> >& boundarySegment )
562  {
563  const std::size_t numVx = vertices.size();
564 
565  const unsigned int faceTopoId = (elementType == tetra) ?
566  Dune::Impl::SimplexTopology< dimension-1 >::type::id :
567  Dune::Impl::CubeTopology< dimension-1 >::type::id ;
568  GeometryType type( faceTopoId, dimension-1 );
569 
570  // we need double here because of the structure of BoundarySegment
571  // and BoundarySegmentWrapper which have double as coordinate type
572  typedef FieldVector< double, dimensionworld > CoordType;
573  std::vector< CoordType > coords( numVx );
574  for( std::size_t i = 0; i < numVx; ++i )
575  {
576  // adapt vertex index for 2d grids
577  const std::size_t vtx = (dimension == 2 ? (elementType == tetra ? vertices[ i ] + 1 : 2 * vertices[ i ]) : vertices[ i ]);
578 
579  // if this assertions is thrown vertices were not inserted at first
580  alugrid_assert ( vertices_.size() > vtx );
581 
582  // get global coordinate and copy it
583  std::copy_n( position( vtx ).begin(), dimensionworld, coords[ i ].begin() );
584  }
585 
586  std::unique_ptr< BoundarySegmentWrapperType > prj( new BoundarySegmentWrapperType( type, coords, boundarySegment ) );
587 
588  // consistency check
589  for( std::size_t i = 0; i < numVx; ++i )
590  {
591  CoordType global = (*prj)( coords [ i ] );
592  if( (global - coords[ i ]).two_norm() > 1e-6 )
593  DUNE_THROW( GridError, "BoundarySegment does not map face vertices to face vertices." );
594  }
595 
596  FaceType faceId = makeFace( vertices );
597 
598  boundaryIds_.insert( makeBndPair( faceId, 1 ) );
599 
600  std::sort( faceId.begin(), faceId.end() );
601  if( boundaryProjections_.find( faceId ) != boundaryProjections_.end() )
602  DUNE_THROW( GridError, "Only one boundary projection can be attached to a face." );
603 
604  boundaryProjections_[ faceId ] = prj.release();
605 
606  insertionOrder_.insert( std::make_pair( faceId, insertionOrder_.size() ) );
607  }
608 
609 
610  template< class ALUGrid >
611  inline void ALU3dGridFactory< ALUGrid >
612  ::generateFace ( const SubEntity &subEntity, FaceType &face ) const
613  {
614  generateFace( elements_[ subEntity.first ], subEntity.second, face );
615  }
616 
617 } // end namespace Dune
618 
619 #if COMPILE_ALUGRID_INLINE
620  #include "gridfactory.cc"
621 #endif
622 #endif
#define ALU3DSPACE
Definition: alu3dinclude.hh:24
Provides base classes for ALUGrid.
#define alugrid_assert(EX)
Definition: alugrid_assert.hh:20
Definition: alu3dinclude.hh:50
static const int ProcessorBoundary_t
Definition: alu3dinclude.hh:54
Definition: alu3dinclude.hh:80
ALU3dGridElementType
Definition: topology.hh:12
@ hexa
Definition: topology.hh:12
@ tetra
Definition: topology.hh:12
unstructured parallel implementation of the DUNE grid interface
Definition: alugrid.hh:31
BaseType::ctype ctype
Definition: alugrid.hh:48
@ dimension
Definition: alugrid.hh:47
@ dimensionworld
Definition: alugrid.hh:47
Factory class for ALUGrids.
Definition: gridfactory.hh:30
unsigned int VertexId
Definition: gridfactory.hh:58
static const unsigned int dimensionworld
Definition: gridfactory.hh:42
virtual void insertVertex(const VertexInputType &pos)
insert a vertex into the coarse grid
Definition: gridfactory.cc:40
Grid::MPICommunicatorType MPICommunicatorType
Definition: gridfactory.hh:44
virtual unsigned int insertionIndex(const typename Codim< 0 >::Entity &entity) const
Definition: gridfactory.hh:280
virtual bool wasInserted(const typename Grid::LeafIntersection &intersection) const
Definition: gridfactory.hh:314
virtual void insertBoundarySegment(const std::vector< VertexId > &vertices, const shared_ptr< BoundarySegment< dimension, dimensionworld > > &boundarySegment)
insert a shaped boundary segment into the macro grid
Definition: gridfactory.hh:560
void insertFaceTransformation(const WorldMatrix &matrix, const WorldVector &shift)
add a face transformation (for periodic identification)
Definition: gridfactory.cc:236
virtual unsigned int insertionIndex(const typename Codim< dimension >::Entity &entity) const
Definition: gridfactory.hh:287
ALU3dGridFactory(const MPICommunicatorType &communicator=Grid::defaultCommunicator(), bool removeGeneratedFile=true)
default constructor
Definition: gridfactory.hh:487
ALUGrid Grid
Definition: gridfactory.hh:35
const std::vector< unsigned int > & ordering() const
Definition: gridfactory.hh:319
Transformation::WorldMatrix WorldMatrix
type of matrix from world coordinates to world coordinates
Definition: gridfactory.hh:66
virtual void insertBoundarySegment(const std::vector< VertexId > &vertices)
insert a boundary segment into the macro grid
Definition: gridfactory.hh:531
ALUGridTransformation< ctype, dimensionworld > Transformation
Definition: gridfactory.hh:61
virtual void insertBoundary(const GeometryType &geometry, const std::vector< VertexId > &faceVertices, int boundaryId=1)
insert a boundary element into the coarse grid
Definition: gridfactory.cc:154
void insertProcessBorder(int element, int face)
Definition: gridfactory.hh:215
DuneBoundaryProjection< dimensionworld > DuneBoundaryProjectionType
type of boundary projection class
Definition: gridfactory.hh:47
Transformation::WorldVector WorldVector
type of vector for world coordinates
Definition: gridfactory.hh:64
virtual void insertProcessBorder(const std::vector< VertexId > &vertices)
Definition: gridfactory.hh:548
ALU3dGridFactory(const std::string &filename, const MPICommunicatorType &communicator=Grid::defaultCommunicator())
constructor taking filename for temporary outfile
Definition: gridfactory.hh:502
static const unsigned int dimension
Definition: gridfactory.hh:41
virtual void insertElement(const GeometryType &geometry, const std::vector< VertexId > &vertices)
insert an element into the coarse grid
Definition: gridfactory.cc:104
Grid * createGrid()
finalize the grid creation and hand over the grid
Definition: gridfactory.cc:335
virtual void insertBoundaryProjection(const GeometryType &type, const std::vector< VertexId > &vertices, const DuneBoundaryProjectionType *projection)
insert a boundary projection into the macro grid
Definition: gridfactory.cc:216
virtual unsigned int insertionIndex(const typename Grid::LevelIntersection &intersection) const
Definition: gridfactory.hh:299
unsigned int GlobalIdType
Definition: gridfactory.hh:59
virtual ~ALU3dGridFactory()
Destructor.
Definition: gridfactory.cc:34
virtual bool wasInserted(const typename Grid::LevelIntersection &intersection) const
Definition: gridfactory.hh:309
Grid::ctype ctype
Definition: gridfactory.hh:37
ALU3dGridFactory(const bool verbose, const MPICommunicatorType &communicator)
constructor taking verbose flag
Definition: gridfactory.hh:517
static const ALU3dGridElementType elementType
Definition: gridfactory.hh:39
virtual unsigned int insertionIndex(const typename Grid::LeafIntersection &intersection) const
Definition: gridfactory.hh:304
Definition: alugrid/3d/grid.hh:135
Definition: gridfactory.hh:51
Grid::template Codim< codim >::EntityPointer EntityPointer
Definition: gridfactory.hh:54
Grid::template Codim< codim >::Entity Entity
Definition: gridfactory.hh:52
Definition: gridfactory.hh:399
Specialization of the generic GridFactory for ALUGrid.
Definition: gridfactory.hh:436
GridFactory(const MPICommunicatorType &communicator=Grid::defaultCommunicator())
Default constructor.
Definition: gridfactory.hh:446
GridFactory(const std::string &filename, const MPICommunicatorType &communicator=Grid::defaultCommunicator())
constructor taking filename
Definition: gridfactory.hh:451
BaseType::MPICommunicatorType MPICommunicatorType
Definition: gridfactory.hh:443
Definition: gridfactory.hh:458
BaseType::MPICommunicatorType MPICommunicatorType
Definition: gridfactory.hh:471
ReferenceGridFactory()
Default constructor.
Definition: gridfactory.hh:474
Definition: topology.hh:15
Definition: topology.hh:40
Definition: topology.hh:151
static int dune2aluVertex(int index)
Maps vertex index from Dune onto ALU3dGrid reference face.
Definition: topology.hh:310
Definition: hsfc.hh:167
CurveType
Definition: hsfc.hh:175
Definition: transformation.hh:12
FieldMatrix< ctype, dimension, dimension > WorldMatrix
Definition: transformation.hh:16
FieldVector< ctype, dimension > WorldVector
Definition: transformation.hh:15