From 58aef1c3141f73dae2058398e898fd23d4a05646 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Dennis=20Gl=C3=A4ser?= <dennis.glaeser@iws.uni-stuttgart.de>
Date: Fri, 27 Oct 2023 13:32:38 +0200
Subject: [PATCH] [mpfa][fvelemgeom] expose vertex/facet corners
---
.../cellcentered/mpfa/fvelementgeometry.hh | 41 +++++++++++++++++++
.../mpfa/omethod/scvgeometryhelper.hh | 6 +--
.../cellcentered/mpfa/subcontrolvolumeface.hh | 3 ++
.../facet/cellcentered/mpfa/couplingmapper.hh | 2 +-
4 files changed, 48 insertions(+), 4 deletions(-)
diff --git a/dumux/discretization/cellcentered/mpfa/fvelementgeometry.hh b/dumux/discretization/cellcentered/mpfa/fvelementgeometry.hh
index 1f02570882..e1c7531bcd 100644
--- a/dumux/discretization/cellcentered/mpfa/fvelementgeometry.hh
+++ b/dumux/discretization/cellcentered/mpfa/fvelementgeometry.hh
@@ -59,6 +59,31 @@ typename SubControlVolumeFace::Traits::Geometry makeScvfGeometry(const GridGeome
DUNE_THROW(Dune::InvalidStateException, "Could not construct scvf geometry");
}
+template<typename GridGeometry, typename SubControlVolumeFace>
+auto getVertexCorner(const GridGeometry& gridGeometry, const SubControlVolumeFace& scvf)
+{
+ static constexpr int dim = GridGeometry::GridView::dimension;
+
+ const auto& facetInfo = scvf.facetInfo();
+ const auto element = gridGeometry.element(facetInfo.elementIndex);
+ const auto elemGeo = element.geometry();
+ const auto refElement = referenceElement(elemGeo);
+ return elemGeo.global(refElement.position(
+ refElement.subEntity(facetInfo.facetIndex, 1, facetInfo.facetCornerIndex, dim),
+ dim
+ ));
+}
+
+template<typename GridGeometry, typename SubControlVolumeFace>
+auto getFacetCorner(const GridGeometry& gridGeometry, const SubControlVolumeFace& scvf)
+{
+ const auto& facetInfo = scvf.facetInfo();
+ const auto element = gridGeometry.element(facetInfo.elementIndex);
+ const auto elemGeo = element.geometry();
+ const auto refElement = referenceElement(elemGeo);
+ return elemGeo.global(refElement.position(facetInfo.facetIndex, 1));
+}
+
} // namespace Detail::Mpfa
#endif // DOXYGEN
@@ -225,6 +250,14 @@ public:
typename SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace& scvf) const
{ return Detail::Mpfa::makeScvfGeometry(gridGeometry(), scvf); }
+ //! Return the position of the scvf corner that coincides with an element vertex
+ typename SubControlVolumeFace::Traits::GlobalPosition vertexCorner(const SubControlVolumeFace& scvf) const
+ { return Detail::Mpfa::getVertexCorner(gridGeometry(), scvf); }
+
+ //! Return the corner of the scvf that is inside the facet the scvf is embedded in
+ typename SubControlVolumeFace::Traits::GlobalPosition facetCorner(const SubControlVolumeFace& scvf) const
+ { return Detail::Mpfa::getFacetCorner(gridGeometry(), scvf); }
+
private:
std::optional<Element> element_;
@@ -383,6 +416,14 @@ public:
typename SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace& scvf) const
{ return Detail::Mpfa::makeScvfGeometry(gridGeometry(), scvf); }
+ //! Return the position of the scvf corner that coincides with an element vertex
+ typename SubControlVolumeFace::Traits::GlobalPosition vertexCorner(const SubControlVolumeFace& scvf) const
+ { return Detail::Mpfa::getVertexCorner(gridGeometry(), scvf); }
+
+ //! Return the corner of the scvf that is inside the facet the scvf is embedded in
+ typename SubControlVolumeFace::Traits::GlobalPosition facetCorner(const SubControlVolumeFace& scvf) const
+ { return Detail::Mpfa::getFacetCorner(gridGeometry(), scvf); }
+
private:
//! Binding of an element preparing the geometries of the whole stencil
diff --git a/dumux/discretization/cellcentered/mpfa/omethod/scvgeometryhelper.hh b/dumux/discretization/cellcentered/mpfa/omethod/scvgeometryhelper.hh
index 4d772a3982..4db457755c 100644
--- a/dumux/discretization/cellcentered/mpfa/omethod/scvgeometryhelper.hh
+++ b/dumux/discretization/cellcentered/mpfa/omethod/scvgeometryhelper.hh
@@ -71,9 +71,9 @@ public:
typename MLGTraits::template CornerStorage<dim, dimWorld>::Type corners;
corners[0] = fvGeometry.scv( scv.gridScvIndex() ).center();
- corners[1] = firstGridScvf.facetCorner();
- corners[2] = secondGridScvf.facetCorner();
- corners[3] = secondGridScvf.vertexCorner();
+ corners[1] = fvGeometry.facetCorner(firstGridScvf);
+ corners[2] = fvGeometry.facetCorner(secondGridScvf);
+ corners[3] = fvGeometry.vertexCorner(secondGridScvf);
using std::swap;
typename LocalScvType::LocalBasis basis;
diff --git a/dumux/discretization/cellcentered/mpfa/subcontrolvolumeface.hh b/dumux/discretization/cellcentered/mpfa/subcontrolvolumeface.hh
index 427617e4ae..7249250cad 100644
--- a/dumux/discretization/cellcentered/mpfa/subcontrolvolumeface.hh
+++ b/dumux/discretization/cellcentered/mpfa/subcontrolvolumeface.hh
@@ -194,14 +194,17 @@ public:
{ return outsideScvIndices_; }
//! Returns the number of corners
+ [[deprecated("Will be removed after 3.8. Use fvGeometry.geometry(scvf).corners().")]]
std::size_t corners() const
{ return corners_.size(); }
//! Returns the global position of the vertex the scvf is connected to
+ [[deprecated("Will be removed after 3.8. Use fvGeometry.vertexCorner(scvf)")]]
const GlobalPosition& vertexCorner() const
{ return corners_.back(); }
//! Returns the global position of the center of the element facet this scvf is embedded in
+ [[deprecated("Will be removed after 3.8. Use fvGeometry.facetCorner(scvf)")]]
const GlobalPosition& facetCorner() const
{ return corners_[0]; }
diff --git a/dumux/multidomain/facet/cellcentered/mpfa/couplingmapper.hh b/dumux/multidomain/facet/cellcentered/mpfa/couplingmapper.hh
index 9c9d612d1d..54fc7ca8de 100644
--- a/dumux/multidomain/facet/cellcentered/mpfa/couplingmapper.hh
+++ b/dumux/multidomain/facet/cellcentered/mpfa/couplingmapper.hh
@@ -127,7 +127,7 @@ public:
else
{
const auto eps = lowDimGeometry.volume()*1e-8;
- const auto diffVec = lowDimGeometry.center()-scvf.facetCorner();
+ const auto diffVec = lowDimGeometry.center()-fvGeometry.facetCorner(scvf);
if ( Dune::FloatCmp::eq<GlobalPosition, Dune::FloatCmp::CmpStyle::absolute>(diffVec, GlobalPosition(0.0), eps) )
embeddedScvfIndices.push_back(scvf.index());
--
GitLab