Line | Branch | Exec | Source |
---|---|---|---|
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 | // | ||
7 | /*! | ||
8 | * \file | ||
9 | * \ingroup BoxDFMModel | ||
10 | * \brief the sub control volume for the box discrete fracture scheme | ||
11 | */ | ||
12 | |||
13 | #ifndef DUMUX_POROUSMEDIUMFLOW_BOXDFM_SUBCONTROLVOLUME_HH | ||
14 | #define DUMUX_POROUSMEDIUMFLOW_BOXDFM_SUBCONTROLVOLUME_HH | ||
15 | |||
16 | #include <dune/common/reservedvector.hh> | ||
17 | #include <dune/geometry/type.hh> | ||
18 | #include <dune/geometry/multilineargeometry.hh> | ||
19 | |||
20 | #include <dumux/discretization/subcontrolvolumebase.hh> | ||
21 | #include <dumux/porousmediumflow/boxdfm/geometryhelper.hh> | ||
22 | #include <dumux/common/math.hh> | ||
23 | #include <dumux/geometry/volume.hh> | ||
24 | |||
25 | namespace Dumux { | ||
26 | |||
27 | /*! | ||
28 | * \ingroup BoxDFMModel | ||
29 | * \brief Default traits class to be used for the sub-control volumes | ||
30 | * for the box discrete fracture scheme | ||
31 | * | ||
32 | * \tparam GV the type of the grid view | ||
33 | * | ||
34 | * \note We define new traits for the box-dfm sub-control volume face | ||
35 | * as we use a different type of container for storing the scvf corners! | ||
36 | */ | ||
37 | template<class GridView> | ||
38 | struct BoxDfmDefaultScvGeometryTraits | ||
39 | { | ||
40 | using Grid = typename GridView::Grid; | ||
41 | |||
42 | static const int dim = Grid::dimension; | ||
43 | static const int dimWorld = Grid::dimensionworld; | ||
44 | using GridIndexType = typename Grid::LeafGridView::IndexSet::IndexType; | ||
45 | using LocalIndexType = unsigned int; | ||
46 | using Scalar = typename Grid::ctype; | ||
47 | using GeometryTraits = BoxDfmMLGeometryTraits<Scalar>; | ||
48 | using Geometry = Dune::MultiLinearGeometry<Scalar, dim, dimWorld, GeometryTraits>; | ||
49 | using CornerStorage = typename GeometryTraits::template CornerStorage<dim, dimWorld>::Type; | ||
50 | using GlobalPosition = typename CornerStorage::value_type; | ||
51 | }; | ||
52 | |||
53 | /*! | ||
54 | * \ingroup BoxDFMModel | ||
55 | * \brief the sub control volume for the box discrete fracture scheme | ||
56 | * | ||
57 | * \tparam GV the type of the grid view | ||
58 | * \tparam T the scvf geometry traits | ||
59 | */ | ||
60 | template<class GV, | ||
61 | class T = BoxDfmDefaultScvGeometryTraits<GV> > | ||
62 | class BoxDfmSubControlVolume | ||
63 | : public SubControlVolumeBase<BoxDfmSubControlVolume<GV, T>, T> | ||
64 | { | ||
65 | using ThisType = BoxDfmSubControlVolume<GV, T>; | ||
66 | using ParentType = SubControlVolumeBase<ThisType, T>; | ||
67 | using Geometry = typename T::Geometry; | ||
68 | using GridIndexType = typename T::GridIndexType; | ||
69 | using LocalIndexType = typename T::LocalIndexType; | ||
70 | using Scalar = typename T::Scalar; | ||
71 | using GlobalPosition = typename T::GlobalPosition; | ||
72 | using CornerStorage = typename T::CornerStorage; | ||
73 | enum { dim = Geometry::mydimension }; | ||
74 | |||
75 | static_assert(dim == 2 || dim == 3, "Box-Dfm sub-control volume only implemented in 2d or 3d"); | ||
76 | |||
77 | public: | ||
78 | //! State the traits public and thus export all types | ||
79 | using Traits = T; | ||
80 | |||
81 | //! The default constructor | ||
82 | 26094996 | BoxDfmSubControlVolume() = default; | |
83 | |||
84 | // the constructor for standard scvs | ||
85 | template<class GeometryHelper> | ||
86 | 8761392 | BoxDfmSubControlVolume(const GeometryHelper& geometryHelper, | |
87 | LocalIndexType scvIdx, | ||
88 | GridIndexType elementIndex, | ||
89 | GridIndexType dofIndex) | ||
90 | : isFractureScv_(false) | ||
91 | , corners_(geometryHelper.getScvCorners(scvIdx)) | ||
92 | , center_(0.0) | ||
93 | 17522784 | , volume_(Dumux::convexPolytopeVolume<T::dim>( | |
94 | Dune::GeometryTypes::cube(T::dim), | ||
95 | 92005056 | [&](unsigned int i){ return corners_[i]; }) | |
96 | ) | ||
97 | , elementIndex_(elementIndex) | ||
98 | , vIdxLocal_(scvIdx) | ||
99 | , elemLocalScvIdx_(scvIdx) | ||
100 | , dofIndex_(dofIndex) | ||
101 | 17522784 | , facetIdx_(0) | |
102 | { | ||
103 | // compute center point | ||
104 |
2/2✓ Branch 0 taken 39428352 times.
✓ Branch 1 taken 8761392 times.
|
65712528 | for (const auto& corner : corners_) |
105 | 39428352 | center_ += corner; | |
106 | 8761392 | center_ /= corners_.size(); | |
107 | 8761392 | } | |
108 | |||
109 | /*! | ||
110 | * \brief Constructor for fracture scvs | ||
111 | * | ||
112 | * The corner computation is the same as for boundary scvfs. | ||
113 | * Also, the scvf area of a boundary scvf is equal to the scv | ||
114 | * volume (unscaled by the aperture) Thus, we reuse functionality here. | ||
115 | * In order to get the right dimensions later, one must provide appropriate | ||
116 | * extrusion factors in the problem corresponding to the fracture aperture. * | ||
117 | */ | ||
118 | template<class GeometryHelper, class Intersection> | ||
119 | 600528 | BoxDfmSubControlVolume(const GeometryHelper& geometryHelper, | |
120 | const Intersection& intersection, | ||
121 | const typename Intersection::Geometry& isGeometry, | ||
122 | LocalIndexType indexInIntersection, | ||
123 | LocalIndexType vIdxLocal, | ||
124 | LocalIndexType elemLocalScvIdx, | ||
125 | LocalIndexType elemLocalFacetIdx, | ||
126 | GridIndexType elementIndex, | ||
127 | GridIndexType dofIndex) | ||
128 | : isFractureScv_(true) | ||
129 | , corners_() | ||
130 | , center_(0.0) | ||
131 | , volume_(0.0) | ||
132 | , elementIndex_(elementIndex) | ||
133 | , vIdxLocal_(vIdxLocal) | ||
134 | , elemLocalScvIdx_(elemLocalScvIdx) | ||
135 | , dofIndex_(dofIndex) | ||
136 | 1201056 | , facetIdx_(elemLocalFacetIdx) | |
137 | { | ||
138 | // copy corners | ||
139 |
2/2✓ Branch 0 taken 4284 times.
✓ Branch 1 taken 254400 times.
|
900792 | auto corners = geometryHelper.getBoundaryScvfCorners(intersection.indexInInside(), indexInIntersection); |
140 | 600528 | const auto numCorners = corners.size(); | |
141 | 600528 | corners_.resize(numCorners); | |
142 |
2/2✓ Branch 0 taken 1367376 times.
✓ Branch 1 taken 600528 times.
|
1967904 | for (unsigned int i = 0; i < numCorners; ++i) |
143 | 4102128 | corners_[i] = corners[i]; | |
144 | |||
145 | // compute volume and scv center | ||
146 | 600528 | volume_ = Dumux::convexPolytopeVolume<T::dim-1>( | |
147 | Dune::GeometryTypes::cube(T::dim-1), | ||
148 | 665280 | [&](unsigned int i){ return corners_[i]; } | |
149 | ); | ||
150 |
2/2✓ Branch 0 taken 1367376 times.
✓ Branch 1 taken 600528 times.
|
3168960 | for (const auto& corner : corners_) |
151 | 1367376 | center_ += corner; | |
152 | 600528 | center_ /= corners_.size(); | |
153 | 600528 | } | |
154 | |||
155 | //! The center of the sub control volume | ||
156 | const GlobalPosition& center() const | ||
157 | 132538736 | { return center_; } | |
158 | |||
159 | //! The volume of the sub control volume | ||
160 | ✗ | Scalar volume() const | |
161 | ✗ | { return volume_; } | |
162 | |||
163 | //! The element-local vertex index this scv is connected to | ||
164 | ✗ | LocalIndexType localDofIndex() const | |
165 | ✗ | { return vIdxLocal_; } | |
166 | |||
167 | //! The element-local index of this scv | ||
168 | ✗ | LocalIndexType indexInElement() const | |
169 | ✗ | { return elemLocalScvIdx_; } | |
170 | |||
171 | //! The element-local facet index for which a fracture scv was created | ||
172 | ✗ | LocalIndexType facetIndexInElement() const | |
173 |
6/8✓ Branch 0 taken 648016 times.
✓ Branch 1 taken 13656 times.
✓ Branch 2 taken 648016 times.
✓ Branch 3 taken 13656 times.
✓ Branch 4 taken 13656 times.
✓ Branch 5 taken 648016 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
2057224 | { assert(isFractureScv_); return facetIdx_; } |
174 | |||
175 | //! The index of the dof this scv is embedded in | ||
176 | ✗ | GridIndexType dofIndex() const | |
177 | ✗ | { return dofIndex_; } | |
178 | |||
179 | // The position of the dof this scv is embedded in (list is defined such that first entry is vertex itself) | ||
180 | const GlobalPosition& dofPosition() const | ||
181 |
4/12✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 159266 times.
✓ Branch 5 taken 121716 times.
✓ Branch 6 taken 159266 times.
✓ Branch 7 taken 121716 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
1928520 | { return corners_[0]; } |
182 | |||
183 | //! The global index of the element this scv is embedded in | ||
184 | GridIndexType elementIndex() const | ||
185 | { return elementIndex_; } | ||
186 | |||
187 | //! Return true if this scv is part of the fracture domain | ||
188 | ✗ | bool isOnFracture() const | |
189 | ✗ | { return isFractureScv_; } | |
190 | |||
191 | //! The geometry of the sub control volume | ||
192 | // e.g. for integration | ||
193 | [[deprecated("Will be removed after 3.7. Use fvGeometry.geometry(scv).")]] | ||
194 | Geometry geometry() const | ||
195 | { | ||
196 | if (isFractureScv_) | ||
197 | DUNE_THROW(Dune::InvalidStateException, "The geometry object cannot be defined for fracture scvs " | ||
198 | "because the number of known corners is insufficient. " | ||
199 | "You can do this manually by extract the corners from this scv " | ||
200 | "and extrude them by the corresponding aperture. "); | ||
201 | |||
202 | return Geometry(Dune::GeometryTypes::cube(dim), corners_); | ||
203 | } | ||
204 | |||
205 | private: | ||
206 | bool isFractureScv_; | ||
207 | CornerStorage corners_; | ||
208 | GlobalPosition center_; | ||
209 | Scalar volume_; | ||
210 | GridIndexType elementIndex_; | ||
211 | LocalIndexType vIdxLocal_; | ||
212 | LocalIndexType elemLocalScvIdx_; | ||
213 | GridIndexType dofIndex_; | ||
214 | |||
215 | // for fracture scvs only! | ||
216 | LocalIndexType facetIdx_; | ||
217 | }; | ||
218 | |||
219 | } // end namespace | ||
220 | |||
221 | #endif | ||
222 |