All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
fvbaseboundarycontext.hh
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
28 #ifndef EWOMS_FV_BASE_BOUNDARY_CONTEXT_HH
29 #define EWOMS_FV_BASE_BOUNDARY_CONTEXT_HH
30 
31 #include "fvbaseproperties.hh"
32 
33 #include <opm/common/Unused.hpp>
34 
35 #include <dune/common/fvector.hh>
36 
37 namespace Ewoms {
38 
44 template<class TypeTag>
46 {
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;
55 
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;
60 
61  enum { dim = GridView::dimension };
62  enum { dimWorld = GridView::dimensionworld };
63 
64  typedef typename GridView::ctype CoordScalar;
65  typedef Dune::FieldVector<CoordScalar, dim> GlobalPosition;
66  typedef Dune::FieldVector<Scalar, dimWorld> Vector;
67 
68 public:
72  explicit FvBaseBoundaryContext(const ElementContext& elemCtx)
73  : elemCtx_(elemCtx)
74  , intersectionIt_(gridView().ibegin(element()))
75  { }
76 
80  const Problem& problem() const
81  { return elemCtx_.problem(); }
82 
86  const Model& model() const
87  { return elemCtx_.model(); }
88 
92  const GridView& gridView() const
93  { return elemCtx_.gridView(); }
94 
98  const Element& element() const
99  { return elemCtx_.element(); }
100 
104  const ElementContext& elementContext() const
105  { return elemCtx_; }
106 
110  const GradientCalculator& gradientCalculator() const
111  { return elemCtx_.gradientCalculator(); }
112 
116  size_t numDof(unsigned timeIdx) const
117  { return elemCtx_.numDof(timeIdx); }
118 
122  size_t numPrimaryDof(unsigned timeIdx) const
123  { return elemCtx_.numPrimaryDof(timeIdx); }
124 
128  size_t numInteriorFaces(unsigned timeIdx) const
129  { return elemCtx_.numInteriorFaces(timeIdx); }
130 
134  size_t numBoundaryFaces(unsigned timeIdx) const
135  { return elemCtx_.stencil(timeIdx).numBoundaryFaces(); }
136 
140  const Stencil& stencil(unsigned timeIdx) const
141  { return elemCtx_.stencil(timeIdx); }
142 
149  Vector normal(unsigned boundaryFaceIdx, unsigned timeIdx) const
150  {
151  auto tmp = stencil(timeIdx).boundaryFace[boundaryFaceIdx].normal;
152  tmp /= tmp.two_norm();
153  return tmp;
154  }
155 
159  Scalar boundarySegmentArea(unsigned boundaryFaceIdx, unsigned timeIdx) const
160  { return elemCtx_.stencil(timeIdx).boundaryFace(boundaryFaceIdx).area(); }
161 
168  const GlobalPosition& pos(unsigned boundaryFaceIdx, unsigned timeIdx) const
169  { return stencil(timeIdx).boundaryFace(boundaryFaceIdx).integrationPos(); }
170 
177  const GlobalPosition& cvCenter(unsigned boundaryFaceIdx, unsigned timeIdx) const
178  {
179  unsigned scvIdx = stencil(timeIdx).boundaryFace(boundaryFaceIdx).interiorIndex();
180  return stencil(timeIdx).subControlVolume(scvIdx).globalPos();
181  }
182 
187  unsigned focusDofIndex() const
188  { return elemCtx_.focusDofIndex(); }
189 
197  unsigned interiorScvIndex(unsigned boundaryFaceIdx, unsigned timeIdx) const
198  { return stencil(timeIdx).boundaryFace(boundaryFaceIdx).interiorIndex(); }
199 
207  unsigned globalSpaceIndex(unsigned boundaryFaceIdx, unsigned timeIdx) const
208  { return elemCtx_.globalSpaceIndex(interiorScvIndex(boundaryFaceIdx, timeIdx), timeIdx); }
209 
217  const IntensiveQuantities& intensiveQuantities(unsigned boundaryFaceIdx, unsigned timeIdx) const
218  {
219  unsigned interiorScvIdx = this->interiorScvIndex(boundaryFaceIdx, timeIdx);
220  return elemCtx_.intensiveQuantities(interiorScvIdx, timeIdx);
221  }
222 
229  const ExtensiveQuantities& extensiveQuantities(unsigned boundaryFaceIdx, unsigned timeIdx) const
230  { return elemCtx_.boundaryExtensiveQuantities(boundaryFaceIdx, timeIdx); }
231 
241  const Intersection& intersection(unsigned boundaryFaceIdx OPM_UNUSED) const
242  { return *intersectionIt_; }
243 
252  IntersectionIterator& intersectionIt()
253  { return intersectionIt_; }
254 
255 protected:
256  const ElementContext& elemCtx_;
257  IntersectionIterator intersectionIt_;
258 };
259 
260 } // namespace Ewoms
261 
262 #endif
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&#39;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