3 #ifndef DUNE_FOAMGRID_FACTORY_HH
4 #define DUNE_FOAMGRID_FACTORY_HH
13 #include <unordered_map>
16 #include <dune/common/deprecated.hh>
17 #include <dune/common/version.hh>
19 #define DUNE_FUNCTION_HH_SILENCE_DEPRECATION
20 #include <dune/common/function.hh>
21 #include <dune/common/fvector.hh>
22 #include <dune/common/hash.hh>
24 #if DUNE_VERSION_EQUAL(DUNE_COMMON, 2, 7)
25 #include <dune/common/to_unique_ptr.hh>
28 #include <dune/grid/common/gridfactory.hh>
35 template <
int dimgr
id,
int dimworld,
class ct>
37 :
public GridFactoryInterface<FoamGrid<dimgrid, dimworld, ct> >
44 #if DUNE_VERSION_LT(DUNE_COMMON, 2, 7)
46 #elif DUNE_VERSION_LT(DUNE_COMMON, 2, 8)
47 using GridPtrType = ToUniquePtr<FoamGrid<dimgrid, dimworld, ctype>>;
49 using GridPtrType = std::unique_ptr<FoamGrid<dimgrid, dimworld,ctype>>;
59 grid_->entityImps_.resize(1);
75 grid_->entityImps_.resize(1);
85 void insertVertex(
const FieldVector<ctype,dimworld>& pos)
override {
86 std::get<0>(
grid_->entityImps_[0]).emplace_back(
89 grid_->getNextFreeId()
99 return entity.impl().target_->leafIndex_;
107 return vertex.impl().target_->leafIndex_;
115 return intersection.boundarySegmentIndex();
128 std::vector<FoamGridEntityImp<0, dimgrid, dimworld, ctype>*>
vertexArray_;
134 template <
int dimgr
id,
int dimworld,
class ct>
135 class GridFactory<
FoamGrid<dimgrid, dimworld, ct> >
141 template <
int dimworld,
class ct>
170 boundarySegmentIndices_[ vertices[0] ] = this->boundarySegmentCounter_++;
177 const std::shared_ptr<BoundarySegment<dimgrid, dimworld> >& boundarySegment)
override
179 DUNE_THROW(Dune::NotImplemented,
"Parameterized boundary segments are not implemented");
186 if ( !intersection.boundary() || intersection.inside().level() != 0 )
189 const auto& vertex = intersection.inside().template subEntity<1>(intersection.indexInInside());
190 const auto& it = boundarySegmentIndices_.find( this->insertionIndex(vertex) );
191 return (it != boundarySegmentIndices_.end());
199 const std::vector<unsigned int>& vertices)
override
201 assert(type.isLine());
203 std::get<1>(this->grid_->entityImps_[0]).emplace_back(
204 this->vertexArray_[vertices[0]],
205 this->vertexArray_[vertices[1]],
207 this->grid_->getNextFreeId()
212 #if DUNE_VERSION_LTE(DUNE_COMMON, 2, 8)
213 DUNE_NO_DEPRECATED_BEGIN
219 [[deprecated(
"Signature with VirtualFunction is deprecated and will be removed after Dune 2.8. Use signature with std::function.")]]
220 void insertElement(
const GeometryType& type,
221 const std::vector<unsigned int>& vertices,
222 const std::shared_ptr<VirtualFunction<FieldVector<ctype,dimgrid>,FieldVector<ctype,dimworld> > >& elementParametrization)
override
224 assert(type.isLine());
225 EntityImp<1> newElement(this->vertexArray_[vertices[0]],
226 this->vertexArray_[vertices[1]],
228 this->grid_->getNextFreeId());
230 newElement.elementParametrization_ =
231 [elementParametrization](
const FieldVector<ctype,dimgrid>& x){
232 FieldVector<ctype,dimworld> y;
233 elementParametrization->evaluate(x, y);
237 std::get<dimgrid>(this->grid_->entityImps_[0]).push_back(newElement);
239 DUNE_NO_DEPRECATED_END
248 const std::vector<unsigned int>& vertices,
249 std::function<FieldVector<ctype,dimworld>(FieldVector<ctype,dimgrid>)> elementParametrization)
250 #if DUNE_VERSION_GT(DUNE_COMMON, 2, 7)
254 assert(type.isLine());
255 EntityImp<1> newElement(this->vertexArray_[vertices[0]],
256 this->vertexArray_[vertices[1]],
258 this->grid_->getNextFreeId());
260 newElement.elementParametrization_ = elementParametrization;
261 std::get<dimgrid>(this->grid_->entityImps_[0]).push_back(newElement);
272 if (this->grid_==
nullptr)
276 for (
auto& element : std::get<dimgrid>(this->grid_->entityImps_[0]))
278 element.vertex_[0]->elements_.push_back(&element);
279 element.vertex_[1]->elements_.push_back(&element);
283 this->grid_->setIndices();
289 for (
auto& facet : std::get<0>(this->grid_->entityImps_[0]))
291 if (facet.elements_.size() == 1)
293 const auto it = boundarySegmentIndices_.find( facet.vertex_[0]->leafIndex_ );
294 if (it != boundarySegmentIndices_.end())
295 facet.boundarySegmentIndex_ = it->second;
297 facet.boundarySegmentIndex_ = this->boundarySegmentCounter_++;
306 tmp->numBoundarySegments_ = this->boundarySegmentCounter_;
307 this->grid_ =
nullptr;
308 return GridPtrType(tmp);
312 std::unordered_map<unsigned int, unsigned int> boundarySegmentIndices_;
317 template <
int dimworld,
class ct>
346 std::array<unsigned int, 2> vertexIndices {{ vertices[0], vertices[1] }};
349 if ( vertexIndices[0] > vertexIndices[1] )
350 std::swap( vertexIndices[0], vertexIndices[1] );
352 boundarySegmentIndices_[ vertexIndices ] = this->boundarySegmentCounter_++;
359 const std::shared_ptr<BoundarySegment<dimgrid, dimworld> >& boundarySegment)
override
361 DUNE_THROW(Dune::NotImplemented,
"Parameterized boundary segments are not implemented");
368 if ( !intersection.boundary() || intersection.inside().level() != 0 )
372 const auto refElement = ReferenceElements<ctype, dimgrid>::general(intersection.inside().type());
374 const int subIdx0 = refElement.subEntity(intersection.indexInInside(), 1, 0, dimgrid);
375 const auto vertex0 = intersection.inside().template subEntity<2>( subIdx0 );
376 const int subIdx1 = refElement.subEntity(intersection.indexInInside(), 1, 1, dimgrid);
377 const auto vertex1 = intersection.inside().template subEntity<2>( subIdx1 );
379 std::array<unsigned int, 2> vertexIndices {{
380 this->insertionIndex( vertex0 ),
381 this->insertionIndex( vertex1 )
385 if ( vertexIndices[0] > vertexIndices[1] )
386 std::swap( vertexIndices[0], vertexIndices[1] );
388 const auto& it = boundarySegmentIndices_.find( vertexIndices );
389 return (it != boundarySegmentIndices_.end());
397 const std::vector<unsigned int>& vertices)
override {
399 assert(type.isTriangle());
402 newElement.vertex_[0] = this->vertexArray_[vertices[0]];
403 newElement.vertex_[1] = this->vertexArray_[vertices[1]];
404 newElement.vertex_[2] = this->vertexArray_[vertices[2]];
406 std::get<dimgrid>(this->grid_->entityImps_[0]).push_back(newElement);
410 #if DUNE_VERSION_LTE(DUNE_COMMON, 2, 8)
411 DUNE_NO_DEPRECATED_BEGIN
417 [[deprecated(
"Signature with VirtualFunction is deprecated and will be removed after Dune 2.8. Use signature with std::function.")]]
418 void insertElement(
const GeometryType& type,
419 const std::vector<unsigned int>& vertices,
420 const std::shared_ptr<VirtualFunction<FieldVector<ctype,dimgrid>,FieldVector<ctype,dimworld> > >& elementParametrization)
override
422 assert(type.isTriangle());
423 EntityImp<dimgrid> newElement(0, this->grid_->getNextFreeId());
424 newElement.vertex_[0] = this->vertexArray_[vertices[0]];
425 newElement.vertex_[1] = this->vertexArray_[vertices[1]];
426 newElement.vertex_[2] = this->vertexArray_[vertices[2]];
428 newElement.elementParametrization_ =
429 [elementParametrization](
const FieldVector<ctype,dimgrid>& x){
430 FieldVector<ctype,dimworld> y;
431 elementParametrization->evaluate(x, y);
435 std::get<dimgrid>(this->grid_->entityImps_[0]).push_back(newElement);
437 DUNE_NO_DEPRECATED_END
446 const std::vector<unsigned int>& vertices,
447 std::function<FieldVector<ctype,dimworld>(FieldVector<ctype,dimgrid>)> elementParametrization)
448 #if DUNE_VERSION_GT(DUNE_COMMON, 2, 7)
452 assert(type.isTriangle());
454 newElement.vertex_[0] = this->vertexArray_[vertices[0]];
455 newElement.vertex_[1] = this->vertexArray_[vertices[1]];
456 newElement.vertex_[2] = this->vertexArray_[vertices[2]];
457 newElement.elementParametrization_ = elementParametrization;
458 std::get<dimgrid>(this->grid_->entityImps_[0]).push_back(newElement);
468 if (this->grid_==
nullptr)
477 HashPair<const EntityImp<0>*>> edgeMap;
478 for (
auto& element : std::get<dimgrid>(this->grid_->entityImps_[0]))
480 const auto refElement = ReferenceElements<ctype, dimgrid>::general(element.type());
483 for (std::size_t i=0; i<element.facet_.size(); ++i)
488 auto v0 = element.vertex_[refElement.subEntity(i, 1, 0, 2)];
489 auto v1 = element.vertex_[refElement.subEntity(i, 1, 1, 2)];
496 auto e = edgeMap.find(std::make_pair(v0, v1));
497 if (e != edgeMap.end())
501 std::get<1>(this->grid_->entityImps_[0]).emplace_back(v0, v1, 0, this->grid_->getNextFreeId());
503 auto newEdge = &*std::get<1>(this->grid_->entityImps_[0]).rbegin();
504 edgeMap.insert(std::make_pair(std::make_pair(v0, v1), newEdge));
509 element.facet_[i] = edge;
512 edge->elements_.push_back(&element);
517 this->grid_->setIndices();
524 for (
auto& facet : std::get<1>(this->grid_->entityImps_[0]))
526 if (facet.elements_.size() == 1)
528 std::array<unsigned int, 2> vertexIndices {{ facet.vertex_[0]->leafIndex_, facet.vertex_[1]->leafIndex_ }};
530 if ( vertexIndices[0] > vertexIndices[1] )
531 std::swap( vertexIndices[0], vertexIndices[1] );
533 const auto it = boundarySegmentIndices_.find( vertexIndices );
534 if (it != boundarySegmentIndices_.end())
535 facet.boundarySegmentIndex_ = it->second;
537 facet.boundarySegmentIndex_ = this->boundarySegmentCounter_++;
546 tmp->numBoundarySegments_ = this->boundarySegmentCounter_;
547 this->grid_ =
nullptr;
548 return GridPtrType(tmp);
552 template<
class T,
class U = T>
554 std::size_t operator() (
const std::pair<T, U>& a)
const {
555 std::size_t seed = 0;
556 hash_combine(seed, a.first);
557 hash_combine(seed, a.second);
562 struct HashUIntArray {
563 std::size_t operator() (
const std::array<unsigned int, 2>& a)
const {
564 return hash_range(a.begin(), a.end());
568 std::unordered_map<std::array<unsigned int, 2>,
unsigned int, HashUIntArray> boundarySegmentIndices_;
Specialization of the generic GridFactory for FoamGrid<dimgrid, dimworld>
Definition: foamgridfactory.hh:38
~GridFactoryBase() override
Destructor.
Definition: foamgridfactory.hh:79
void insertVertex(const FieldVector< ctype, dimworld > &pos) override
Insert a vertex into the coarse grid.
Definition: foamgridfactory.hh:85
unsigned int insertionIndex(const typename FoamGrid< dimgrid, dimworld, ctype >::LeafIntersection &intersection) const override
Obtain a boundary's insertion index.
Definition: foamgridfactory.hh:113
GridFactoryBase(FoamGrid< dimgrid, dimworld, ctype > *grid)
Constructor for a given grid object.
Definition: foamgridfactory.hh:71
std::vector< FoamGridEntityImp< 0, dimgrid, dimworld, ctype > * > vertexArray_
Array containing all vertices.
Definition: foamgridfactory.hh:128
unsigned int boundarySegmentCounter_
Counter that creates the boundary segment indices.
Definition: foamgridfactory.hh:131
bool factoryOwnsGrid_
Definition: foamgridfactory.hh:125
FoamGrid< dimgrid, dimworld, ctype > * grid_
Definition: foamgridfactory.hh:121
unsigned int insertionIndex(const typename FoamGrid< dimgrid, dimworld, ctype >::Traits::template Codim< dimgrid >::Entity &vertex) const override
Obtain a vertex' insertion index.
Definition: foamgridfactory.hh:105
GridFactoryBase()
Default constructor.
Definition: foamgridfactory.hh:55
std::unique_ptr< FoamGrid< dimgrid, dimworld, ctype > > GridPtrType
Definition: foamgridfactory.hh:49
unsigned int insertionIndex(const typename FoamGrid< dimgrid, dimworld, ctype >::Traits::template Codim< 0 >::Entity &entity) const override
Obtain an element's insertion index.
Definition: foamgridfactory.hh:97
void insertElement(const GeometryType &type, const std::vector< unsigned int > &vertices, std::function< FieldVector< ctype, dimworld >(FieldVector< ctype, dimgrid >)> elementParametrization)
Insert a parametrized element into the coarse grid.
Definition: foamgridfactory.hh:247
GridPtrType createGrid() override
Finalize grid creation and hand over the grid.
Definition: foamgridfactory.hh:268
void insertBoundarySegment(const std::vector< unsigned int > &vertices, const std::shared_ptr< BoundarySegment< dimgrid, dimworld > > &boundarySegment) override
Insert a boundary segment and the boundary segment geometry This influences the ordering of the bound...
Definition: foamgridfactory.hh:176
bool wasInserted(const typename FoamGrid< dimgrid, dimworld, ctype >::LeafIntersection &intersection) const override
Return true if leaf intersection was inserted as boundary segment.
Definition: foamgridfactory.hh:184
void insertBoundarySegment(const std::vector< unsigned int > &vertices) override
Insert a boundary segment. This is only needed if you want to control the numbering of the boundary s...
Definition: foamgridfactory.hh:168
void insertElement(const GeometryType &type, const std::vector< unsigned int > &vertices) override
Insert an element into the coarse grid.
Definition: foamgridfactory.hh:198
GridFactory(FoamGrid< 1, dimworld, ctype > *grid)
Definition: foamgridfactory.hh:161
GridFactory()
Definition: foamgridfactory.hh:159
GridFactory()
Definition: foamgridfactory.hh:335
void insertBoundarySegment(const std::vector< unsigned int > &vertices) override
Insert a boundary segment. This is only needed if you want to control the numbering of the boundary s...
Definition: foamgridfactory.hh:344
void insertBoundarySegment(const std::vector< unsigned int > &vertices, const std::shared_ptr< BoundarySegment< dimgrid, dimworld > > &boundarySegment) override
Insert a boundary segment (== a line) and the boundary segment geometry This influences the ordering ...
Definition: foamgridfactory.hh:358
bool wasInserted(const typename FoamGrid< dimgrid, dimworld, ctype >::LeafIntersection &intersection) const override
Return true if leaf intersection was inserted as boundary segment.
Definition: foamgridfactory.hh:366
GridPtrType createGrid() override
Finalize grid creation and hand over the grid The receiver takes responsibility of the memory allocat...
Definition: foamgridfactory.hh:464
void insertElement(const GeometryType &type, const std::vector< unsigned int > &vertices, std::function< FieldVector< ctype, dimworld >(FieldVector< ctype, dimgrid >)> elementParametrization)
Insert a parametrized element into the coarse grid.
Definition: foamgridfactory.hh:445
void insertElement(const GeometryType &type, const std::vector< unsigned int > &vertices) override
Insert an element into the coarse grid.
Definition: foamgridfactory.hh:396
GridFactory(FoamGrid< 2, dimworld, ctype > *grid)
Definition: foamgridfactory.hh:337
The actual entity implementation.
Definition: foamgridvertex.hh:47