1 #ifndef DUNE_ALU3DGRID_FACTORY_HH
2 #define DUNE_ALU3DGRID_FACTORY_HH
10 #include <dune/common/shared_ptr.hh>
11 #include <dune/common/parallel/mpihelper.hh>
12 #include <dune/common/version.hh>
14 #include <dune/geometry/referenceelements.hh>
16 #include <dune/grid/common/gridfactory.hh>
17 #include <dune/grid/common/boundaryprojection.hh>
27 template<
class ALUGr
id >
29 :
public GridFactoryInterface< ALUGrid >
32 typedef GridFactoryInterface< ALUGrid > BaseType;
53 #if !DUNE_VERSION_NEWER(DUNE_GRID,2,5)
70 "ALU3dGridFactory supports only grids containing "
71 "tetrahedrons or hexahedrons exclusively." );
73 typedef Dune::BoundarySegmentWrapper< dimension, dimensionworld > BoundarySegmentWrapperType;
83 typedef FieldVector< ctype, dimensionworld > VertexInputType;
87 typedef FieldVector< ctype, 3 > VertexType;
89 typedef std::vector< unsigned int > ElementType;
90 typedef std::array< unsigned int, numFaceCorners > FaceType;
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;
102 typedef std::map< FaceType, const DuneBoundaryProjectionType* > BoundaryProjectionMap;
103 typedef std::vector< const DuneBoundaryProjectionType* > BoundaryProjectionVector;
105 typedef std::vector< Transformation > FaceTransformationVector;
107 static void copy (
const std::initializer_list< unsigned int > &vertices, FaceType &faceId )
109 std::copy_n( vertices.begin(), faceId.size(), faceId.begin() );
112 static FaceType makeFace (
const std::vector< unsigned int > &vertices )
114 if( vertices.size() != (
dimension == 2 ? 2 : numFaceCorners) )
115 DUNE_THROW( GridError,
"Wrong number of face vertices passed: " << vertices.size() <<
"." );
121 copy( { 0, vertices[ 1 ]+1, vertices[ 0 ]+1 }, faceId );
123 copy( { 2*vertices[ 0 ], 2*vertices[ 1 ], 2*vertices[ 0 ]+1, 2*vertices[ 1 ]+1 }, faceId );
126 std::copy_n( vertices.begin(), numFaceCorners, faceId.begin() );
130 static BndPair makeBndPair (
const FaceType &face,
int id )
133 for(
unsigned int i = 0; i < numFaceCorners; ++i )
136 bndPair.first[ j ] = face[ i ];
144 virtual Grid* createGridObj( BoundaryProjectionVector* bndProjections,
const std::string& name )
const
146 return new Grid( communicator_, globalProjection_, bndProjections , name, realGrid_ );
156 bool removeGeneratedFile =
true );
169 virtual void insertVertex (
const VertexInputType &pos );
188 const std::vector< VertexId > &vertices );
202 insertBoundary (
const GeometryType &geometry,
const std::vector< VertexId > &faceVertices,
int boundaryId = 1 );
212 virtual void insertBoundary (
int element,
int face,
int boundaryId = 1 );
230 const std::vector< VertexId > &vertices,
250 const shared_ptr<BoundarySegment<dimension,dimensionworld> >& boundarySegment ) ;
275 Grid *
createGrid (
const bool addMissingBoundaries,
const std::string dgfName =
"" );
277 Grid *
createGrid (
const bool addMissingBoundaries,
bool temporary,
const std::string dgfName =
"" );
282 alugrid_assert( Grid::getRealImplementation( entity ).getIndex() <
int(ordering_.size()) );
283 return ordering_[ Grid::getRealImplementation( entity ).getIndex() ];
291 return Grid::getRealImplementation( entity ).getIndex()/2;
294 return Grid::getRealImplementation( entity ).getIndex() - 1;
296 return Grid::getRealImplementation( entity ).getIndex();
299 virtual unsigned int insertionIndex (
const typename Grid::LevelIntersection &intersection )
const
301 return boundaryInsertionIndex( intersection.inside(), intersection.indexInInside() );
304 virtual unsigned int insertionIndex (
const typename Grid::LeafIntersection &intersection )
const
306 return boundaryInsertionIndex( intersection.inside(), intersection.indexInInside() );
309 virtual bool wasInserted (
const typename Grid::LevelIntersection &intersection )
const
311 return intersection.boundary() && (
insertionIndex(intersection) < std::numeric_limits<unsigned int>::max());
314 virtual bool wasInserted (
const typename Grid::LeafIntersection &intersection )
const
316 return intersection.boundary() && (
insertionIndex(intersection) < std::numeric_limits<unsigned int>::max());
319 const std::vector<unsigned int>&
ordering ()
const {
return ordering_; }
322 unsigned int boundaryInsertionIndex (
const typename Codim< 0 >::Entity &entity,
int face )
const
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 ) ) );
330 FaceType faceId = makeFace( vertices );
331 std::sort( faceId.begin(), faceId.end() );
333 const auto pos = insertionOrder_.find( faceId );
334 return (pos != insertionOrder_.end() ? pos->second : std::numeric_limits< unsigned int >::max());
337 void doInsertVertex (
const VertexInputType &pos,
const GlobalIdType globalId );
338 void doInsertBoundary (
int element,
int face,
int boundaryId );
343 return vertices_[ id ].second;
346 const VertexType &position (
const VertexId &
id )
const
349 return vertices_[ id ].first;
352 const VertexInputType inputPosition (
const VertexId &
id )
const
355 VertexType vertex = vertices_[ id ].first;
356 VertexInputType iVtx(0.);
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 );
372 void sortElements(
const VertexVector& vertices,
const ElementVector& elements, std::vector< unsigned int >&
ordering );
376 VertexVector vertices_;
377 ElementVector elements_;
378 BoundaryIdMap boundaryIds_,insertionOrder_;
379 PeriodicBoundaryVector periodicBoundaries_;
381 BoundaryProjectionMap boundaryProjections_;
382 FaceTransformationVector faceTransformations_;
383 unsigned int numFacesInserted_;
385 const bool allowGridGeneration_;
386 bool foundGlobalIndex_ ;
391 std::vector< unsigned int > ordering_;
396 template<
class ALUGr
id >
398 :
public std::binary_function< FaceType, FaceType, bool >
400 bool operator() (
const FaceType &a,
const FaceType &b )
const
402 for(
unsigned int i = 0; i < numFaceCorners; ++i )
404 if( a[ i ] != b[ i ] )
405 return (a[ i ] < b[ i ]);
412 template<
class ALUGr
id >
416 if( elementType ==
tetra )
418 if( !geometry.isSimplex() )
419 DUNE_THROW( GridError,
"Only simplex geometries can be inserted into "
420 "ALUGrid< 3, 3, simplex, refrule >." );
424 if( !geometry.isCube() )
425 DUNE_THROW( GridError,
"Only cube geometries can be inserted into "
426 "ALUGrid< 3, 3, cube, refrule >." );
433 template<
int dim,
int dimw, ALUGr
idElementType eltype, ALUGr
idRefinementType refinementtype ,
class Comm >
434 class GridFactory<
ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
435 :
public ALU3dGridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
437 typedef GridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
ThisType;
453 :
BaseType( filename, communicator )
457 template<
class Gr
id >
461 template<
int dim,
int dimw, ALUGr
idElementType eltype, ALUGr
idRefinementType refinementtype ,
class Comm >
463 :
public ALU3dGridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
484 template<
class ALUGr
id >
488 bool removeGeneratedFile )
490 globalProjection_ ( 0 ),
491 numFacesInserted_ ( 0 ),
493 allowGridGeneration_( rank_ == 0 ),
494 foundGlobalIndex_( false ),
495 communicator_( communicator ),
499 template<
class ALUGr
id >
505 globalProjection_ ( 0 ),
506 numFacesInserted_ ( 0 ),
508 allowGridGeneration_( rank_ == 0 ),
509 foundGlobalIndex_( false ),
510 communicator_( communicator ),
514 template<
class ALUGr
id >
520 globalProjection_ ( 0 ),
521 numFacesInserted_ ( 0 ),
522 realGrid_( realGrid ),
523 allowGridGeneration_( true ),
524 foundGlobalIndex_( false ),
525 communicator_( communicator ),
529 template<
class ALUGr
id >
533 FaceType faceId = makeFace( vertices );
535 boundaryIds_.insert( makeBndPair( faceId, 1 ) );
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." );
541 boundaryProjections_[ faceId ] =
nullptr;
543 insertionOrder_.insert( std::make_pair( faceId, insertionOrder_.size() ) );
546 template<
class ALUGr
id >
550 FaceType faceId = makeFace( vertices );
554 std::sort( faceId.begin(), faceId.end() );
555 boundaryProjections_[ faceId ] =
nullptr;
558 template<
class ALUGr
id >
561 const shared_ptr<BoundarySegment<dimension,dimensionworld> >& boundarySegment )
563 const std::size_t numVx = vertices.size();
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 );
572 typedef FieldVector< double, dimensionworld > CoordType;
573 std::vector< CoordType > coords( numVx );
574 for( std::size_t i = 0; i < numVx; ++i )
577 const std::size_t vtx = (dimension == 2 ? (elementType ==
tetra ? vertices[ i ] + 1 : 2 * vertices[ i ]) : vertices[ i ]);
583 std::copy_n( position( vtx ).begin(), dimensionworld, coords[ i ].begin() );
586 std::unique_ptr< BoundarySegmentWrapperType > prj(
new BoundarySegmentWrapperType( type, coords, boundarySegment ) );
589 for( std::size_t i = 0; i < numVx; ++i )
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." );
596 FaceType faceId = makeFace( vertices );
598 boundaryIds_.insert( makeBndPair( faceId, 1 ) );
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." );
604 boundaryProjections_[ faceId ] = prj.release();
606 insertionOrder_.insert( std::make_pair( faceId, insertionOrder_.size() ) );
610 template<
class ALUGr
id >
612 ::generateFace (
const SubEntity &subEntity, FaceType &face )
const
614 generateFace( elements_[ subEntity.first ], subEntity.second, face );
619 #if COMPILE_ALUGRID_INLINE
#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::Grid Grid
Definition: gridfactory.hh:441
BaseType::MPICommunicatorType MPICommunicatorType
Definition: gridfactory.hh:443
Definition: gridfactory.hh:458
Definition: gridfactory.hh:464
BaseType::MPICommunicatorType MPICommunicatorType
Definition: gridfactory.hh:471
ReferenceGridFactory()
Default constructor.
Definition: gridfactory.hh:474
BaseType::Grid Grid
Definition: gridfactory.hh:469
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
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