28 #ifndef EWOMS_FV_BASE_BOUNDARY_CONTEXT_HH
29 #define EWOMS_FV_BASE_BOUNDARY_CONTEXT_HH
33 #include <opm/common/Unused.hpp>
35 #include <dune/common/fvector.hh>
44 template<
class TypeTag>
47 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
48 typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
49 typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
50 typedef typename GET_PROP_TYPE(TypeTag, Stencil) Stencil;
51 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
52 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
53 typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
54 typedef typename GET_PROP_TYPE(TypeTag, GradientCalculator) GradientCalculator;
56 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
57 typedef typename GridView::template Codim<0>::Entity Element;
58 typedef typename GridView::IntersectionIterator IntersectionIterator;
59 typedef typename GridView::Intersection Intersection;
61 enum { dim = GridView::dimension };
62 enum { dimWorld = GridView::dimensionworld };
64 typedef typename GridView::ctype CoordScalar;
65 typedef Dune::FieldVector<CoordScalar, dim> GlobalPosition;
66 typedef Dune::FieldVector<Scalar, dimWorld> Vector;
81 {
return elemCtx_.problem(); }
87 {
return elemCtx_.model(); }
93 {
return elemCtx_.gridView(); }
99 {
return elemCtx_.element(); }
111 {
return elemCtx_.gradientCalculator(); }
117 {
return elemCtx_.numDof(timeIdx); }
123 {
return elemCtx_.numPrimaryDof(timeIdx); }
129 {
return elemCtx_.numInteriorFaces(timeIdx); }
135 {
return elemCtx_.stencil(timeIdx).numBoundaryFaces(); }
140 const Stencil&
stencil(
unsigned timeIdx)
const
141 {
return elemCtx_.stencil(timeIdx); }
149 Vector
normal(
unsigned boundaryFaceIdx,
unsigned timeIdx)
const
151 auto tmp =
stencil(timeIdx).boundaryFace[boundaryFaceIdx].normal;
152 tmp /= tmp.two_norm();
160 {
return elemCtx_.stencil(timeIdx).boundaryFace(boundaryFaceIdx).area(); }
168 const GlobalPosition&
pos(
unsigned boundaryFaceIdx,
unsigned timeIdx)
const
169 {
return stencil(timeIdx).boundaryFace(boundaryFaceIdx).integrationPos(); }
177 const GlobalPosition&
cvCenter(
unsigned boundaryFaceIdx,
unsigned timeIdx)
const
179 unsigned scvIdx =
stencil(timeIdx).boundaryFace(boundaryFaceIdx).interiorIndex();
180 return stencil(timeIdx).subControlVolume(scvIdx).globalPos();
188 {
return elemCtx_.focusDofIndex(); }
198 {
return stencil(timeIdx).boundaryFace(boundaryFaceIdx).interiorIndex(); }
208 {
return elemCtx_.globalSpaceIndex(
interiorScvIndex(boundaryFaceIdx, timeIdx), timeIdx); }
220 return elemCtx_.intensiveQuantities(interiorScvIdx, timeIdx);
230 {
return elemCtx_.boundaryExtensiveQuantities(boundaryFaceIdx, timeIdx); }
241 const Intersection&
intersection(
unsigned boundaryFaceIdx OPM_UNUSED)
const
242 {
return *intersectionIt_; }
253 {
return intersectionIt_; }
256 const ElementContext& elemCtx_;
257 IntersectionIterator intersectionIt_;
size_t numPrimaryDof(unsigned timeIdx) const
Definition: fvbaseboundarycontext.hh:122
IntersectionIterator & intersectionIt()
Return the intersection for the neumann segment.
Definition: fvbaseboundarycontext.hh:252
unsigned focusDofIndex() const
Return the local sub-control volume index upon which the linearization is currently focused...
Definition: fvbaseboundarycontext.hh:187
const GridView & gridView() const
Definition: fvbaseboundarycontext.hh:92
const Problem & problem() const
Definition: fvbaseboundarycontext.hh:80
unsigned globalSpaceIndex(unsigned boundaryFaceIdx, unsigned timeIdx) const
Return the global space index of the sub-control volume at the interior of a boundary segment...
Definition: fvbaseboundarycontext.hh:207
Represents all quantities which available on boundary segments.
Definition: fvbaseboundarycontext.hh:45
const GlobalPosition & pos(unsigned boundaryFaceIdx, unsigned timeIdx) const
Return the position of a local entity in global coordinates.
Definition: fvbaseboundarycontext.hh:168
const Element & element() const
Definition: fvbaseboundarycontext.hh:98
unsigned interiorScvIndex(unsigned boundaryFaceIdx, unsigned timeIdx) const
Return the local sub-control volume index of the interior of a boundary segment.
Definition: fvbaseboundarycontext.hh:197
Scalar boundarySegmentArea(unsigned boundaryFaceIdx, unsigned timeIdx) const
Returns the area [m^2] of a given boudary segment.
Definition: fvbaseboundarycontext.hh:159
const ExtensiveQuantities & extensiveQuantities(unsigned boundaryFaceIdx, unsigned timeIdx) const
Return the extensive quantities for a given boundary face.
Definition: fvbaseboundarycontext.hh:229
Declare the properties used by the infrastructure code of the finite volume discretizations.
size_t numBoundaryFaces(unsigned timeIdx) const
Return the number of boundary segments of the current element.
Definition: fvbaseboundarycontext.hh:134
FvBaseBoundaryContext(const ElementContext &elemCtx)
The constructor.
Definition: fvbaseboundarycontext.hh:72
const Model & model() const
Definition: fvbaseboundarycontext.hh:86
size_t numInteriorFaces(unsigned timeIdx) const
Definition: fvbaseboundarycontext.hh:128
const GlobalPosition & cvCenter(unsigned boundaryFaceIdx, unsigned timeIdx) const
Return the position of a control volume's center in global coordinates.
Definition: fvbaseboundarycontext.hh:177
const Stencil & stencil(unsigned timeIdx) const
Definition: fvbaseboundarycontext.hh:140
const GradientCalculator & gradientCalculator() const
Returns a reference to the current gradient calculator.
Definition: fvbaseboundarycontext.hh:110
Vector normal(unsigned boundaryFaceIdx, unsigned timeIdx) const
Returns the outer unit normal of the boundary segment.
Definition: fvbaseboundarycontext.hh:149
size_t numDof(unsigned timeIdx) const
Definition: fvbaseboundarycontext.hh:116
const ElementContext & elementContext() const
Returns a reference to the element context object.
Definition: fvbaseboundarycontext.hh:104
const Intersection & intersection(unsigned boundaryFaceIdx OPM_UNUSED) const
Return the intersection for the neumann segment.
Definition: fvbaseboundarycontext.hh:241
const IntensiveQuantities & intensiveQuantities(unsigned boundaryFaceIdx, unsigned timeIdx) const
Return the intensive quantities for the finite volume in the interiour of a boundary segment...
Definition: fvbaseboundarycontext.hh:217