3.2-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
discretization/staggered/subcontrolvolumeface.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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program 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 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program 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 this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_DISCRETIZATION_STAGGERED_SUBCONTROLVOLUMEFACE_HH
25#define DUMUX_DISCRETIZATION_STAGGERED_SUBCONTROLVOLUMEFACE_HH
26
27#include <utility>
28
29#include <dune/common/fvector.hh>
30#include <dune/geometry/type.hh>
31
34
35#include <typeinfo>
36
37namespace Dumux {
38
43template<class GridView>
45{
46 using Element = typename GridView::template Codim<0>::Entity;
47 using Intersection = typename GridView::Intersection;
48 static constexpr int codimIntersection = 1;
49
50public:
51
52 BaseStaggeredGeometryHelper(const Element& element, const GridView& gridView)
53 : element_(element)
54 , gridView_(gridView)
55 { }
56
60 template<class IntersectionMapper>
61 void updateLocalFace(const IntersectionMapper& intersectionMapper, const Intersection& intersection)
62 {
63 intersection_ = intersection;
64 }
65
69 int dofIndex() const
70 {
71 //TODO: use proper intersection mapper!
72 const auto inIdx = intersection_.indexInInside();
73 return gridView_.indexSet().subIndex(intersection_.inside(), inIdx, codimIntersection);
74 }
75
79 int localFaceIndex() const
80 {
81 return intersection_.indexInInside();
82 }
83
84private:
85 Intersection intersection_;
86 const Element element_;
87 const GridView gridView_;
88};
89
90
97template<class GridView>
99{
100 using Geometry = typename GridView::template Codim<1>::Geometry;
103 using Scalar = typename GridView::ctype;
104 using GlobalPosition = Dune::FieldVector<Scalar, GridView::dimensionworld>;
105};
106
112template<class GV,
115: public SubControlVolumeFaceBase<StaggeredSubControlVolumeFace<GV, T>, T>
116{
119 using Geometry = typename T::Geometry;
120 using GridIndexType = typename T::GridIndexType;
121 using LocalIndexType = typename T::LocalIndexType;
122
123 using Scalar = typename T::Scalar;
124 static const int dim = Geometry::mydimension;
125 static const int dimworld = Geometry::coorddimension;
126
127public:
128 using GlobalPosition = typename T::GlobalPosition;
129
131 using Traits = T;
132
133 // the default constructor
135
137 template <class Intersection, class GeometryHelper>
138 StaggeredSubControlVolumeFace(const Intersection& is,
139 const typename Intersection::Geometry& isGeometry,
140 GridIndexType scvfIndex,
141 const std::vector<GridIndexType>& scvIndices,
142 const GeometryHelper& geometryHelper)
143 : ParentType()
144 , geomType_(isGeometry.type())
145 , area_(isGeometry.volume())
146 , center_(isGeometry.center())
147 , unitOuterNormal_(is.centerUnitOuterNormal())
148 , scvfIndex_(scvfIndex)
149 , scvIndices_(scvIndices)
150 , boundary_(is.boundary())
151 {
152 corners_.resize(isGeometry.corners());
153 for (int i = 0; i < isGeometry.corners(); ++i)
154 corners_[i] = isGeometry.corner(i);
155
156 dofIdx_ = geometryHelper.dofIndex();
157 localFaceIdx_ = geometryHelper.localFaceIndex();
158 }
159
161 const GlobalPosition& center() const
162 {
163 return center_;
164 }
165
168 {
169 return center_;
170 }
171
174 {
175 // Return center for now
176 return center_;
177 }
178
180 Scalar area() const
181 {
182 return area_;
183 }
184
186 bool boundary() const
187 {
188 return boundary_;
189 }
190
193 {
194 return unitOuterNormal_;
195 }
196
198 GridIndexType insideScvIdx() const
199 {
200 return scvIndices_[0];
201 }
202
204 GridIndexType outsideScvIdx() const
205 {
206 return scvIndices_[1];
207 }
208
210 GridIndexType index() const
211 {
212 return scvfIndex_;
213 }
214
216 const GlobalPosition& corner(unsigned int localIdx) const
217 {
218 assert(localIdx < corners_.size() && "provided index exceeds the number of corners");
219 return corners_[localIdx];
220 }
221
223 const Geometry geometry() const
224 {
225 return Geometry(geomType_, corners_);
226 }
227
229 GridIndexType dofIndex() const
230 {
231 return dofIdx_;
232 }
233
235 LocalIndexType localFaceIdx() const
236 {
237 return localFaceIdx_;
238 }
239
240private:
241 Dune::GeometryType geomType_;
242 std::vector<GlobalPosition> corners_;
243 Scalar area_;
244 GlobalPosition center_;
245 GlobalPosition unitOuterNormal_;
246 GridIndexType scvfIndex_;
247 std::vector<GridIndexType> scvIndices_;
248 bool boundary_;
249
250 GridIndexType dofIdx_;
251 LocalIndexType localFaceIdx_;
252};
253
254} // end namespace Dumux
255
256#endif
Defines the index types used for grid and local indices.
Base class for a sub control volume face.
Definition: adapt.hh:29
typename GridView::IndexSet::IndexType GridIndex
Definition: indextraits.hh:39
unsigned int LocalIndex
Definition: indextraits.hh:40
defines an intersection mapper for mapping of global DOFs assigned to faces which also works for adap...
Definition: intersectionmapper.hh:185
Base class for a staggered grid geometry helper.
Definition: discretization/staggered/subcontrolvolumeface.hh:45
void updateLocalFace(const IntersectionMapper &intersectionMapper, const Intersection &intersection)
Updates the current face, i.e. sets the correct intersection.
Definition: discretization/staggered/subcontrolvolumeface.hh:61
BaseStaggeredGeometryHelper(const Element &element, const GridView &gridView)
Definition: discretization/staggered/subcontrolvolumeface.hh:52
int dofIndex() const
Returns the global dofIdx of the intersection itself.
Definition: discretization/staggered/subcontrolvolumeface.hh:69
int localFaceIndex() const
Returns the local index of the face (i.e. the intersection)
Definition: discretization/staggered/subcontrolvolumeface.hh:79
Default traits class to be used for the sub-control volume faces for the staggered finite volume sche...
Definition: discretization/staggered/subcontrolvolumeface.hh:99
Dune::FieldVector< Scalar, GridView::dimensionworld > GlobalPosition
Definition: discretization/staggered/subcontrolvolumeface.hh:104
typename GridView::template Codim< 1 >::Geometry Geometry
Definition: discretization/staggered/subcontrolvolumeface.hh:100
typename IndexTraits< GridView >::GridIndex GridIndexType
Definition: discretization/staggered/subcontrolvolumeface.hh:101
typename GridView::ctype Scalar
Definition: discretization/staggered/subcontrolvolumeface.hh:103
typename IndexTraits< GridView >::LocalIndex LocalIndexType
Definition: discretization/staggered/subcontrolvolumeface.hh:102
Class for a sub control volume face in the staggered method, i.e a part of the boundary of a sub cont...
Definition: discretization/staggered/subcontrolvolumeface.hh:116
const GlobalPosition & ipGlobal() const
The integration point for flux evaluations in global coordinates.
Definition: discretization/staggered/subcontrolvolumeface.hh:173
GridIndexType index() const
The global index of this sub control volume face.
Definition: discretization/staggered/subcontrolvolumeface.hh:210
LocalIndexType localFaceIdx() const
The local index of this sub control volume face.
Definition: discretization/staggered/subcontrolvolumeface.hh:235
Scalar area() const
The area of the sub control volume face.
Definition: discretization/staggered/subcontrolvolumeface.hh:180
T Traits
state the traits public and thus export all types
Definition: discretization/staggered/subcontrolvolumeface.hh:131
StaggeredSubControlVolumeFace(const Intersection &is, const typename Intersection::Geometry &isGeometry, GridIndexType scvfIndex, const std::vector< GridIndexType > &scvIndices, const GeometryHelper &geometryHelper)
Constructor with intersection.
Definition: discretization/staggered/subcontrolvolumeface.hh:138
GridIndexType insideScvIdx() const
Index of the inside sub control volume.
Definition: discretization/staggered/subcontrolvolumeface.hh:198
const GlobalPosition & dofPosition() const
The position of the dof living on the face.
Definition: discretization/staggered/subcontrolvolumeface.hh:167
const GlobalPosition & unitOuterNormal() const
The unit outer normal vector.
Definition: discretization/staggered/subcontrolvolumeface.hh:192
const GlobalPosition & center() const
The center of the sub control volume face.
Definition: discretization/staggered/subcontrolvolumeface.hh:161
bool boundary() const
Returns bolean if the sub control volume face is on the boundary.
Definition: discretization/staggered/subcontrolvolumeface.hh:186
const Geometry geometry() const
The geometry of the sub control volume face.
Definition: discretization/staggered/subcontrolvolumeface.hh:223
typename T::GlobalPosition GlobalPosition
Definition: discretization/staggered/subcontrolvolumeface.hh:128
GridIndexType dofIndex() const
The global index of the dof living on this face.
Definition: discretization/staggered/subcontrolvolumeface.hh:229
GridIndexType outsideScvIdx() const
Index of the outside sub control volume.
Definition: discretization/staggered/subcontrolvolumeface.hh:204
const GlobalPosition & corner(unsigned int localIdx) const
The positions of the corners.
Definition: discretization/staggered/subcontrolvolumeface.hh:216
Base class for a sub control volume face, i.e a part of the boundary of a sub control volume we compu...
Definition: subcontrolvolumefacebase.hh:41