3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
porousmediumflow/boxdfm/fvgridgeometry.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 *****************************************************************************/
29#ifndef DUMUX_POROUSMEDIUMFLOW_BOXDFM_GRID_FVGEOMETRY_HH
30#define DUMUX_POROUSMEDIUMFLOW_BOXDFM_GRID_FVGEOMETRY_HH
31
32#include <unordered_map>
33
34#include <dune/geometry/referenceelements.hh>
35#include <dune/localfunctions/lagrange/pqkfactory.hh>
36#include <dune/geometry/multilineargeometry.hh>
37#include <dune/grid/common/mcmgmapper.hh>
38
43
44#include "fvelementgeometry.hh"
45#include "geometryhelper.hh"
46#include "subcontrolvolume.hh"
48
49namespace Dumux {
50
59template<class GridView, class MapperTraits = DefaultMapperTraits<GridView>>
61: public MapperTraits
62{
65
66 template<class GridGeometry, bool enableCache>
68
69 // Mapper type for mapping edges
70 using FacetMapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
71};
72
81template<class Scalar,
82 class GridView,
83 bool enableGridGeometryCache = false,
86
97template<class Scalar, class GV, class Traits>
98class BoxDfmFVGridGeometry<Scalar, GV, true, Traits>
99: public BaseGridGeometry<GV, Traits>
100{
103 using GridIndexType = typename GV::IndexSet::IndexType;
104
105 using Element = typename GV::template Codim<0>::Entity;
106 using CoordScalar = typename GV::ctype;
107 static const int dim = GV::dimension;
108 static const int dimWorld = GV::dimensionworld;
109 static_assert(dim == 2 || dim == 3, "The box-dfm GridGeometry is only implemented in 2 or 3 dimensions.");
110
111 using ReferenceElements = typename Dune::ReferenceElements<CoordScalar, dim>;
112 using FaceReferenceElements = typename Dune::ReferenceElements<CoordScalar, dim-1>;
113
114 using GeometryHelper = BoxDfmGeometryHelper<GV, dim,
115 typename Traits::SubControlVolume,
116 typename Traits::SubControlVolumeFace>;
117
118public:
120 static constexpr DiscretizationMethod discMethod = DiscretizationMethod::box;
121
123 using LocalView = typename Traits::template LocalView<ThisType, true>;
125 using SubControlVolume = typename Traits::SubControlVolume;
127 using SubControlVolumeFace = typename Traits::SubControlVolumeFace;
129 using DofMapper = typename Traits::VertexMapper;
131 using FeCache = Dune::PQkLocalFiniteElementCache<CoordScalar, Scalar, dim, 1>;
133 using GridView = GV;
134
137 : ParentType(gridView) {}
138
141 const DofMapper& dofMapper() const
142 { return this->vertexMapper(); }
143
145 std::size_t numScv() const
146 { return numScv_; }
147
149 std::size_t numScvf() const
150 { return numScvf_; }
151
154 std::size_t numBoundaryScvf() const
155 { return numBoundaryScvf_; }
156
158 std::size_t numDofs() const
159 { return this->gridView().size(dim); }
160
162 template< class FractureGridAdapter >
163 void update(const FractureGridAdapter& fractureGridAdapter)
164 {
165 ParentType::update();
166
167 scvs_.clear();
168 scvfs_.clear();
169
170 auto numElements = this->gridView().size(0);
171 scvs_.resize(numElements);
172 scvfs_.resize(numElements);
173
174 boundaryDofIndices_.assign(numDofs(), false);
175 fractureDofIndices_.assign(this->gridView.size(dim), false);
176
177 numScv_ = 0;
178 numScvf_ = 0;
179 numBoundaryScvf_ = 0;
180 // Build the SCV and SCV faces
181 for (const auto& element : elements(this->gridView()))
182 {
183 // fill the element map with seeds
184 auto eIdx = this->elementMapper().index(element);
185
186 // count
187 numScv_ += element.subEntities(dim);
188 numScvf_ += element.subEntities(dim-1);
189
190 // get the element geometry
191 auto elementGeometry = element.geometry();
192 const auto referenceElement = ReferenceElements::general(elementGeometry.type());
193
194 // instantiate the geometry helper
195 GeometryHelper geometryHelper(elementGeometry);
196
197 // construct the sub control volumes
198 scvs_[eIdx].resize(elementGeometry.corners());
199 using LocalIndexType = typename SubControlVolumeFace::Traits::LocalIndexType;
200 for (LocalIndexType scvLocalIdx = 0; scvLocalIdx < elementGeometry.corners(); ++scvLocalIdx)
201 {
202 const auto dofIdxGlobal = this->vertexMapper().subIndex(element, scvLocalIdx, dim);
203
204 scvs_[eIdx][scvLocalIdx] = SubControlVolume(geometryHelper,
205 scvLocalIdx,
206 eIdx,
207 dofIdxGlobal);
208 }
209
210 // construct the sub control volume faces
211 LocalIndexType scvfLocalIdx = 0;
212 scvfs_[eIdx].resize(element.subEntities(dim-1));
213 for (; scvfLocalIdx < element.subEntities(dim-1); ++scvfLocalIdx)
214 {
215 // find the global and local scv indices this scvf is belonging to
216 std::vector<LocalIndexType> localScvIndices({static_cast<LocalIndexType>(referenceElement.subEntity(scvfLocalIdx, dim-1, 0, dim)),
217 static_cast<LocalIndexType>(referenceElement.subEntity(scvfLocalIdx, dim-1, 1, dim))});
218
219 scvfs_[eIdx][scvfLocalIdx] = SubControlVolumeFace(geometryHelper,
220 element,
221 elementGeometry,
222 scvfLocalIdx,
223 std::move(localScvIndices));
224 }
225
226 // construct the ...
227 // ... sub-control volume faces on the domain boundary
228 // ... sub-control volumes on fracture facets
229 // ... sub-control volume faces on fracture facets
230 // NOTE We do not construct fracture scvfs on boundaries here!
231 // That means specifying Neumann fluxes on only the fractures is not possible
232 // However, it is difficult to interpret this here as a fracture ending on the boundary
233 // could also be connected to a facet which is both boundary and fracture at the same time!
234 // In that case, the fracture boundary scvf wouldn't make sense. In order to do it properly
235 // we would have to find only those fractures that are at the boundary and aren't connected
236 // to a fracture which is a boundary.
237 LocalIndexType scvLocalIdx = element.subEntities(dim);
238 for (const auto& intersection : intersections(this->gridView(), element))
239 {
240 // first, obtain all vertex indices on this intersection
241 const auto& isGeometry = intersection.geometry();
242 const auto numCorners = isGeometry.corners();
243 const auto idxInInside = intersection.indexInInside();
244
245 std::vector<GridIndexType> isVertexIndices(numCorners);
246 for (unsigned int vIdxLocal = 0; vIdxLocal < numCorners; ++vIdxLocal)
247 isVertexIndices[vIdxLocal] = this->vertexMapper().subIndex(element,
248 referenceElement.subEntity(idxInInside, 1, vIdxLocal, dim),
249 dim);
250 // maybe add boundary scvf
251 if (intersection.boundary() && !intersection.neighbor())
252 {
253 numScvf_ += isGeometry.corners();
254 numBoundaryScvf_ += isGeometry.corners();
255
256 for (unsigned int isScvfLocalIdx = 0; isScvfLocalIdx < numCorners; ++isScvfLocalIdx)
257 {
258 // find the scvs this scvf is belonging to
259 const LocalIndexType insideScvIdx = static_cast<LocalIndexType>(referenceElement.subEntity(idxInInside, 1, isScvfLocalIdx, dim));
260 std::vector<LocalIndexType> localScvIndices = {insideScvIdx, insideScvIdx};
261 scvfs_[eIdx].emplace_back(geometryHelper,
262 intersection,
263 isGeometry,
264 isScvfLocalIdx,
265 scvfLocalIdx++,
266 std::move(localScvIndices));
267 }
268
269 // add all vertices on the intersection to the set of
270 // boundary vertices
271 const auto numFaceVerts = referenceElement.size(idxInInside, 1, dim);
272 for (int localVIdx = 0; localVIdx < numFaceVerts; ++localVIdx)
273 {
274 const auto vIdx = referenceElement.subEntity(idxInInside, 1, localVIdx, dim);
275 const auto vIdxGlobal = this->vertexMapper().subIndex(element, vIdx, dim);
276 boundaryDofIndices_[vIdxGlobal] = true;
277 }
278 }
279 // make sure we have no periodic boundaries
280 else if (intersection.boundary() && intersection.neighbor())
281 DUNE_THROW(Dune::InvalidStateException, "Periodic boundaries are not supported by the box-dfm scheme");
282
283 // maybe add fracture scvs & scvfs
284 if (fractureGridAdapter.composeFacetElement(isVertexIndices))
285 {
286 for (auto vIdx : isVertexIndices)
287 fractureDofIndices_[vIdx] = true;
288
289 // add fracture scv for each vertex of intersection
290 numScv_ += numCorners;
291 const auto curNumScvs = scvs_[eIdx].size();
292 scvs_[eIdx].reserve(curNumScvs+numCorners);
293 for (unsigned int vIdxLocal = 0; vIdxLocal < numCorners; ++vIdxLocal)
294 scvs_[eIdx].emplace_back(geometryHelper,
295 intersection,
296 isGeometry,
297 vIdxLocal,
298 static_cast<LocalIndexType>(referenceElement.subEntity(idxInInside, 1, vIdxLocal, dim)),
299 scvLocalIdx++,
300 idxInInside,
301 eIdx,
302 isVertexIndices[vIdxLocal]);
303
304 // add fracture scvf for each edge of the intersection in 3d
305 if (dim == 3)
306 {
307 const auto& faceReferenceElement = FaceReferenceElements::general(isGeometry.type());
308 for (unsigned int edgeIdx = 0; edgeIdx < faceReferenceElement.size(1); ++edgeIdx)
309 {
310 // inside/outside scv indices in face local node numbering
311 std::vector<LocalIndexType> localScvIndices({static_cast<LocalIndexType>(faceReferenceElement.subEntity(edgeIdx, 1, 0, dim-1)),
312 static_cast<LocalIndexType>(faceReferenceElement.subEntity(edgeIdx, 1, 1, dim-1))});
313
314 // add offset to get the right scv indices
315 std::for_each( localScvIndices.begin(),
316 localScvIndices.end(),
317 [curNumScvs] (auto& elemLocalIdx) { elemLocalIdx += curNumScvs; } );
318
319 // add scvf
320 numScvf_++;
321 scvfs_[eIdx].emplace_back(geometryHelper,
322 intersection,
323 isGeometry,
324 edgeIdx,
325 scvfLocalIdx++,
326 std::move(localScvIndices),
327 intersection.boundary());
328 }
329 }
330
331 // dim == 2, intersection is an edge, make 1 scvf
332 else
333 {
334 // inside/outside scv indices in face local node numbering
335 std::vector<LocalIndexType> localScvIndices({0, 1});
336
337 // add offset such that the fracture scvs above are addressed
338 std::for_each( localScvIndices.begin(),
339 localScvIndices.end(),
340 [curNumScvs] (auto& elemLocalIdx) { elemLocalIdx += curNumScvs; } );
341
342 // add scvf
343 numScvf_++;
344 scvfs_[eIdx].emplace_back(geometryHelper,
345 intersection,
346 isGeometry,
347 /*idxOnIntersection*/0,
348 scvfLocalIdx++,
349 std::move(localScvIndices),
350 intersection.boundary());
351 }
352 }
353 }
354 }
355 }
356
358 const FeCache& feCache() const { return feCache_; }
360 const std::vector<SubControlVolume>& scvs(GridIndexType eIdx) const { return scvs_[eIdx]; }
362 const std::vector<SubControlVolumeFace>& scvfs(GridIndexType eIdx) const { return scvfs_[eIdx]; }
364 bool dofOnBoundary(unsigned int dofIdx) const { return boundaryDofIndices_[dofIdx]; }
366 bool dofOnFracture(unsigned int dofIdx) const { return fractureDofIndices_[dofIdx]; }
368 bool dofOnPeriodicBoundary(std::size_t dofIdx) const { return false; }
369
371 std::size_t periodicallyMappedDof(std::size_t dofIdx) const
372 { DUNE_THROW(Dune::InvalidStateException, "Periodic boundaries are not supported by the box-dfm scheme"); }
373
375 std::unordered_map<std::size_t, std::size_t> periodicVertexMap() const
376 { return std::unordered_map<std::size_t, std::size_t>(); }
377
378private:
379 const FeCache feCache_;
380
381 std::vector<std::vector<SubControlVolume>> scvs_;
382 std::vector<std::vector<SubControlVolumeFace>> scvfs_;
383
384 // TODO do we need those?
385 std::size_t numScv_;
386 std::size_t numScvf_;
387 std::size_t numBoundaryScvf_;
388
389 // vertices on the boundary & fracture facets
390 std::vector<bool> boundaryDofIndices_;
391 std::vector<bool> fractureDofIndices_;
392};
393
401template<class Scalar, class GV, class Traits>
402class BoxDfmFVGridGeometry<Scalar, GV, false, Traits>
403: public BaseGridGeometry<GV, Traits>
404{
407 using GridIndexType = typename GV::IndexSet::IndexType;
408
409 static const int dim = GV::dimension;
410 static const int dimWorld = GV::dimensionworld;
411
412 using Element = typename GV::template Codim<0>::Entity;
413 using Intersection = typename GV::Intersection;
414 using CoordScalar = typename GV::ctype;
415
416 using ReferenceElements = typename Dune::ReferenceElements<CoordScalar, dim>;
417 using FaceReferenceElements = typename Dune::ReferenceElements<CoordScalar, dim-1>;
418
419public:
421 static constexpr DiscretizationMethod discMethod = DiscretizationMethod::box;
422
424 using LocalView = typename Traits::template LocalView<ThisType, false>;
426 using SubControlVolume = typename Traits::SubControlVolume;
428 using SubControlVolumeFace = typename Traits::SubControlVolumeFace;
430 using DofMapper = typename Traits::VertexMapper;
432 using FeCache = Dune::PQkLocalFiniteElementCache<CoordScalar, Scalar, dim, 1>;
434 using GridView = GV;
435
438 : ParentType(gridView)
439 , facetMapper_(gridView, Dune::mcmgLayout(Dune::template Codim<1>()))
440 {}
441
444 const DofMapper& dofMapper() const
445 { return this->vertexMapper(); }
446
448 std::size_t numScv() const
449 { return numScv_; }
450
452 std::size_t numScvf() const
453 { return numScvf_; }
454
457 std::size_t numBoundaryScvf() const
458 { return numBoundaryScvf_; }
459
461 std::size_t numDofs() const
462 { return this->gridView().size(dim); }
463
465 template< class FractureGridAdapter >
466 void update(const FractureGridAdapter& fractureGridAdapter)
467 {
468 ParentType::update();
469
470 boundaryDofIndices_.assign(numDofs(), false);
471 fractureDofIndices_.assign(numDofs(), false);
472 facetOnFracture_.assign(this->gridView().size(1), false);
473
474 // save global data on the grid's scvs and scvfs
475 // TODO do we need those information?
476 numScv_ = 0;
477 numScvf_ = 0;
478 numBoundaryScvf_ = 0;
479 for (const auto& element : elements(this->gridView()))
480 {
481 numScv_ += element.subEntities(dim);
482 numScvf_ += element.subEntities(dim-1);
483
484 const auto elementGeometry = element.geometry();
485 const auto referenceElement = ReferenceElements::general(elementGeometry.type());
486
487 // store the sub control volume face indices on the domain boundary
488 for (const auto& intersection : intersections(this->gridView(), element))
489 {
490 // first, obtain all vertex indices on this intersection
491 const auto& isGeometry = intersection.geometry();
492 const auto numCorners = isGeometry.corners();
493 const auto idxInInside = intersection.indexInInside();
494
495 std::vector<GridIndexType> isVertexIndices(numCorners);
496 for (unsigned int vIdxLocal = 0; vIdxLocal < numCorners; ++vIdxLocal)
497 isVertexIndices[vIdxLocal] = this->vertexMapper().subIndex(element,
498 referenceElement.subEntity(idxInInside, 1, vIdxLocal, dim),
499 dim);
500
501 if (intersection.boundary() && !intersection.neighbor())
502 {
503 numScvf_ += numCorners;
504 numBoundaryScvf_ += numCorners;
505
506 // add all vertices on the intersection to the set of
507 // boundary vertices
508 const auto fIdx = intersection.indexInInside();
509 const auto numFaceVerts = referenceElement.size(fIdx, 1, dim);
510 for (int localVIdx = 0; localVIdx < numFaceVerts; ++localVIdx)
511 {
512 const auto vIdx = referenceElement.subEntity(fIdx, 1, localVIdx, dim);
513 const auto vIdxGlobal = this->vertexMapper().subIndex(element, vIdx, dim);
514 boundaryDofIndices_[vIdxGlobal] = true;
515 }
516 }
517
518 // make sure we have no periodic boundaries
519 else if (intersection.boundary() && intersection.neighbor())
520 DUNE_THROW(Dune::InvalidStateException, "Periodic boundaries are not supported by the box-dfm scheme");
521
522 // maybe add fracture scvs & scvfs
523 if (fractureGridAdapter.composeFacetElement(isVertexIndices))
524 {
525 facetOnFracture_[facetMapper_.subIndex(element, idxInInside, 1)] = true;
526 for (auto vIdx : isVertexIndices)
527 fractureDofIndices_[vIdx] = true;
528
529 const auto isGeometry = intersection.geometry();
530 numScv_ += isGeometry.corners();
531 numScvf_ += dim == 3 ? FaceReferenceElements::general(isGeometry.type()).size(1) : 1;
532 }
533 }
534 }
535 }
536
538 const FeCache& feCache() const { return feCache_; }
540 bool dofOnBoundary(unsigned int dofIdx) const { return boundaryDofIndices_[dofIdx]; }
542 bool dofOnFracture(unsigned int dofIdx) const { return fractureDofIndices_[dofIdx]; }
544 bool dofOnPeriodicBoundary(std::size_t dofIdx) const { return false; }
545
547 bool isOnFracture(const Element& element, const Intersection& intersection) const
548 { return facetOnFracture_[facetMapper_.subIndex(element, intersection.indexInInside(), 1)]; }
549
551 std::size_t periodicallyMappedDof(std::size_t dofIdx) const
552 { DUNE_THROW(Dune::InvalidStateException, "Periodic boundaries are not supported by the box-dfm scheme"); }
553
555 std::unordered_map<std::size_t, std::size_t> periodicVertexMap() const
556 { return std::unordered_map<std::size_t, std::size_t>(); }
557
558private:
559
560 const FeCache feCache_;
561
562 // Information on the global number of geometries
563 // TODO do we need those information?
564 std::size_t numScv_;
565 std::size_t numScvf_;
566 std::size_t numBoundaryScvf_;
567
568 // vertices on the boundary & fracture facets
569 std::vector<bool> boundaryDofIndices_;
570 std::vector<bool> fractureDofIndices_;
571
572 // facet mapper and markers which facets lie on interior boundaries
573 typename Traits::FacetMapper facetMapper_;
574 std::vector<bool> facetOnFracture_;
575};
576
577} // end namespace Dumux
578
579#endif
Defines the default element and vertex mapper types.
Helper class constructing the dual grid finite volume geometries for the box discretizazion method.
Base class for grid geometries.
The available discretization methods in Dumux.
Helper class constructing the dual grid finite volume geometries for the box discrete fracture model.
DiscretizationMethod
The available discretization methods in Dumux.
Definition: method.hh:37
Dune::IteratorRange< typename MultiDomainGlue< DomainGridView, TargetGridView, DomainMapper, TargetMapper >::Intersections::const_iterator > intersections(const MultiDomainGlue< DomainGridView, TargetGridView, DomainMapper, TargetMapper > &glue)
Range generator to iterate with range-based for loops over all intersections as follows: for (const a...
Definition: glue.hh:62
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
Definition: common/properties/model.hh:34
Base class for all finite volume grid geometries.
Definition: basegridgeometry.hh:49
GV GridView
export the grid view type
Definition: basegridgeometry.hh:64
Base class for the finite volume geometry vector for box discrete fracture model.
Definition: porousmediumflow/boxdfm/fvelementgeometry.hh:52
The default traits for the box finite volume grid geometry.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:62
Dune::MultipleCodimMultipleGeomTypeMapper< GridView > FacetMapper
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:70
Base class for the finite volume geometry vector for box schemes.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:85
Base class for the finite volume geometry vector for box schemes that consider extra connectivity bet...
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:100
typename Traits::SubControlVolume SubControlVolume
Export the type of sub control volume.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:125
std::unordered_map< std::size_t, std::size_t > periodicVertexMap() const
Returns the map between dofs across periodic boundaries.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:375
const std::vector< SubControlVolumeFace > & scvfs(GridIndexType eIdx) const
Get the local scvfs for an element.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:362
std::size_t numDofs() const
The total number of degrees of freedom.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:158
typename Traits::VertexMapper DofMapper
Export dof mapper type.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:129
BoxDfmFVGridGeometry(const GridView gridView)
Constructor.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:136
Dune::PQkLocalFiniteElementCache< CoordScalar, Scalar, dim, 1 > FeCache
Export the finite element cache type.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:131
const std::vector< SubControlVolume > & scvs(GridIndexType eIdx) const
Get the local scvs for an element.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:360
typename Traits::template LocalView< ThisType, true > LocalView
Export the type of the fv element geometry (the local view type)
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:123
const DofMapper & dofMapper() const
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:141
std::size_t numScvf() const
The total number of sun control volume faces.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:149
std::size_t numScv() const
The total number of sub control volumes.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:145
std::size_t periodicallyMappedDof(std::size_t dofIdx) const
The index of the vertex / d.o.f. on the other side of the periodic boundary.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:371
std::size_t numBoundaryScvf() const
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:154
bool dofOnBoundary(unsigned int dofIdx) const
If a vertex / d.o.f. is on the boundary.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:364
void update(const FractureGridAdapter &fractureGridAdapter)
Update all fvElementGeometries (do this again after grid adaption)
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:163
const FeCache & feCache() const
The finite element cache for creating local FE bases.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:358
typename Traits::SubControlVolumeFace SubControlVolumeFace
Export the type of sub control volume.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:127
bool dofOnPeriodicBoundary(std::size_t dofIdx) const
Periodic boundaries are not supported for the box-dfm scheme.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:368
bool dofOnFracture(unsigned int dofIdx) const
If a vertex / d.o.f. is on a fracture.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:366
Base class for the finite volume geometry vector for box schemes This builds up the sub control volum...
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:404
bool dofOnBoundary(unsigned int dofIdx) const
If a vertex / d.o.f. is on the boundary.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:540
bool dofOnPeriodicBoundary(std::size_t dofIdx) const
Periodic boundaries are not supported for the box-dfm scheme.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:544
std::size_t numScv() const
The total number of sub control volumes.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:448
bool isOnFracture(const Element &element, const Intersection &intersection) const
Returns true if an intersection coincides with a fracture element.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:547
std::unordered_map< std::size_t, std::size_t > periodicVertexMap() const
Returns the map between dofs across periodic boundaries.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:555
std::size_t numDofs() const
The total number of degrees of freedom.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:461
BoxDfmFVGridGeometry(const GridView gridView)
Constructor.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:437
std::size_t numBoundaryScvf() const
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:457
typename Traits::VertexMapper DofMapper
export dof mapper type
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:430
typename Traits::template LocalView< ThisType, false > LocalView
export the type of the fv element geometry (the local view type)
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:424
std::size_t numScvf() const
The total number of sun control volume faces.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:452
Dune::PQkLocalFiniteElementCache< CoordScalar, Scalar, dim, 1 > FeCache
export the finite element cache type
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:432
const FeCache & feCache() const
The finite element cache for creating local FE bases.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:538
const DofMapper & dofMapper() const
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:444
typename Traits::SubControlVolumeFace SubControlVolumeFace
export the type of sub control volume
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:428
void update(const FractureGridAdapter &fractureGridAdapter)
update all fvElementGeometries (do this again after grid adaption)
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:466
typename Traits::SubControlVolume SubControlVolume
export the type of sub control volume
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:426
bool dofOnFracture(unsigned int dofIdx) const
If a vertex / d.o.f. is on a fracture.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:542
std::size_t periodicallyMappedDof(std::size_t dofIdx) const
The index of the vertex / d.o.f. on the other side of the periodic boundary.
Definition: porousmediumflow/boxdfm/fvgridgeometry.hh:551
Create sub control volumes and sub control volume face geometries.
Definition: geometryhelper.hh:35
the sub control volume for the box discrete fracture scheme
Definition: porousmediumflow/boxdfm/subcontrolvolume.hh:96
Class for a sub control volume face in the box discrete fracture method, i.e a part of the boundary o...
Definition: porousmediumflow/boxdfm/subcontrolvolumeface.hh:100
Base class for the local finite volume geometry for the box discrete fracture model.
the sub control volume for the box discrete fracture scheme
The sub control volume face class for the box discrete fracture model.