3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
discretization/cellcentered/mpfa/omethod/interactionvolume.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_CC_MPFA_O_INTERACTIONVOLUME_HH
25#define DUMUX_DISCRETIZATION_CC_MPFA_O_INTERACTIONVOLUME_HH
26
27#include <type_traits>
28
29#include <dune/common/dynmatrix.hh>
30#include <dune/common/dynvector.hh>
31#include <dune/common/fvector.hh>
32#include <dune/common/reservedvector.hh>
33
37
38#include "localassembler.hh"
41#include "scvgeometryhelper.hh"
42
43namespace Dumux {
44
46template< class Traits > class CCMpfaOInteractionVolume;
47
58template< class NodalIndexSet, class Scalar >
60{
61private:
62 using GridIndexType = typename NodalIndexSet::GridIndexType;
63 using LocalIndexType = typename NodalIndexSet::LocalIndexType;
64
65 static constexpr int dim = NodalIndexSet::Traits::GridView::dimension;
66 static constexpr int dimWorld = NodalIndexSet::Traits::GridView::dimensionworld;
67
68 using DimVector = Dune::FieldVector<Scalar, dim>;
69 using FaceOmegas = typename std::conditional< (dim<dimWorld),
70 std::vector<DimVector>,
71 Dune::ReservedVector<DimVector, 2> >::type;
72
74 struct MVTraits
75 {
76 using OmegaStorage = std::vector< FaceOmegas >;
77
78 using AMatrix = Dune::DynamicMatrix< Scalar >;
79 using BMatrix = Dune::DynamicMatrix< Scalar >;
80 using CMatrix = Dune::DynamicMatrix< Scalar >;
81 using DMatrix = Dune::DynamicMatrix< Scalar >;
82 using TMatrix = Dune::DynamicMatrix< Scalar >;
83 using CellVector = Dune::DynamicVector< Scalar >;
84 using FaceVector = Dune::DynamicVector< Scalar >;
85 };
86
87public:
89 using GridView = typename NodalIndexSet::Traits::GridView;
99 using MatVecTraits = MVTraits;
100
102 template<class Problem, class FVElementGeometry, class ElemVolVars>
104};
105
112template< class Traits >
114: public CCMpfaInteractionVolumeBase< Traits >
115{
116 using GridView = typename Traits::GridView;
117 using Element = typename GridView::template Codim<0>::Entity;
118
119 using IndexSet = typename Traits::IndexSet;
120 using GridIndexType = typename IndexSet::GridIndexType;
121 using LocalIndexType = typename IndexSet::LocalIndexType;
122 using Stencil = typename IndexSet::NodalGridStencilType;
123
124 using LocalScvType = typename Traits::LocalScvType;
125 using LocalScvfType = typename Traits::LocalScvfType;
126 using LocalFaceData = typename Traits::LocalFaceData;
127
128public:
132 {
133 GridIndexType volVarIndex_;
134 public:
136 DirichletData(const GridIndexType index) : volVarIndex_(index) {}
137
139 GridIndexType volVarIndex() const { return volVarIndex_; }
140 };
141
144
146 template< class Problem, class FVElementGeometry >
147 void bind(const IndexSet& indexSet,
148 const Problem& problem,
149 const FVElementGeometry& fvGeometry)
150 {
151 // for the o-scheme, the stencil is equal to the scv
152 // index set of the dual grid's nodal index set
153 stencil_ = &indexSet.nodalIndexSet().gridScvIndices();
154
155 // number of interaction-volume-local scvs(=node-local for o-scheme) and scvfs
156 numFaces_ = indexSet.numFaces();
157 const auto numLocalScvs = indexSet.numScvs();
158 const auto numGlobalScvfs = indexSet.nodalIndexSet().numScvfs();
159
160 // reserve memory for local entities
161 elements_.clear(); elements_.reserve(numLocalScvs);
162 scvs_.clear(); scvs_.reserve(numLocalScvs);
163 scvfs_.clear(); scvfs_.reserve(numFaces_);
164 localFaceData_.clear(); localFaceData_.reserve(numGlobalScvfs);
165 dirichletData_.clear(); dirichletData_.reserve(numFaces_);
166
167 // set up stuff related to sub-control volumes
168 for (LocalIndexType scvIdxLocal = 0; scvIdxLocal < numLocalScvs; scvIdxLocal++)
169 {
170 elements_.emplace_back(fvGeometry.gridGeometry().element( stencil()[scvIdxLocal] ));
171 scvs_.emplace_back(fvGeometry.gridGeometry().mpfaHelper(),
172 fvGeometry,
173 fvGeometry.scv( stencil()[scvIdxLocal] ),
174 scvIdxLocal,
175 indexSet);
176 }
177
178 // keep track of the number of unknowns etc
179 numUnknowns_ = 0;
180 numKnowns_ = numLocalScvs;
181
182 // set up quantitites related to sub-control volume faces
183 for (LocalIndexType faceIdxLocal = 0; faceIdxLocal < numFaces_; ++faceIdxLocal)
184 {
185 const auto& scvf = fvGeometry.scvf(indexSet.gridScvfIndex(faceIdxLocal));
186
187 // the neighboring scvs in local indices (order: 0 - inside scv, 1..n - outside scvs)
188 const auto& neighborScvIndicesLocal = indexSet.neighboringLocalScvIndices(faceIdxLocal);
189 const auto numNeighborScvs = neighborScvIndicesLocal.size();
190 localFaceData_.emplace_back(faceIdxLocal, neighborScvIndicesLocal[0], scvf.index());
191
192 // create iv-local scvf object
193 if (scvf.boundary())
194 {
195 if (problem.boundaryTypes(elements_[neighborScvIndicesLocal[0]], scvf).hasOnlyDirichlet())
196 {
197 scvfs_.emplace_back(scvf, neighborScvIndicesLocal, numKnowns_++, /*isDirichlet*/true);
198 dirichletData_.emplace_back(scvf.outsideScvIdx());
199 }
200 else
201 scvfs_.emplace_back(scvf, neighborScvIndicesLocal, numUnknowns_++, /*isDirichlet*/false);
202 }
203 else
204 {
205 scvfs_.emplace_back(scvf, neighborScvIndicesLocal, numUnknowns_++, /*isDirichlet*/false);
206
207 // add local face data objects for the outside faces
208 for (LocalIndexType i = 1; i < numNeighborScvs; ++i)
209 {
210 // loop over scvfs in outside scv until we find the one coinciding with current scvf
211 const auto outsideLocalScvIdx = neighborScvIndicesLocal[i];
212 const auto& flipScvfIndex = fvGeometry.gridGeometry().flipScvfIndexSet()[scvf.index()][i-1];
213 const auto& flipScvf = fvGeometry.scvf(flipScvfIndex);
214 localFaceData_.emplace_back(faceIdxLocal, // iv-local scvf idx
215 outsideLocalScvIdx, // iv-local scv index
216 i-1, // scvf-local index in outside faces
217 flipScvf.index()); // global scvf index
218 }
219 }
220 }
221 }
222
224 std::size_t numFaces() const
225 { return numFaces_; }
226
228 std::size_t numUnknowns() const
229 { return numUnknowns_; }
230
232 std::size_t numKnowns() const
233 { return numKnowns_; }
234
236 std::size_t numScvs() const
237 { return scvs_.size(); }
238
240 const Stencil& stencil() const
241 { return *stencil_; }
242
244 const Element& element(LocalIndexType ivLocalScvIdx) const
245 { return elements_[ivLocalScvIdx]; }
246
248 const LocalScvfType& localScvf(LocalIndexType ivLocalScvfIdx) const
249 { return scvfs_[ivLocalScvfIdx]; }
250
252 const LocalScvType& localScv(LocalIndexType ivLocalScvIdx) const
253 { return scvs_[ivLocalScvIdx]; }
254
256 const std::vector<LocalFaceData>& localFaceData() const
257 { return localFaceData_; }
258
260 const std::vector<DirichletData>& dirichletData() const
261 { return dirichletData_; }
262
264 template< class FVElementGeometry >
265 auto getScvGeometry(LocalIndexType ivLocalScvIdx, const FVElementGeometry& fvGeometry) const
266 { return CCMpfaOScvGeometryHelper<LocalScvType>::computeScvGeometry(ivLocalScvIdx, *this, fvGeometry); }
267
269 template< class NI >
270 static constexpr std::size_t numIVAtVertex(const NI& nodalIndexSet)
271 { return 1; }
272
275 template< class IvIndexSetContainer,
276 class ScvfIndexMap,
277 class NodalIndexSet,
278 class FlipScvfIndexSet >
279 static void addIVIndexSets(IvIndexSetContainer& ivIndexSetContainer,
280 ScvfIndexMap& scvfIndexMap,
281 const NodalIndexSet& nodalIndexSet,
282 const FlipScvfIndexSet& flipScvfIndexSet)
283 {
284 // the global index of the iv index set that is about to be created
285 const auto curGlobalIndex = ivIndexSetContainer.size();
286
287 // make the one index set for this node
288 ivIndexSetContainer.emplace_back(nodalIndexSet, flipScvfIndexSet);
289
290 // store the index mapping
291 for (const auto scvfIdx : nodalIndexSet.gridScvfIndices())
292 scvfIndexMap[scvfIdx] = curGlobalIndex;
293 }
294
295private:
296 // pointer to cell stencil (in iv index set)
297 const Stencil* stencil_;
298
299 // Variables defining the local scope
300 std::vector<Element> elements_;
301 std::vector<LocalScvType> scvs_;
302 std::vector<LocalScvfType> scvfs_;
303 std::vector<LocalFaceData> localFaceData_;
304 std::vector<DirichletData> dirichletData_;
305
306 // sizes involved in the local system equations
307 std::size_t numFaces_;
308 std::size_t numUnknowns_;
309 std::size_t numKnowns_;
310};
311
312} // end namespace Dumux
313
314#endif
Class for the index set within an interaction volume of the mpfa-o scheme.
Helper class providing functionality to compute the geometry of the interaction-volume local sub-cont...
Base class for interaction volumes of mpfa methods.
Data structure holding interaction volume-local information for a grid subb-control volume face embed...
The available mpfa schemes in Dumux.
MpfaMethods
The available mpfa schemes in Dumux.
Definition: methods.hh:34
Definition: adapt.hh:29
Base class for the interaction volumes of mpfa methods. It defines the interface and actual implement...
Definition: interactionvolumebase.hh:68
General implementation of a data structure holding interaction volume-local information for a grid su...
Definition: localfacedata.hh:42
Forward declaration of the o-method's interaction volume.
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:115
const LocalScvfType & localScvf(LocalIndexType ivLocalScvfIdx) const
returns the local scvf entity corresponding to a given iv-local scvf idx
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:248
std::size_t numUnknowns() const
returns the number of intermediate unknowns within this interaction volume
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:228
static constexpr std::size_t numIVAtVertex(const NI &nodalIndexSet)
returns the number of interaction volumes living around a vertex
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:270
void bind(const IndexSet &indexSet, const Problem &problem, const FVElementGeometry &fvGeometry)
Sets up the local scope for a given iv index set.
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:147
static constexpr MpfaMethods MpfaMethod
publicly state the mpfa-scheme this interaction volume is associated with
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:143
static void addIVIndexSets(IvIndexSetContainer &ivIndexSetContainer, ScvfIndexMap &scvfIndexMap, const NodalIndexSet &nodalIndexSet, const FlipScvfIndexSet &flipScvfIndexSet)
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:279
std::size_t numKnowns() const
returns the number of (in this context) known solution values within this interaction volume
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:232
const Element & element(LocalIndexType ivLocalScvIdx) const
returns the grid element corresponding to a given iv-local scv idx
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:244
const LocalScvType & localScv(LocalIndexType ivLocalScvIdx) const
returns the local scv entity corresponding to a given iv-local scv idx
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:252
const std::vector< LocalFaceData > & localFaceData() const
returns a reference to the container with the local face data
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:256
const Stencil & stencil() const
returns the cell-stencil of this interaction volume
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:240
const std::vector< DirichletData > & dirichletData() const
returns a reference to the information container on Dirichlet BCs within this iv
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:260
std::size_t numScvs() const
returns the number of scvs embedded in this interaction volume
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:236
std::size_t numFaces() const
returns the number of primary scvfs of this interaction volume
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:224
auto getScvGeometry(LocalIndexType ivLocalScvIdx, const FVElementGeometry &fvGeometry) const
returns the geometry of the i-th local scv
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:265
The default interaction volume traits class for the mpfa-o method. This uses dynamic types types for ...
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:60
MVTraits MatVecTraits
export the matrix/vector traits to be used by the iv
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:99
typename NodalIndexSet::Traits::GridView GridView
export the type of grid view
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:89
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:132
GridIndexType volVarIndex() const
Return corresponding vol var index.
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:139
DirichletData(const GridIndexType index)
Constructor.
Definition: discretization/cellcentered/mpfa/omethod/interactionvolume.hh:136
The interaction volume index set class for the mpfa-o scheme.
Definition: interactionvolumeindexset.hh:41
Specialization of the interaction volume-local assembler class for the schemes using an mpfa-o type a...
Definition: discretization/cellcentered/mpfa/omethod/localassembler.hh:48
Class for the interaction volume-local sub-control volume used in the mpfa-o scheme.
Definition: discretization/cellcentered/mpfa/omethod/localsubcontrolentities.hh:43
Class for the interaction volume-local sub-control volume face used in the mpfa-o scheme.
Definition: discretization/cellcentered/mpfa/omethod/localsubcontrolentities.hh:138
static ScvGeometry computeScvGeometry(LocalIndexType ivLocalScvIdx, const InteractionVolume &iv, const FVElementGeometry &fvGeometry)
returns the geometry of the i-th local scv
Definition: scvgeometryhelper.hh:73
Classes for sub control entities of the mpfa-o method in the context of facet coupling.