1#ifndef DUNE_ALU3DGRID_FACTORY_HH
2#define DUNE_ALU3DGRID_FACTORY_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;
69 "ALU3dGridFactory supports only grids containing "
70 "tetrahedrons or hexahedrons exclusively." );
73 typedef DuneBoundaryProjection< dimensionworld > DuneBoundaryProjectionType;
75 typedef Dune::BoundarySegmentWrapper< dimension, dimensionworld > BoundarySegmentWrapperType;
86 typedef FieldVector< ctype, dimensionworld > VertexInputType;
90 typedef FieldVector< ctype, 3 > VertexType;
92 typedef std::vector< unsigned int > ElementType;
93 typedef std::array< unsigned int, numFaceCorners > FaceType;
97 typedef std::vector< std::pair< VertexType, GlobalIdType > > VertexVector;
98 typedef std::vector< ElementType > ElementVector;
99 typedef std::pair< FaceType, int > BndPair ;
100 typedef std::map < FaceType, int > BoundaryIdMap;
101 typedef std::vector< std::pair< BndPair, BndPair > > PeriodicBoundaryVector;
102 typedef std::pair< unsigned int, int > SubEntity;
103 typedef std::map< FaceType, SubEntity, FaceLess > FaceMap;
105 typedef std::map< FaceType, const DuneBoundaryProjectionType* > BoundaryProjectionMap;
106 typedef std::vector< const DuneBoundaryProjectionType* > BoundaryProjectionVector;
108 typedef std::vector< Transformation > FaceTransformationVector;
110 static void copy (
const std::initializer_list< unsigned int > &vertices, FaceType &faceId )
112 std::copy_n( vertices.begin(), faceId.size(), faceId.begin() );
115 static FaceType makeFace (
const std::vector< unsigned int > &vertices )
117 if( vertices.size() != (
dimension == 2 ? 2 : numFaceCorners) )
118 DUNE_THROW( GridError,
"Wrong number of face vertices passed: " << vertices.size() <<
"." );
124 copy( { 0, vertices[ 1 ]+1, vertices[ 0 ]+1 }, faceId );
126 copy( { 2*vertices[ 0 ], 2*vertices[ 1 ], 2*vertices[ 0 ]+1, 2*vertices[ 1 ]+1 }, faceId );
129 std::copy_n( vertices.begin(), numFaceCorners, faceId.begin() );
133 static BndPair makeBndPair (
const FaceType &face,
const int id )
136 for(
unsigned int i = 0; i < numFaceCorners; ++i )
139 bndPair.first[ j ] = face[ i ];
145 void markLongestEdge( std::vector< bool >& elementOrientation,
const bool resortElements =
true ) ;
146 void markLongestEdge();
150 virtual Grid* createGridObj(
const std::string& name )
const
152 ALU3DSPACE ProjectVertexPtrPair pv = std::make_pair( globalProjection_, surfaceProjection_ );
153 return new Grid( communicator_, pv, name, realGrid_ );
163 bool removeGeneratedFile =
true );
176 virtual void insertVertex (
const VertexInputType &pos );
195 const std::vector< VertexId > &vertices );
209 insertBoundary (
const GeometryType &geometry,
const std::vector< VertexId > &faceVertices,
int boundaryId = 1 );
219 virtual void insertBoundary (
int element,
int face,
int boundaryId = 1 );
237 const std::vector< VertexId > &vertices,
238 const DuneBoundaryProjectionType *projection );
257 const std::shared_ptr<BoundarySegment<dimension,dimensionworld> >& boundarySegment ) ;
284 GridPtrType createGrid (
const bool addMissingBoundaries,
bool temporary,
const std::string dgfName =
"" );
289 alugrid_assert( entity.impl().getIndex() <
int(ordering_.size()) );
290 return ordering_[ entity.impl().getIndex() ];
298 return entity.impl().getIndex()/2;
301 return entity.impl().getIndex() - 1;
303 return entity.impl().getIndex();
306 virtual unsigned int insertionIndex (
const typename Grid::LevelIntersection &intersection )
const
308 return boundaryInsertionIndex( intersection.inside(), intersection.indexInInside() );
311 virtual unsigned int insertionIndex (
const typename Grid::LeafIntersection &intersection )
const
313 return boundaryInsertionIndex( intersection.inside(), intersection.indexInInside() );
316 virtual bool wasInserted (
const typename Grid::LevelIntersection &intersection )
const
318 return intersection.boundary() && (
insertionIndex(intersection) < std::numeric_limits<unsigned int>::max());
321 virtual bool wasInserted (
const typename Grid::LeafIntersection &intersection )
const
323 return intersection.boundary() && (
insertionIndex(intersection) < std::numeric_limits<unsigned int>::max());
326 const std::vector<unsigned int>&
ordering ()
const {
return ordering_; }
342 unsigned int boundaryInsertionIndex (
const typename Codim< 0 >::Entity &entity,
int face )
const
344 const auto& refElem = Dune::ReferenceElements< double, dimension >::general( entity.type() );
345 const int vxSize = refElem.size( face, 1,
dimension );
346 std::vector< unsigned int > vertices( vxSize );
347 for(
int i = 0; i < vxSize; ++i )
348 vertices[ i ] =
insertionIndex( entity.template subEntity< dimension >( refElem.subEntity( face, 1, i,
dimension ) ) );
350 FaceType faceId = makeFace( vertices );
351 std::sort( faceId.begin(), faceId.end() );
353 const auto pos = insertionOrder_.find( faceId );
354 return (pos != insertionOrder_.end() ? pos->second : std::numeric_limits< unsigned int >::max());
357 void doInsertVertex (
const VertexInputType &pos,
const GlobalIdType globalId );
358 void doInsertBoundary (
int element,
int face,
int boundaryId );
363 return vertices_[ id ].second;
366 const VertexType &position (
const VertexId &
id )
const
369 return vertices_[ id ].first;
372 const VertexInputType inputPosition (
const VertexId &
id )
const
375 VertexType vertex = vertices_[ id ].first;
376 VertexInputType iVtx(0.);
382 void assertGeometryType(
const GeometryType &geometry );
383 static void generateFace (
const ElementType &element,
const int f, FaceType &face );
384 void generateFace (
const SubEntity &subEntity, FaceType &face )
const;
385 void correctElementOrientation ();
386 bool identifyFaces (
const Transformation &transformation,
const FaceType &key1,
const FaceType &key2,
const int defaultId );
387 void searchPeriodicNeighbor ( FaceMap &faceMap,
typename FaceMap::iterator &pos,
const int defaultId );
388 void reinsertBoundary (
const FaceMap &faceMap,
const typename FaceMap::const_iterator &pos,
const int id );
389 void recreateBoundaryIds (
const int defaultId = 1 );
392 void sortElements(
const VertexVector& vertices,
const ElementVector& elements, std::vector< unsigned int >&
ordering );
396 VertexVector vertices_;
397 ElementVector elements_;
398 BoundaryIdMap boundaryIds_,insertionOrder_;
399 PeriodicBoundaryVector periodicBoundaries_;
400 ALU3DSPACE ProjectVertexPtr globalProjection_ ;
401 ALU3DSPACE ProjectVertexPtr surfaceProjection_ ;
402 BoundaryProjectionMap boundaryProjections_;
403 FaceTransformationVector faceTransformations_;
404 unsigned int numFacesInserted_;
406 const bool allowGridGeneration_;
407 bool foundGlobalIndex_ ;
411 typename SpaceFillingCurveOrderingType :: CurveType curveType_;
412 std::vector< unsigned int > ordering_;
414 bool markLongestEdge_;
419 template<
class ALUGr
id >
421 :
public std::binary_function< FaceType, FaceType, bool >
423 bool operator() (
const FaceType &a,
const FaceType &b )
const
425 for(
unsigned int i = 0; i < numFaceCorners; ++i )
427 if( a[ i ] != b[ i ] )
428 return (a[ i ] < b[ i ]);
435 template<
class ALUGr
id >
439 if( elementType ==
tetra )
441 if( !geometry.isSimplex() )
442 DUNE_THROW( GridError,
"Only simplex geometries can be inserted into "
443 "ALUGrid< 3, 3, simplex, refrule >." << geometry );
447 if( !geometry.isCube() )
448 DUNE_THROW( GridError,
"Only cube geometries can be inserted into "
449 "ALUGrid< 3, 3, cube, refrule >." );
456 template<
int dim,
int dimw, ALUGr
idElementType eltype, ALUGr
idRefinementType refinementtype ,
class Comm >
457 class GridFactory<
ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
458 :
public ALU3dGridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
460 typedef GridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
ThisType;
474 template <
class MPIComm>
482 :
BaseType( filename, communicator )
486 template <
class MPIComm>
493 template<
class Gr
id >
497 template<
int dim,
int dimw, ALUGr
idElementType eltype, ALUGr
idRefinementType refinementtype ,
class Comm >
499 :
public ALU3dGridFactory< ALUGrid< dim, dimw, eltype, refinementtype, Comm > >
520 template<
class ALUGr
id >
522 ALU3dGridFactory< ALUGrid >
524 bool removeGeneratedFile )
526 globalProjection_ ( 0 ),
527 surfaceProjection_ ( 0 ),
528 numFacesInserted_ ( 0 ),
530 allowGridGeneration_( rank_ == 0 ),
531 foundGlobalIndex_( false ),
532 communicator_( communicator ),
534 markLongestEdge_(
ALUGrid::dimension == 2 )
536 BoundarySegmentWrapperType::registerFactory();
540 template<
class ALUGr
id >
546 globalProjection_ ( 0 ),
547 surfaceProjection_ ( 0 ),
548 numFacesInserted_ ( 0 ),
550 allowGridGeneration_( rank_ == 0 ),
551 foundGlobalIndex_( false ),
552 communicator_( communicator ),
554 markLongestEdge_(
ALUGrid::dimension == 2 )
556 BoundarySegmentWrapperType::registerFactory();
560 template<
class ALUGr
id >
566 globalProjection_ ( 0 ),
567 surfaceProjection_ ( 0 ),
568 numFacesInserted_ ( 0 ),
569 realGrid_( realGrid ),
570 allowGridGeneration_( true ),
571 foundGlobalIndex_( false ),
572 communicator_( communicator ),
574 markLongestEdge_(
ALUGrid::dimension == 2 )
576 BoundarySegmentWrapperType::registerFactory();
580 template<
class ALUGr
id >
584 FaceType faceId = makeFace( vertices );
586 boundaryIds_.insert( makeBndPair( faceId, 1 ) );
588 std::sort( faceId.begin(), faceId.end() );
589 if( boundaryProjections_.find( faceId ) != boundaryProjections_.end() )
590 DUNE_THROW( GridError,
"Only one boundary projection can be attached to a face." );
592 boundaryProjections_[ faceId ] =
nullptr;
594 insertionOrder_.insert( std::make_pair( faceId, insertionOrder_.size() ) );
597 template<
class ALUGr
id >
601 FaceType faceId = makeFace( vertices );
603 boundaryIds_.insert( makeBndPair( faceId,
ALU3DSPACE ProcessorBoundary_t ) );
605 std::sort( faceId.begin(), faceId.end() );
606 boundaryProjections_[ faceId ] =
nullptr;
609 template<
class ALUGr
id >
612 const std::shared_ptr<BoundarySegment<dimension,dimensionworld> >& boundarySegment )
614 const std::size_t numVx = vertices.size();
616 GeometryType type = (elementType ==
tetra) ?
617 GeometryTypes::simplex(dimension-1) :
618 GeometryTypes::cube(dimension-1);
622 typedef FieldVector< double, dimensionworld > CoordType;
623 std::vector< CoordType > coords( numVx );
624 for( std::size_t i = 0; i < numVx; ++i )
627 const std::size_t vtx = (dimension == 2 ? (elementType ==
tetra ? vertices[ i ] + 1 : 2 * vertices[ i ]) : vertices[ i ]);
633 std::copy_n( position( vtx ).begin(), dimensionworld, coords[ i ].begin() );
636 std::unique_ptr< BoundarySegmentWrapperType > prj(
new BoundarySegmentWrapperType( type, coords, boundarySegment ) );
639 for( std::size_t i = 0; i < numVx; ++i )
641 CoordType global = (*prj)( coords [ i ] );
642 if( (global - coords[ i ]).two_norm() > 1e-6 )
643 DUNE_THROW( GridError,
"BoundarySegment does not map face vertices to face vertices." );
646 FaceType faceId = makeFace( vertices );
648 boundaryIds_.insert( makeBndPair( faceId, 1 ) );
650 std::sort( faceId.begin(), faceId.end() );
651 if( boundaryProjections_.find( faceId ) != boundaryProjections_.end() )
652 DUNE_THROW( GridError,
"Only one boundary projection can be attached to a face." );
654 boundaryProjections_[ faceId ] = prj.release();
656 insertionOrder_.insert( std::make_pair( faceId, insertionOrder_.size() ) );
660 template<
class ALUGr
id >
662 ::generateFace (
const SubEntity &subEntity, FaceType &face )
const
664 generateFace( elements_[ subEntity.first ], subEntity.second, face );
669#if COMPILE_ALUGRID_INLINE
#define ALU3DSPACE
Definition: alu3dinclude.hh:7
Provides base classes for ALUGrid.
#define alugrid_assert(EX)
Definition: alugrid_assert.hh:20
Definition: alu3dinclude.hh:33
Definition: alu3dinclude.hh:63
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:46
@ dimension
Definition: alugrid.hh:44
@ dimensionworld
Definition: alugrid.hh:44
BaseType::MPICommunicatorType MPICommunicatorType
Definition: alugrid.hh:39
Factory class for ALUGrids.
Definition: gridfactory.hh:30
Communication comm() const
Return the Communication used by the grid factory.
Definition: gridfactory.hh:336
unsigned int VertexId
Definition: gridfactory.hh:52
GridPtrType createGrid()
finalize the grid creation and hand over the grid
Definition: gridfactory.cc:476
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:34
Grid::MPICommunicatorType MPICommunicatorType
Definition: gridfactory.hh:44
virtual unsigned int insertionIndex(const typename Codim< 0 >::Entity &entity) const
Definition: gridfactory.hh:287
virtual bool wasInserted(const typename Grid::LeafIntersection &intersection) const
Definition: gridfactory.hh:321
void insertFaceTransformation(const WorldMatrix &matrix, const WorldVector &shift)
add a face transformation (for periodic identification)
Definition: gridfactory.cc:242
virtual unsigned int insertionIndex(const typename Codim< dimension >::Entity &entity) const
Definition: gridfactory.hh:294
ALU3dGridFactory(const MPICommunicatorType &communicator=Grid::defaultCommunicator(), bool removeGeneratedFile=true)
default constructor
Definition: gridfactory.hh:523
ALUGrid Grid
Definition: gridfactory.hh:35
const std::vector< unsigned int > & ordering() const
Definition: gridfactory.hh:326
Grid::CollectiveCommunication Communication
Definition: gridfactory.hh:62
Transformation::WorldMatrix WorldMatrix
type of matrix from world coordinates to world coordinates
Definition: gridfactory.hh:60
virtual void insertBoundarySegment(const std::vector< VertexId > &vertices)
insert a boundary segment into the macro grid
Definition: gridfactory.hh:582
ALUGridTransformation< ctype, dimensionworld > Transformation
Definition: gridfactory.hh:55
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:148
void insertProcessBorder(int element, int face)
Definition: gridfactory.hh:222
Transformation::WorldVector WorldVector
type of vector for world coordinates
Definition: gridfactory.hh:58
virtual void insertProcessBorder(const std::vector< VertexId > &vertices)
Definition: gridfactory.hh:599
ALU3dGridFactory(const std::string &filename, const MPICommunicatorType &communicator=Grid::defaultCommunicator())
constructor taking filename for temporary outfile
Definition: gridfactory.hh:543
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:98
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:221
virtual unsigned int insertionIndex(const typename Grid::LevelIntersection &intersection) const
Definition: gridfactory.hh:306
unsigned int GlobalIdType
Definition: gridfactory.hh:53
virtual ~ALU3dGridFactory()
Destructor.
Definition: gridfactory.cc:28
virtual bool wasInserted(const typename Grid::LevelIntersection &intersection) const
Definition: gridfactory.hh:316
Grid::ctype ctype
Definition: gridfactory.hh:37
void setLongestEdgeFlag(bool flag=true)
set longest edge marking for biscetion grids (default is off)
Definition: gridfactory.hh:330
ALU3dGridFactory(const bool verbose, const MPICommunicatorType &communicator)
constructor taking verbose flag
Definition: gridfactory.hh:563
static const ALU3dGridElementType elementType
Definition: gridfactory.hh:39
virtual void insertBoundarySegment(const std::vector< VertexId > &vertices, const std::shared_ptr< BoundarySegment< dimension, dimensionworld > > &boundarySegment)
insert a shaped boundary segment into the macro grid
Definition: gridfactory.hh:611
virtual unsigned int insertionIndex(const typename Grid::LeafIntersection &intersection) const
Definition: gridfactory.hh:311
decltype(std::declval< Dune::GridFactoryInterface< Grid > * >() ->createGrid() GridPtrType)
Definition: gridfactory.hh:65
Definition: 3d/grid.hh:117
Definition: gridfactory.hh:48
Grid::template Codim< codim >::Entity Entity
Definition: gridfactory.hh:49
Definition: gridfactory.hh:422
Specialization of the generic GridFactory for ALUGrid.
Definition: gridfactory.hh:459
GridFactory(const MPICommunicatorType &communicator=Grid::defaultCommunicator())
Default constructor.
Definition: gridfactory.hh:469
GridFactory(const std::string &filename, const MPIComm &)
constructor taking filename and ignoring MPIComm
Definition: gridfactory.hh:487
GridFactory(const MPIComm &)
Default constructor ignoring MPIComm.
Definition: gridfactory.hh:475
GridFactory(const std::string &filename, const MPICommunicatorType &communicator=Grid::defaultCommunicator())
constructor taking filename
Definition: gridfactory.hh:480
BaseType::Grid Grid
Definition: gridfactory.hh:464
BaseType::MPICommunicatorType MPICommunicatorType
Definition: gridfactory.hh:466
Definition: gridfactory.hh:494
Definition: gridfactory.hh:500
BaseType::MPICommunicatorType MPICommunicatorType
Definition: gridfactory.hh:507
ReferenceGridFactory()
Default constructor.
Definition: gridfactory.hh:510
BaseType::Grid Grid
Definition: gridfactory.hh:505
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.cc:497
ALUGrid boundary projection implementation DuneBndProjection has to fulfil the DuneBoundaryProjection...
Definition: bndprojection.hh:14
static void registerFactory()
Definition: bndprojection.hh:91
Definition: transformation.hh:12
FieldMatrix< ctype, dimension, dimension > WorldMatrix
Definition: transformation.hh:16
FieldVector< ctype, dimension > WorldVector
Definition: transformation.hh:15