15#ifndef DUMUX_FACETCOUPLING_BOX_GRID_FVGEOMETRY_HH
16#define DUMUX_FACETCOUPLING_BOX_GRID_FVGEOMETRY_HH
21#include <dune/grid/common/mcmgmapper.hh>
22#include <dune/localfunctions/lagrange/lagrangelfecache.hh>
38template<
class GV,
class T>
53template<
class Gr
idView>
60 template<
class Gr
idGeometry,
bool enableCache>
64 using ElementMapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
68 using FacetMapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
80 bool enableGridGeometryCache =
false,
91template<
class Scalar,
class GV,
class Traits>
100 using Element =
typename GV::template Codim<0>::Entity;
101 using CoordScalar =
typename GV::ctype;
102 static const int dim = GV::dimension;
103 static const int dimWorld = GV::dimensionworld;
108 static constexpr DiscretizationMethod discMethod{};
111 using LocalView =
typename Traits::template LocalView<ThisType, true>;
113 using SubControlVolume =
typename Traits::SubControlVolume;
115 using SubControlVolumeFace =
typename Traits::SubControlVolumeFace;
119 using DofMapper =
typename Traits::VertexMapper;
121 using FeCache = Dune::LagrangeLocalFiniteElementCache<CoordScalar, Scalar, dim, 1>;
128 template<
class FacetGr
idView,
class CodimOneGr
idAdapter>
130 const FacetGridView& facetGridView,
132 bool verbose =
false)
133 : ParentType(gridView)
135 update_(facetGridView, codimOneGridAdapter, verbose);
139 const DofMapper& dofMapper()
const
140 {
return this->vertexMapper(); }
143 std::size_t numScv()
const
147 std::size_t numScvf()
const
152 std::size_t numBoundaryScvf()
const
153 {
return numBoundaryScvf_; }
156 std::size_t numDofs()
const
157 {
return this->vertexMapper().size(); }
172 template<
class FacetGr
idView,
class CodimOneGr
idAdapter>
173 void update(
const GridView& gridView,
174 const FacetGridView& facetGridView,
175 const CodimOneGridAdapter& codimOneGridAdapter,
176 bool verbose =
false)
178 ParentType::update(gridView);
179 update_(facetGridView, codimOneGridAdapter, verbose);
183 template<
class FacetGr
idView,
class CodimOneGr
idAdapter>
184 void update(GridView&& gridView,
185 const FacetGridView& facetGridView,
186 const CodimOneGridAdapter& codimOneGridAdapter,
187 bool verbose =
false)
189 ParentType::update(std::move(gridView));
190 update_(facetGridView, codimOneGridAdapter, verbose);
194 const FeCache& feCache()
const
198 const std::vector<SubControlVolume>& scvs(GridIndexType eIdx)
const
199 {
return scvs_[eIdx]; }
202 const std::vector<SubControlVolumeFace>& scvfs(GridIndexType eIdx)
const
203 {
return scvfs_[eIdx]; }
206 bool dofOnBoundary(GridIndexType dofIdx)
const
207 {
return boundaryDofIndices_[dofIdx]; }
210 bool dofOnInteriorBoundary(GridIndexType dofIdx)
const
211 {
return interiorBoundaryDofIndices_[dofIdx]; }
214 bool dofOnPeriodicBoundary(GridIndexType dofIdx)
const
218 GridIndexType periodicallyMappedDof(GridIndexType dofIdx)
const
219 { DUNE_THROW(Dune::InvalidStateException,
"Periodic boundaries are not supported by the box facet coupling scheme"); }
223 template<
class FacetGr
idView,
class CodimOneGr
idAdapter>
224 void update_(
const FacetGridView& facetGridView,
225 const CodimOneGridAdapter& codimOneGridAdapter,
226 bool verbose =
false)
229 this->vertexMapper().enrich(facetGridView, codimOneGridAdapter, verbose);
232 const auto numDof = numDofs();
233 const auto numElements = this->gridView().size(0);
236 scvs_.resize(numElements);
237 scvfs_.resize(numElements);
238 boundaryDofIndices_.assign(numDof,
false);
239 interiorBoundaryDofIndices_.assign(numDof,
false);
244 numBoundaryScvf_ = 0;
245 for (
const auto& element : elements(this->gridView()))
247 auto eIdx = this->elementMapper().index(element);
250 numScv_ +=
element.subEntities(dim);
251 numScvf_ +=
element.subEntities(dim-1);
254 auto elementGeometry =
element.geometry();
255 const auto refElement = referenceElement(elementGeometry);
258 GeometryHelper geometryHelper(elementGeometry);
262 scvs_[eIdx].reserve(elementGeometry.corners());
263 for (LocalIndexType scvLocalIdx = 0; scvLocalIdx < elementGeometry.corners(); ++scvLocalIdx)
264 scvs_[eIdx].emplace_back(geometryHelper.getScvCorners(scvLocalIdx),
267 this->vertexMapper().subIndex(element, scvLocalIdx, dim));
270 LocalIndexType scvfLocalIdx = 0;
271 scvfs_[eIdx].clear();
272 scvfs_[eIdx].reserve(
element.subEntities(dim-1));
273 for (; scvfLocalIdx <
element.subEntities(dim-1); ++scvfLocalIdx)
276 std::vector<LocalIndexType> localScvIndices({
static_cast<LocalIndexType
>(refElement.subEntity(scvfLocalIdx, dim-1, 0, dim)),
277 static_cast<LocalIndexType
>(refElement.subEntity(scvfLocalIdx, dim-1, 1, dim))});
280 scvfs_[eIdx].emplace_back(geometryHelper,
284 std::move(localScvIndices));
289 std::vector<unsigned int> handledFacets;
290 for (
const auto& intersection : intersections(this->gridView(), element))
292 if (std::count(handledFacets.begin(), handledFacets.end(), intersection.indexInInside()))
295 handledFacets.push_back(intersection.indexInInside());
298 const auto isGeometry = intersection.geometry();
299 const auto numFaceCorners = isGeometry.corners();
300 const auto idxInInside = intersection.indexInInside();
301 const auto boundary = intersection.boundary();
303 std::vector<LocalIndexType> vIndicesLocal(numFaceCorners);
304 for (
int i = 0; i < numFaceCorners; ++i)
305 vIndicesLocal[i] =
static_cast<LocalIndexType
>(refElement.subEntity(idxInInside, 1, i, dim));
307 std::vector<LocalIndexType> gridVertexIndices(numFaceCorners);
308 for (
int i = 0; i < numFaceCorners; ++i)
309 gridVertexIndices[i] = this->vertexMapper().vertexIndex(element, vIndicesLocal[i], dim);
312 const bool isOnFacet = codimOneGridAdapter.composeFacetElement(gridVertexIndices);
315 if (boundary && intersection.neighbor())
316 DUNE_THROW(Dune::InvalidStateException,
"Periodic boundaries are not supported by the box facet coupling scheme");
318 if (isOnFacet || boundary)
321 numScvf_ += numFaceCorners;
322 numBoundaryScvf_ += int(boundary)*numFaceCorners;
324 for (
unsigned int isScvfLocalIdx = 0; isScvfLocalIdx < numFaceCorners; ++isScvfLocalIdx)
327 std::vector<LocalIndexType> localScvIndices = {vIndicesLocal[isScvfLocalIdx], vIndicesLocal[isScvfLocalIdx]};
330 scvfs_[eIdx].emplace_back(geometryHelper,
335 std::move(localScvIndices),
340 const auto dofIndex = this->vertexMapper().subIndex(element, vIndicesLocal[isScvfLocalIdx], dim);
341 if (boundary) boundaryDofIndices_[ dofIndex ] = boundary && !isOnFacet;
342 if (isOnFacet) interiorBoundaryDofIndices_[ dofIndex ] = isOnFacet;
352 const FeCache feCache_;
354 std::vector<std::vector<SubControlVolume>> scvs_;
355 std::vector<std::vector<SubControlVolumeFace>> scvfs_;
359 std::size_t numScvf_;
360 std::size_t numBoundaryScvf_;
363 std::vector<bool> boundaryDofIndices_;
364 std::vector<bool> interiorBoundaryDofIndices_;
374template<
class Scalar,
class GV,
class Traits>
378 using ThisType = BoxFacetCouplingFVGridGeometry<Scalar, GV, false, Traits>;
379 using ParentType = BaseGridGeometry<GV, Traits>;
383 static const int dim = GV::dimension;
384 static const int dimWorld = GV::dimensionworld;
386 using Element =
typename GV::template Codim<0>::Entity;
387 using Intersection =
typename GV::Intersection;
388 using CoordScalar =
typename GV::ctype;
393 static constexpr DiscretizationMethod discMethod{};
396 using LocalView =
typename Traits::template LocalView<ThisType, false>;
398 using SubControlVolume =
typename Traits::SubControlVolume;
400 using SubControlVolumeFace =
typename Traits::SubControlVolumeFace;
404 using DofMapper =
typename Traits::VertexMapper;
406 using FeCache = Dune::LagrangeLocalFiniteElementCache<CoordScalar, Scalar, dim, 1>;
413 template<
class FacetGr
idView,
class CodimOneGr
idAdapter>
414 BoxFacetCouplingFVGridGeometry(
const GridView& gridView,
415 const FacetGridView& facetGridView,
416 const CodimOneGridAdapter& codimOneGridAdapter,
417 bool verbose =
false)
418 : ParentType(gridView)
419 , facetMapper_(gridView, Dune::mcmgLayout(Dune::template Codim<1>()))
421 update_(facetGridView, codimOneGridAdapter, verbose);
426 const DofMapper& dofMapper()
const
427 {
return this->vertexMapper(); }
430 std::size_t numScv()
const
434 std::size_t numScvf()
const
439 std::size_t numBoundaryScvf()
const
440 {
return numBoundaryScvf_; }
443 std::size_t numDofs()
const
444 {
return this->vertexMapper().size(); }
459 template<
class FacetGr
idView,
class CodimOneGr
idAdapter>
460 void update(
const GridView& gridView,
461 const FacetGridView& facetGridView,
462 const CodimOneGridAdapter& codimOneGridAdapter,
463 bool verbose =
false)
465 ParentType::update(gridView);
466 updateFacetMapper_();
467 update_(facetGridView, codimOneGridAdapter, verbose);
471 template<
class FacetGr
idView,
class CodimOneGr
idAdapter>
472 void update(GridView&& gridView,
473 const FacetGridView& facetGridView,
474 const CodimOneGridAdapter& codimOneGridAdapter,
475 bool verbose =
false)
477 ParentType::update(std::move(gridView));
478 updateFacetMapper_();
479 update_(facetGridView, codimOneGridAdapter, verbose);
483 const FeCache& feCache()
const
487 bool dofOnBoundary(
unsigned int dofIdx)
const
488 {
return boundaryDofIndices_[dofIdx]; }
491 bool dofOnInteriorBoundary(
unsigned int dofIdx)
const
492 {
return interiorBoundaryDofIndices_[dofIdx]; }
495 bool isOnInteriorBoundary(
const Element& element,
const Intersection& intersection)
const
496 {
return facetIsOnInteriorBoundary_[ facetMapper_.subIndex(element, intersection.indexInInside(), 1) ]; }
499 bool dofOnPeriodicBoundary(GridIndexType dofIdx)
const
503 GridIndexType periodicallyMappedDof(GridIndexType dofIdx)
const
504 { DUNE_THROW(Dune::InvalidStateException,
"Periodic boundaries are not supported by the facet coupling scheme"); }
508 void updateFacetMapper_()
510 facetMapper_.update(this->gridView());
513 template<
class FacetGr
idView,
class CodimOneGr
idAdapter>
514 void update_(
const FacetGridView& facetGridView,
515 const CodimOneGridAdapter& codimOneGridAdapter,
519 this->vertexMapper().enrich(facetGridView, codimOneGridAdapter, verbose);
525 numBoundaryScvf_ = 0;
527 const auto numDof = numDofs();
528 boundaryDofIndices_.assign(numDof,
false);
529 interiorBoundaryDofIndices_.assign(numDof,
false);
530 facetIsOnInteriorBoundary_.assign(this->gridView().size(1),
false);
531 for (
const auto& element : elements(this->gridView()))
533 numScv_ +=
element.subEntities(dim);
534 numScvf_ +=
element.subEntities(dim-1);
536 const auto elementGeometry =
element.geometry();
537 const auto refElement = referenceElement(elementGeometry);
541 std::vector<unsigned int> handledFacets;
542 for (
const auto& intersection : intersections(this->gridView(), element))
544 if (std::count(handledFacets.begin(), handledFacets.end(), intersection.indexInInside()))
547 handledFacets.push_back(intersection.indexInInside());
550 const auto isGeometry = intersection.geometry();
551 const auto numFaceCorners = isGeometry.corners();
552 const auto idxInInside = intersection.indexInInside();
553 const auto boundary = intersection.boundary();
555 std::vector<LocalIndexType> vIndicesLocal(numFaceCorners);
556 for (
int i = 0; i < numFaceCorners; ++i)
557 vIndicesLocal[i] =
static_cast<LocalIndexType
>(refElement.subEntity(idxInInside, 1, i, dim));
559 std::vector<GridIndexType> gridVertexIndices(numFaceCorners);
560 for (
int i = 0; i < numFaceCorners; ++i)
561 gridVertexIndices[i] = this->vertexMapper().vertexIndex(element, vIndicesLocal[i], dim);
564 const bool isOnFacet = codimOneGridAdapter.composeFacetElement(gridVertexIndices);
567 if (boundary && intersection.neighbor())
568 DUNE_THROW(Dune::InvalidStateException,
"Periodic boundaries are not supported by the box facet coupling scheme");
570 if (isOnFacet || boundary)
572 numScvf_ += numFaceCorners;
573 numBoundaryScvf_ += int(boundary)*numFaceCorners;
576 for (
int i = 0; i < numFaceCorners; ++i)
578 const auto dofIndex = this->vertexMapper().subIndex(element, vIndicesLocal[i], dim);
579 if (boundary) boundaryDofIndices_[ dofIndex ] = boundary && !isOnFacet;
582 interiorBoundaryDofIndices_[ dofIndex ] =
true;
583 facetIsOnInteriorBoundary_[ facetMapper_.subIndex(element, idxInInside, 1) ] =
true;
591 const FeCache feCache_;
596 std::size_t numScvf_;
597 std::size_t numBoundaryScvf_;
600 std::vector<bool> boundaryDofIndices_;
601 std::vector<bool> interiorBoundaryDofIndices_;
604 typename Traits::FacetMapper facetMapper_;
605 std::vector<bool> facetIsOnInteriorBoundary_;
Base class for grid geometries.
Helper class constructing the dual grid finite volume geometries for the box discretizazion method.
Base class for all grid geometries.
Definition basegridgeometry.hh:52
Base class for the element-local finite volume geometry for box models in the context of models consi...
Definition multidomain/facet/box/fvelementgeometry.hh:36
Base class for the finite volume geometry vector for box schemes in the context of coupled models whe...
Definition multidomain/facet/box/fvgridgeometry.hh:82
Class for a sub control volume face in the box method, i.e a part of the boundary of a sub control vo...
Definition multidomain/facet/box/subcontrolvolumeface.hh:41
Create sub control volumes and sub control volume face geometries.
Definition boxgeometryhelper.hh:257
the sub control volume for the box scheme
Definition discretization/box/subcontrolvolume.hh:60
Adapter that allows retrieving information on a d-dimensional grid for entities of a (d-1)-dimensiona...
Definition codimonegridadapter.hh:40
A vertex mapper that allows for enrichment of nodes. Indication on where to enrich the nodes is done ...
Definition vertexmapper.hh:126
Traits extracting the public Extrusion type from T Defaults to NoExtrusion if no such type is found.
Definition extrusion.hh:155
the sub control volume for the box scheme
Helper classes to compute the integration elements.
The available discretization methods in Dumux.
Base class for the element-local finite volume geometry for box models in the context of models consi...
Base class for a sub control volume face of the box method in the context of of models considering co...
Distance implementation details.
Definition cvfelocalresidual.hh:25
Dune::Std::detected_or_t< Dumux::BoxGeometryHelper< GV, GV::dimension, typename T::SubControlVolume, typename T::SubControlVolumeFace >, SpecifiesGeometryHelper, T > BoxFacetCouplingGeometryHelper_t
Definition multidomain/facet/box/fvgridgeometry.hh:39
typename T::GeometryHelper SpecifiesGeometryHelper
Definition basegridgeometry.hh:30
CVFE< CVFEMethods::PQ1 > Box
Definition method.hh:94
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition extrusion.hh:166
The default traits for the finite volume grid geometry of the box scheme with coupling occurring acro...
Definition multidomain/facet/box/fvgridgeometry.hh:55
BoxFacetCouplingFVElementGeometry< GridGeometry, enableCache > LocalView
Definition multidomain/facet/box/fvgridgeometry.hh:61
Dune::MultipleCodimMultipleGeomTypeMapper< GridView > FacetMapper
Definition multidomain/facet/box/fvgridgeometry.hh:68
BoxFacetCouplingSubControlVolumeFace< GridView > SubControlVolumeFace
Definition multidomain/facet/box/fvgridgeometry.hh:58
EnrichedVertexDofMapper< GridView > VertexMapper
Definition multidomain/facet/box/fvgridgeometry.hh:66
Dune::MultipleCodimMultipleGeomTypeMapper< GridView > ElementMapper
Definition multidomain/facet/box/fvgridgeometry.hh:64
BoxSubControlVolume< GridView > SubControlVolume
Definition multidomain/facet/box/fvgridgeometry.hh:57
typename GridView::IndexSet::IndexType GridIndex
Definition indextraits.hh:27
unsigned int LocalIndex
Definition indextraits.hh:28
A vertex mapper that allows for enrichment of nodes. Indication on where to enrich the nodes is done ...