3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
porousmediumflow/boxdfm/subcontrolvolume.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 *****************************************************************************/
25#ifndef DUMUX_POROUSMEDIUMFLOW_BOXDFM_SUBCONTROLVOLUME_HH
26#define DUMUX_POROUSMEDIUMFLOW_BOXDFM_SUBCONTROLVOLUME_HH
27
28#include <dune/common/reservedvector.hh>
29#include <dune/geometry/type.hh>
30#include <dune/geometry/multilineargeometry.hh>
31
34#include <dumux/common/math.hh>
35
36namespace Dumux {
37
48template<class GridView>
50{
51 using Grid = typename GridView::Grid;
52
53 static const int dim = Grid::dimension;
54 static const int dimWorld = Grid::dimensionworld;
55
56 template <class ct>
57 struct ScvMLGTraits : public Dune::MultiLinearGeometryTraits<ct>
58 {
59 // we use static vectors to store the corners as we know
60 // the number of corners in advance (2^(dim) corners (1<<(dim))
61 // However, on fracture scvs the number might be smaller (use ReservedVector)
62 template< int mydim, int cdim >
64 {
65 using Type = Dune::ReservedVector< Dune::FieldVector< ct, cdim >, (1<<(dim)) >;
66 };
67
68 // we know all scvfs will have the same geometry type
69 template< int mydim >
71 {
72 static const bool v = true;
73 static const unsigned int topologyId = Dune::GeometryTypes::cube(mydim).id();
74 };
75 };
76
77 using GridIndexType = typename Grid::LeafGridView::IndexSet::IndexType;
78 using LocalIndexType = unsigned int;
79 using Scalar = typename Grid::ctype;
80 using Geometry = Dune::MultiLinearGeometry<Scalar, dim, dimWorld, ScvMLGTraits<Scalar>>;
82 using GlobalPosition = typename CornerStorage::value_type;
83};
84
92template<class GV,
95: public SubControlVolumeBase<BoxDfmSubControlVolume<GV, T>, T>
96{
99 using Geometry = typename T::Geometry;
100 using GridIndexType = typename T::GridIndexType;
101 using LocalIndexType = typename T::LocalIndexType;
102 using Scalar = typename T::Scalar;
103 using GlobalPosition = typename T::GlobalPosition;
104 using CornerStorage = typename T::CornerStorage;
105 enum { dim = Geometry::mydimension };
106
107 static_assert(dim == 2 || dim == 3, "Box-Dfm sub-control volume only implemented in 2d or 3d");
108
109public:
111 using Traits = T;
112
115
116 // the constructor for standard scvs
117 template<class GeometryHelper>
118 BoxDfmSubControlVolume(const GeometryHelper& geometryHelper,
119 LocalIndexType scvIdx,
120 GridIndexType elementIndex,
121 GridIndexType dofIndex)
122 : isFractureScv_(false)
123 , corners_(geometryHelper.getScvCorners(scvIdx))
124 , center_(0.0)
125 , volume_(geometryHelper.scvVolume(corners_))
126 , elementIndex_(elementIndex)
127 , vIdxLocal_(scvIdx)
128 , elemLocalScvIdx_(scvIdx)
129 , dofIndex_(dofIndex)
130 , facetIdx_(0)
131 {
132 // compute center point
133 for (const auto& corner : corners_)
134 center_ += corner;
135 center_ /= corners_.size();
136 }
137
147 template<class GeometryHelper, class Intersection>
148 BoxDfmSubControlVolume(const GeometryHelper& geometryHelper,
149 const Intersection& intersection,
150 const typename Intersection::Geometry& isGeometry,
151 LocalIndexType indexInIntersection,
152 LocalIndexType vIdxLocal,
153 LocalIndexType elemLocalScvIdx,
154 LocalIndexType elemLocalFacetIdx,
155 GridIndexType elementIndex,
156 GridIndexType dofIndex)
157 : isFractureScv_(true)
158 , corners_()
159 , center_(0.0)
160 , volume_(0.0)
161 , elementIndex_(elementIndex)
162 , vIdxLocal_(vIdxLocal)
163 , elemLocalScvIdx_(elemLocalScvIdx)
164 , dofIndex_(dofIndex)
165 , facetIdx_(elemLocalFacetIdx)
166 {
167 // copy corners
168 auto corners = geometryHelper.getBoundaryScvfCorners(intersection, isGeometry, indexInIntersection);
169 const auto numCorners = corners.size();
170 corners_.resize(numCorners);
171 for (unsigned int i = 0; i < numCorners; ++i)
172 corners_[i] = corners[i];
173
174 // compute volume and scv center
175 volume_ = geometryHelper.scvfArea(corners);
176 for (const auto& corner : corners_)
177 center_ += corner;
178 center_ /= corners_.size();
179 }
180
182 const GlobalPosition& center() const
183 { return center_; }
184
186 Scalar volume() const
187 { return volume_; }
188
190 LocalIndexType localDofIndex() const
191 { return vIdxLocal_; }
192
194 LocalIndexType indexInElement() const
195 { return elemLocalScvIdx_; }
196
198 LocalIndexType facetIndexInElement() const
199 { assert(isFractureScv_); return facetIdx_; }
200
202 GridIndexType dofIndex() const
203 { return dofIndex_; }
204
205 // The position of the dof this scv is embedded in (list is defined such that first entry is vertex itself)
206 const GlobalPosition& dofPosition() const
207 { return corners_[0]; }
208
210 GridIndexType elementIndex() const
211 { return elementIndex_; }
212
214 bool isOnFracture() const
215 { return isFractureScv_; }
216
218 // e.g. for integration
219 Geometry geometry() const
220 {
221 if (isFractureScv_)
222 DUNE_THROW(Dune::InvalidStateException, "The geometry object cannot be defined for fracture scvs "
223 "because the number of known corners is insufficient. "
224 "You can do this manually by extract the corners from this scv "
225 "and extrude them by the corresponding aperture. ");
226
227 return Geometry(Dune::GeometryTypes::cube(dim), corners_);
228 }
229
231 const GlobalPosition& corner(LocalIndexType localIdx) const
232 {
233 assert(localIdx < corners_.size() && "provided index exceeds the number of corners");
234 return corners_[localIdx];
235 }
236
237private:
238 bool isFractureScv_;
239 CornerStorage corners_;
240 GlobalPosition center_;
241 Scalar volume_;
242 GridIndexType elementIndex_;
243 LocalIndexType vIdxLocal_;
244 LocalIndexType elemLocalScvIdx_;
245 GridIndexType dofIndex_;
246
247 // for fracture scvs only!
248 LocalIndexType facetIdx_;
249};
250
251} // end namespace
252
253#endif
Define some often used mathematical functions.
Helper class constructing the dual grid finite volume geometries for the box discretizazion method.
Base class for a sub control volume.
Definition: adapt.hh:29
std::size_t numCorners(Shape shape)
Returns the number of corners of a given geometry.
Definition: throatproperties.hh:215
Base class for a sub control volume, i.e a part of the control volume we are making the balance for....
Definition: subcontrolvolumebase.hh:38
Default traits class to be used for the sub-control volumes for the box discrete fracture scheme.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:50
static const int dimWorld
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:54
static const int dim
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:53
typename GridView::Grid Grid
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:51
Dune::MultiLinearGeometry< Scalar, dim, dimWorld, ScvMLGTraits< Scalar > > Geometry
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:80
typename Grid::LeafGridView::IndexSet::IndexType GridIndexType
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:77
typename ScvMLGTraits< Scalar >::template CornerStorage< dim, dimWorld >::Type CornerStorage
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:81
typename CornerStorage::value_type GlobalPosition
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:82
unsigned int LocalIndexType
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:78
typename Grid::ctype Scalar
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:79
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:58
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:64
Dune::ReservedVector< Dune::FieldVector< ct, cdim >,(1<<(dim)) > Type
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:65
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:71
static const unsigned int topologyId
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:73
static const bool v
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:72
the sub control volume for the box discrete fracture scheme
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:96
BoxDfmSubControlVolume()=default
The default constructor.
const GlobalPosition & center() const
The center of the sub control volume.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:182
LocalIndexType facetIndexInElement() const
The element-local facet index for which a fracture scv was created.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:198
LocalIndexType indexInElement() const
The element-local index of this scv.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:194
Geometry geometry() const
The geometry of the sub control volume.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:219
bool isOnFracture() const
Return true if this scv is part of the fracture domain.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:214
const GlobalPosition & corner(LocalIndexType localIdx) const
Return the corner for the given local index.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:231
Scalar volume() const
The volume of the sub control volume.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:186
const GlobalPosition & dofPosition() const
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:206
T Traits
State the traits public and thus export all types.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:111
GridIndexType elementIndex() const
The global index of the element this scv is embedded in.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:210
BoxDfmSubControlVolume(const GeometryHelper &geometryHelper, const Intersection &intersection, const typename Intersection::Geometry &isGeometry, LocalIndexType indexInIntersection, LocalIndexType vIdxLocal, LocalIndexType elemLocalScvIdx, LocalIndexType elemLocalFacetIdx, GridIndexType elementIndex, GridIndexType dofIndex)
Constructor for fracture scvs.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:148
LocalIndexType localDofIndex() const
The element-local vertex index this scv is connected to.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:190
BoxDfmSubControlVolume(const GeometryHelper &geometryHelper, LocalIndexType scvIdx, GridIndexType elementIndex, GridIndexType dofIndex)
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:118
GridIndexType dofIndex() const
The index of the dof this scv is embedded in.
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:202