28#ifndef EWOMS_ECFV_STENCIL_HH
29#define EWOMS_ECFV_STENCIL_HH
33#include <opm/material/common/ConditionalStorage.hpp>
35#include <dune/grid/common/mcmgmapper.hh>
36#include <dune/grid/common/intersectioniterator.hh>
37#include <dune/geometry/type.hh>
38#include <dune/common/fvector.hh>
39#include <dune/common/version.hh>
40#include <opm/common/ErrorMacros.hpp>
41#include <opm/input/eclipse/EclipseState/Grid/FaceDir.hpp>
56template <
class Scalar,
62 enum { dimWorld = GridView::dimensionworld };
64 using CoordScalar =
typename GridView::ctype;
65 using Intersection =
typename GridView::Intersection;
66 using Element =
typename GridView::template Codim<0>::Entity;
68 using ElementMapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
70 using GlobalPosition = Dune::FieldVector<CoordScalar, dimWorld>;
72 using WorldVector = Dune::FieldVector<Scalar, dimWorld>;
75 using Entity = Element ;
76 using Mapper = ElementMapper ;
78 using LocalGeometry =
typename Element::Geometry;
100 void update(
const Element&
element)
110 {
return element_.geometry().center(); }
116 {
return element_.geometry().center(); }
122 {
return element_.geometry().volume(); }
128 {
return element_.geometry(); }
134 {
return element_.geometryInFather(); }
143 template <
bool needIntegrationPos,
bool needNormal>
155 (*normal_) = intersection.centerUnitOuterNormal();
156 const auto& geometry = intersection.geometry();
158 (*integrationPos_) = geometry.center();
159 area_ = geometry.volume();
164 dirId_ = intersection.indexInInside();
186 {
return exteriorIdx_; }
193 {
return *integrationPos_; }
224 using Dir = FaceDir::DirEnum;
226 OPM_THROW(std::runtime_error,
"NNC faces does not have a face id");
240 "Unexpected face id" + std::to_string(dirId_));
250 unsigned short exteriorIdx_;
254 using BoundaryFace = EcfvSubControlVolumeFace<
true,
needFaceNormal>;
264 void updateTopology(
const Element&
element)
270 subControlVolumes_.clear();
271 subControlVolumes_.emplace_back(
element);
273 elements_.emplace_back(
element);
275 interiorFaces_.clear();
276 boundaryFaces_.clear();
279 const auto& intersection = *
isIt;
283 if (intersection.neighbor()) {
284 elements_.emplace_back( intersection.outside() );
285 subControlVolumes_.emplace_back(elements_.back());
286 interiorFaces_.emplace_back(intersection, subControlVolumes_.size() - 1);
289 boundaryFaces_.emplace_back(intersection, - 10000);
294 void updatePrimaryTopology(
const Element&
element)
297 subControlVolumes_.clear();
298 subControlVolumes_.emplace_back(
element);
300 elements_.emplace_back(
element);
303 void update(
const Element&
element)
308 void updateCenterGradients()
324 {
return *gridView_; }
331 {
return subControlVolumes_.size(); }
353 return static_cast<unsigned>(elementMapper_.index(
element(
dofIdx)));
360 {
return elements_[
dofIdx]->partitionType(); }
390 {
return subControlVolumes_[
dofIdx]; }
396 {
return interiorFaces_.size(); }
403 {
return interiorFaces_[
faceIdx]; }
409 {
return boundaryFaces_.size(); }
416 {
return boundaryFaces_[
bfIdx]; }
419 const GridView& gridView_;
420 const ElementMapper& elementMapper_;
422 std::vector<Element> elements_;
423 std::vector<SubControlVolume> subControlVolumes_;
424 std::vector<SubControlVolumeFace> interiorFaces_;
425 std::vector<BoundaryFace> boundaryFaces_;
Represents a face of a sub-control volume.
Definition ecfvstencil.hh:145
const GlobalPosition & integrationPos() const
Returns the global position of the face's integration point.
Definition ecfvstencil.hh:192
unsigned short interiorIndex() const
Returns the local index of the degree of freedom to the face's interior.
Definition ecfvstencil.hh:171
unsigned short exteriorIndex() const
Returns the local index of the degree of freedom to the face's outside.
Definition ecfvstencil.hh:185
const WorldVector & normal() const
Returns the outer unit normal at the face's integration point.
Definition ecfvstencil.hh:199
FaceDir::DirEnum faceDirFromDirId() const
Returns the direction of the face.
Definition ecfvstencil.hh:222
Scalar area() const
Returns the area [m^2] of the face.
Definition ecfvstencil.hh:205
int dirId() const
Returns the direction id of the face w.r.t the cell.
Definition ecfvstencil.hh:214
Represents a sub-control volume.
Definition ecfvstencil.hh:89
const LocalGeometry geometry() const
The geometry of the sub-control volume.
Definition ecfvstencil.hh:127
Scalar volume() const
The volume [m^3] occupied by the sub-control volume.
Definition ecfvstencil.hh:121
const LocalGeometry localGeometry() const
Geometry of the sub-control volume relative to parent.
Definition ecfvstencil.hh:133
decltype(auto) center() const
The center of the sub-control volume.
Definition ecfvstencil.hh:115
decltype(auto) globalPos() const
The global position associated with the sub-control volume.
Definition ecfvstencil.hh:109
Represents the stencil (finite volume geometry) of a single element in the ECFV discretization.
Definition ecfvstencil.hh:61
const SubControlVolumeFace & interiorFace(unsigned faceIdx) const
Returns the face object belonging to a given face index in the interior of the domain.
Definition ecfvstencil.hh:402
const BoundaryFace & boundaryFace(unsigned bfIdx) const
Returns the boundary face object belonging to a given boundary face index.
Definition ecfvstencil.hh:415
const Element & element(unsigned dofIdx) const
Return the element given the index of a degree of freedom.
Definition ecfvstencil.hh:369
size_t numPrimaryDof() const
Returns the number of degrees of freedom which are contained by within the current element.
Definition ecfvstencil.hh:342
const Element & element() const
Return the element to which the stencil refers.
Definition ecfvstencil.hh:316
const Entity & entity(unsigned dofIdx) const
Return the entity given the index of a degree of freedom.
Definition ecfvstencil.hh:380
size_t numBoundaryFaces() const
Returns the number of boundary faces of the stencil.
Definition ecfvstencil.hh:408
unsigned globalSpaceIndex(unsigned dofIdx) const
Return the global space index given the index of a degree of freedom.
Definition ecfvstencil.hh:349
const GridView & gridView() const
Return the grid view of the element to which the stencil refers.
Definition ecfvstencil.hh:323
size_t numInteriorFaces() const
Returns the number of interior faces of the stencil.
Definition ecfvstencil.hh:395
Dune::PartitionType partitionType(unsigned dofIdx) const
Return partition type of a given degree of freedom.
Definition ecfvstencil.hh:359
size_t numDof() const
Returns the number of degrees of freedom which the current element interacts with.
Definition ecfvstencil.hh:330
const SubControlVolume & subControlVolume(unsigned dofIdx) const
Returns the sub-control volume object belonging to a given degree of freedom.
Definition ecfvstencil.hh:389
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Quadrature geometry for quadrilaterals.