28 #ifndef EWOMS_ECFV_STENCIL_HH
29 #define EWOMS_ECFV_STENCIL_HH
34 #include <dune/grid/common/mcmgmapper.hh>
35 #include <dune/grid/common/intersectioniterator.hh>
36 #include <dune/geometry/type.hh>
37 #include <dune/common/fvector.hh>
38 #include <dune/common/version.hh>
53 template <
class Scalar,
55 bool needFaceIntegrationPos =
true,
56 bool needFaceNormal =
true>
59 enum { dimWorld = GridView::dimensionworld };
61 typedef typename GridView::ctype CoordScalar;
62 typedef typename GridView::Intersection Intersection;
63 typedef typename GridView::template Codim<0>::Entity Element;
65 typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView,
66 Dune::MCMGElementLayout> ElementMapper;
68 typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
70 typedef Dune::FieldVector<Scalar, dimWorld> WorldVector;
73 typedef Element Entity;
74 typedef ElementMapper Mapper;
76 typedef typename Element::Geometry LocalGeometry;
98 void update(
const Element& element)
103 const auto&
geometry = element_.geometry();
112 {
return centerPos_; }
118 {
return centerPos_; }
130 {
return element_.geometry(); }
133 GlobalPosition centerPos_;
141 template <
bool needNormal,
bool needIntegrationPos>
150 exteriorIdx_ =
static_cast<unsigned short>(localNeighborIdx);
153 (*normal_) = intersection.centerUnitOuterNormal();
155 const auto& geometry = intersection.geometry();
156 if (needIntegrationPos)
157 (*integrationPos_) = geometry.center();
158 area_ = geometry.volume();
180 {
return exteriorIdx_; }
187 {
return *integrationPos_; }
207 unsigned short exteriorIdx_;
210 typedef EcfvSubControlVolumeFace<needFaceIntegrationPos, needFaceNormal> SubControlVolumeFace;
213 : gridView_(gridView)
214 , elementMapper_(mapper)
217 void updateTopology(
const Element&
element)
219 auto isIt = gridView_.ibegin(element);
220 const auto& endIsIt = gridView_.iend(element);
223 subControlVolumes_.clear();
224 subControlVolumes_.emplace_back(element);
226 elements_.emplace_back(element);
228 interiorFaces_.clear();
229 boundaryFaces_.clear();
231 for (; isIt != endIsIt; ++isIt) {
232 const auto& intersection = *isIt;
236 if (intersection.neighbor()) {
237 elements_.emplace_back( intersection.outside() );
238 subControlVolumes_.emplace_back(elements_.back());
239 interiorFaces_.emplace_back(intersection, subControlVolumes_.size() - 1);
242 boundaryFaces_.emplace_back(intersection, - 10000);
247 void updatePrimaryTopology(
const Element& element)
250 subControlVolumes_.clear();
251 subControlVolumes_.emplace_back(element);
253 elements_.emplace_back(element);
256 void update(
const Element& element)
258 updateTopology(element);
261 void updateCenterGradients()
277 {
return *gridView_; }
284 {
return subControlVolumes_.size(); }
304 assert(0 <= dofIdx && dofIdx <
numDof());
306 return static_cast<unsigned>(elementMapper_.index(
element(dofIdx)));
313 {
return elements_[dofIdx]->partitionType(); }
324 assert(0 <= dofIdx && dofIdx <
numDof());
326 return elements_[dofIdx];
333 const Entity&
entity(
unsigned dofIdx)
const
343 {
return subControlVolumes_[dofIdx]; }
349 {
return interiorFaces_.size(); }
356 {
return interiorFaces_[bfIdx]; }
362 {
return boundaryFaces_.size(); }
369 {
return boundaryFaces_[bfIdx]; }
372 const GridView& gridView_;
373 const ElementMapper& elementMapper_;
375 std::vector<Element> elements_;
376 std::vector<SubControlVolume> subControlVolumes_;
377 std::vector<SubControlVolumeFace> interiorFaces_;
378 std::vector<SubControlVolumeFace> boundaryFaces_;
const GlobalPosition & integrationPos() const
Returns the global position of the face's integration point.
Definition: ecfvstencil.hh:186
size_t numDof() const
Returns the number of degrees of freedom which the current element interacts with.
Definition: ecfvstencil.hh:283
const GlobalPosition & globalPos() const
The global position associated with the sub-control volume.
Definition: ecfvstencil.hh:111
Represents the stencil (finite volume geometry) of a single element in the ECFV discretization.
Definition: ecfvstencil.hh:57
unsigned short interiorIndex() const
Returns the local index of the degree of freedom to the face's interior.
Definition: ecfvstencil.hh:165
Represents a face of a sub-control volume.
Definition: ecfvstencil.hh:142
Scalar area() const
Returns the area [m^2] of the face.
Definition: ecfvstencil.hh:199
Represents a sub-control volume.
Definition: ecfvstencil.hh:86
size_t numInteriorFaces() const
Returns the number of interior faces of the stencil.
Definition: ecfvstencil.hh:348
Dune::PartitionType partitionType(unsigned dofIdx) const
Return partition type of a given degree of freedom.
Definition: ecfvstencil.hh:312
const Entity & entity(unsigned dofIdx) const
Return the entity given the index of a degree of freedom.
Definition: ecfvstencil.hh:333
const Element & element(unsigned dofIdx) const
Return the element given the index of a degree of freedom.
Definition: ecfvstencil.hh:322
const SubControlVolumeFace & interiorFace(unsigned bfIdx) const
Returns the face object belonging to a given face index in the interior of the domain.
Definition: ecfvstencil.hh:355
const GlobalPosition & center() const
The center of the sub-control volume.
Definition: ecfvstencil.hh:117
A simple class which only stores a given member attribute if a boolean condition is true...
size_t numBoundaryFaces() const
Returns the number of boundary faces of the stencil.
Definition: ecfvstencil.hh:361
unsigned short exteriorIndex() const
Returns the local index of the degree of freedom to the face's outside.
Definition: ecfvstencil.hh:179
size_t numPrimaryDof() const
Returns the number of degrees of freedom which are contained by within the current element...
Definition: ecfvstencil.hh:295
const SubControlVolume & subControlVolume(unsigned dofIdx) const
Returns the sub-control volume object belonging to a given degree of freedom.
Definition: ecfvstencil.hh:342
const WorldVector & normal() const
Returns the outer unit normal at the face's integration point.
Definition: ecfvstencil.hh:193
Quadrature geometry for quadrilaterals.
const GridView & gridView() const
Return the grid view of the element to which the stencil refers.
Definition: ecfvstencil.hh:276
const SubControlVolumeFace & boundaryFace(unsigned bfIdx) const
Returns the boundary face object belonging to a given boundary face index.
Definition: ecfvstencil.hh:368
const LocalGeometry geometry() const
The geometry of the sub-control volume.
Definition: ecfvstencil.hh:129
Scalar volume() const
The volume [m^3] occupied by the sub-control volume.
Definition: ecfvstencil.hh:123
unsigned globalSpaceIndex(unsigned dofIdx) const
Return the global space index given the index of a degree of freedom.
Definition: ecfvstencil.hh:302
const Element & element() const
Return the element to which the stencil refers.
Definition: ecfvstencil.hh:269