24 :
public GridFactoryInterface< PolyhedralGrid< dim, dimworld, coord_t > >
29 const static int dimension = Grid::dimension;
30 const static int dimensionworld = Grid::dimensionworld;
33 typedef MPIHelper::MPICommunicator MPICommunicatorType;
34 typedef typename Grid::template Codim<0>::Entity Element;
35 typedef typename Grid::template Codim<dimension>::Entity Vertex;
37 typedef Dune::FieldVector<ctype,dimensionworld> CoordinateType;
38 typedef CoordinateType Coordinate;
40 using UniquePtrType = std::unique_ptr<Grid>;
43 explicit GridFactory (
const MPICommunicatorType& = MPIHelper::getCommunicator() )
49 virtual void insertVertex(
const CoordinateType& pos)
51 nodes_.push_back( pos );
64 const std::vector<unsigned int>& items)
69 std::vector< int > numbers( items.size() );
70 std::copy( items.begin(), items.end(), numbers.begin() );
72 if( type.dim() == dimension-1 )
74 faces_.push_back( numbers );
76 else if( type.dim() == dimension )
80 cells_.push_back( numbers );
84 DUNE_THROW(Dune::NotImplemented,
"insertElement not implemented for type " << type );
94 void insertBoundarySegment(
const std::vector<unsigned int>&)
96 DUNE_THROW(NotImplemented,
"yet");
99 UniquePtrType createGrid()
101 std::vector< CoordinateType >& nodes = nodes_;
102 std::vector< std::vector< int > >& faces = faces_;
103 std::vector< std::vector< int > >& cells = cells_;
107 DUNE_THROW( GridError,
"No cells found for PolyhedralGrid" );
110 const auto sumSize = [] ( std::size_t s,
const std::vector< int > &v ) {
return s + v.size(); };
111 const std::size_t numFaceNodes = std::accumulate( faces.begin(), faces.end(), std::size_t( 0 ), sumSize );
112 const std::size_t numCellFaces = std::accumulate( cells.begin(), cells.end(), std::size_t( 0 ), sumSize );
114 typename Grid::UnstructuredGridPtr ug =
115 Grid::allocateGrid( cells.size(), faces.size(), numFaceNodes, numCellFaces, nodes.size() );
120 std::map< std::vector< int >, std::vector< int > > faceMap;
123 const int nFaces = faces.size();
125 std::fill( ug->face_cells, ug->face_cells + 2*nFaces, -1 );
128 std::vector< int > faceVertices;
129 faceVertices.reserve( 30 );
130 for(
int face = 0; face < nFaces; ++face )
133 faceVertices.clear();
134 ug->face_nodepos[ face ] = facepos;
135 const int nVertices = faces[ face ].size();
136 for(
int vx = 0; vx < nVertices; ++vx, ++facepos )
139 ug->face_nodes[ facepos ] = faces[ face ][ vx ];
140 faceVertices.push_back( faces[ face ][ vx ] );
146 std::sort( faceVertices.begin(), faceVertices.end() );
148 faceMap[ faceVertices ].push_back( face );
149 assert( faceMap[ faceVertices ].size() == 1 );
152 ug->face_nodepos[ nFaces ] = facepos ;
157 const int nCells = cells.size();
159 for(
int cell = 0; cell < nCells; ++cell )
162 ug->cell_facepos[ cell ] = cellpos;
163 const int nFaces = cells[ cell ].size();
164 for(
int f = 0; f < nFaces; ++f, ++cellpos )
166 const int face = cells[ cell ][ f ];
168 ug->cell_faces[ cellpos ] = face;
171 if( ug->face_cells[ 2*face ] == -1 )
173 ug->face_cells[ 2*face ] = cell;
178 ug->face_cells[ 2*face+1 ] = cell;
183 ug->cell_facepos[ nCells ] = cellpos ;
188 const int nNodes = nodes.size();
190 for(
int vx = 0 ; vx < nNodes; ++vx )
192 for(
int d=0; d<dim; ++d, ++nodepos )
193 ug->node_coordinates[ nodepos ] = nodes[ vx ][ d ];
206 if( ug->cell_facetag )
208 std::free( ug->cell_facetag );
209 ug->cell_facetag = nullptr ;
210 for(
int i=0; i<3; ++i ) ug->cartdims[ i ] = 0;
214 Grid::computeGeometry( ug );
218 for(
int face = 0 ; face < ug->number_of_faces; ++face )
220 const int a = ug->face_cells[ 2*face ];
221 const int b = ug->face_cells[ 2*face + 1 ];
225 Coordinate centerDiff( 0 );
226 Coordinate normal( 0 );
228 for(
int d=0; d<dim; ++d )
231 centerDiff[ d ] = ug->cell_centroids[ b*dim + d ] - ug->cell_centroids[ a*dim + d ];
232 normal[ d ] = ug->face_normals[ face*dim + d ];
236 if( centerDiff * normal > 0 )
238 ug->face_cells[ 2*face ] = b;
239 ug->face_cells[ 2*face + 1 ] = a;
244 return UniquePtrType(
new Grid( std::move( ug ) ));
248 std::vector< CoordinateType > nodes_;
249 std::vector< std::vector< int > > faces_;
250 std::vector< std::vector< int > > cells_;