version 3.8.0
Loading...
Searching...
No Matches
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// SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
12#ifndef DUMUX_DISCRETIZATION_STAGGERED_SUBCONTROLVOLUMEFACE_HH
13#define DUMUX_DISCRETIZATION_STAGGERED_SUBCONTROLVOLUMEFACE_HH
14
15#include <utility>
16
17#include <dune/common/fvector.hh>
18#include <dune/geometry/type.hh>
19
22
23#include <typeinfo>
24
25namespace Dumux {
26
31template<class GridView>
33{
34 using Element = typename GridView::template Codim<0>::Entity;
35 using Intersection = typename GridView::Intersection;
36 static constexpr int codimIntersection = 1;
37
38public:
39
40 BaseStaggeredGeometryHelper(const Element& element, const GridView& gridView)
41 : element_(element)
42 , gridView_(gridView)
43 { }
44
48 template<class IntersectionMapper>
49 void updateLocalFace(const IntersectionMapper& intersectionMapper, const Intersection& intersection)
50 {
51 intersection_ = intersection;
52 }
53
57 int dofIndex() const
58 {
59 //TODO: use proper intersection mapper!
60 const auto inIdx = intersection_.indexInInside();
61 return gridView_.indexSet().subIndex(intersection_.inside(), inIdx, codimIntersection);
62 }
63
67 int localFaceIndex() const
68 {
69 return intersection_.indexInInside();
70 }
71
72private:
73 Intersection intersection_;
74 const Element element_;
75 const GridView gridView_;
76};
77
78
85template<class GridView>
87{
90 using Scalar = typename GridView::ctype;
91 using Element = typename GridView::template Codim<0>::Entity;
92 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
93
94 static constexpr int dim = GridView::Grid::dimension;
95 static constexpr int dimWorld = GridView::Grid::dimensionworld;
96 using Geometry = Dune::AxisAlignedCubeGeometry<Scalar, dim-1, dimWorld>;
97};
98
104template<class GV,
107: public SubControlVolumeFaceBase<StaggeredSubControlVolumeFace<GV, T>, T>
108{
111 using Geometry = typename T::Geometry;
112 using GridIndexType = typename T::GridIndexType;
113 using LocalIndexType = typename T::LocalIndexType;
114
115 using Scalar = typename T::Scalar;
116 static const int dim = Geometry::mydimension;
117 static const int dimworld = Geometry::coorddimension;
118
119public:
120 using GlobalPosition = typename T::GlobalPosition;
121
123 using Traits = T;
124
125 // the default constructor
127
129 template <class Intersection, class GeometryHelper>
130 StaggeredSubControlVolumeFace(const Intersection& is,
131 const typename Intersection::Geometry& isGeometry,
132 GridIndexType scvfIndex,
133 const std::vector<GridIndexType>& scvIndices,
134 const GeometryHelper& geometryHelper)
135 : ParentType()
136 , area_(isGeometry.volume())
137 , center_(isGeometry.center())
138 , unitOuterNormal_(is.centerUnitOuterNormal())
139 , scvfIndex_(scvfIndex)
140 , scvIndices_(scvIndices)
141 , boundary_(is.boundary())
142 {
143 dofIdx_ = geometryHelper.dofIndex();
144 localFaceIdx_ = geometryHelper.localFaceIndex();
145 }
146
148 const GlobalPosition& center() const
149 {
150 return center_;
151 }
152
155 {
156 return center_;
157 }
158
161 {
162 // Return center for now
163 return center_;
164 }
165
167 Scalar area() const
168 {
169 return area_;
170 }
171
173 bool boundary() const
174 {
175 return boundary_;
176 }
177
180 {
181 return unitOuterNormal_;
182 }
183
185 GridIndexType insideScvIdx() const
186 {
187 return scvIndices_[0];
188 }
189
191 GridIndexType outsideScvIdx() const
192 {
193 return scvIndices_[1];
194 }
195
197 GridIndexType index() const
198 {
199 return scvfIndex_;
200 }
201
203 GridIndexType dofIndex() const
204 {
205 return dofIdx_;
206 }
207
209 LocalIndexType localFaceIdx() const
210 {
211 return localFaceIdx_;
212 }
213
214private:
215 Scalar area_;
216 GlobalPosition center_;
217 GlobalPosition unitOuterNormal_;
218 GridIndexType scvfIndex_;
219 std::vector<GridIndexType> scvIndices_;
220 bool boundary_;
221
222 GridIndexType dofIdx_;
223 LocalIndexType localFaceIdx_;
224};
225
226} // end namespace Dumux
227
228#endif
Base class for a staggered grid geometry helper.
Definition discretization/staggered/subcontrolvolumeface.hh:33
void updateLocalFace(const IntersectionMapper &intersectionMapper, const Intersection &intersection)
Updates the current face, i.e. sets the correct intersection.
Definition discretization/staggered/subcontrolvolumeface.hh:49
BaseStaggeredGeometryHelper(const Element &element, const GridView &gridView)
Definition discretization/staggered/subcontrolvolumeface.hh:40
int dofIndex() const
Returns the global dofIdx of the intersection itself.
Definition discretization/staggered/subcontrolvolumeface.hh:57
int localFaceIndex() const
Returns the local index of the face (i.e. the intersection)
Definition discretization/staggered/subcontrolvolumeface.hh:67
defines an intersection mapper for mapping of global DOFs assigned to faces which also works for adap...
Definition intersectionmapper.hh:200
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:108
const GlobalPosition & ipGlobal() const
The integration point for flux evaluations in global coordinates.
Definition discretization/staggered/subcontrolvolumeface.hh:160
GridIndexType index() const
The global index of this sub control volume face.
Definition discretization/staggered/subcontrolvolumeface.hh:197
LocalIndexType localFaceIdx() const
The local index of this sub control volume face.
Definition discretization/staggered/subcontrolvolumeface.hh:209
Scalar area() const
The area of the sub control volume face.
Definition discretization/staggered/subcontrolvolumeface.hh:167
T Traits
state the traits public and thus export all types
Definition discretization/staggered/subcontrolvolumeface.hh:123
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:130
GridIndexType insideScvIdx() const
Index of the inside sub control volume.
Definition discretization/staggered/subcontrolvolumeface.hh:185
const GlobalPosition & dofPosition() const
The position of the dof living on the face.
Definition discretization/staggered/subcontrolvolumeface.hh:154
const GlobalPosition & unitOuterNormal() const
The unit outer normal vector.
Definition discretization/staggered/subcontrolvolumeface.hh:179
const GlobalPosition & center() const
The center of the sub control volume face.
Definition discretization/staggered/subcontrolvolumeface.hh:148
bool boundary() const
Returns boolean if the sub control volume face is on the boundary.
Definition discretization/staggered/subcontrolvolumeface.hh:173
typename T::GlobalPosition GlobalPosition
Definition discretization/staggered/subcontrolvolumeface.hh:120
GridIndexType dofIndex() const
The global index of the dof living on this face.
Definition discretization/staggered/subcontrolvolumeface.hh:203
GridIndexType outsideScvIdx() const
Index of the outside sub control volume.
Definition discretization/staggered/subcontrolvolumeface.hh:191
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:29
auto volume(const Geometry &geo, unsigned int integrationOrder=4)
The volume of a given geometry.
Definition volume.hh:159
Defines the index types used for grid and local indices.
Definition adapt.hh:17
typename GridView::IndexSet::IndexType GridIndex
Definition indextraits.hh:27
unsigned int LocalIndex
Definition indextraits.hh:28
Default traits class to be used for the sub-control volume faces for the staggered finite volume sche...
Definition discretization/staggered/subcontrolvolumeface.hh:87
typename GridView::template Codim< 0 >::Entity Element
Definition discretization/staggered/subcontrolvolumeface.hh:91
static constexpr int dim
Definition discretization/staggered/subcontrolvolumeface.hh:94
typename IndexTraits< GridView >::GridIndex GridIndexType
Definition discretization/staggered/subcontrolvolumeface.hh:88
typename GridView::ctype Scalar
Definition discretization/staggered/subcontrolvolumeface.hh:90
typename IndexTraits< GridView >::LocalIndex LocalIndexType
Definition discretization/staggered/subcontrolvolumeface.hh:89
Dune::AxisAlignedCubeGeometry< Scalar, dim-1, dimWorld > Geometry
Definition discretization/staggered/subcontrolvolumeface.hh:96
static constexpr int dimWorld
Definition discretization/staggered/subcontrolvolumeface.hh:95
typename Element::Geometry::GlobalCoordinate GlobalPosition
Definition discretization/staggered/subcontrolvolumeface.hh:92
Base class for a sub control volume face.