Skip to content
Snippets Groups Projects

Feature/fix deprecation warnings

3 files
+ 65
3
Compare changes
  • Side-by-side
  • Inline
Files
3
@@ -5,7 +5,13 @@
* \brief The SubGridGeometry class and its specializations
*/
#include <type_traits>
#include <dune/common/version.hh>
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
#include <dune/common/std/variant.hh>
#else
#include <variant>
#endif
#include <dune/common/fmatrix.hh>
#include <dune/common/typetraits.hh>
@@ -147,7 +153,11 @@ class SubGridLocalGeometry
using MultiLinGeometry = MultiLinearGeometry<ctype, mydim, coorddim>;
// use a variant to store either a HostGridLocalGeometry or a MultiLinearGeometry
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
using GeometryContainer = Std::variant<HostGridLocalGeometry, MultiLinGeometry>;
#else
using GeometryContainer = std::variant<HostGridLocalGeometry, MultiLinGeometry>;
#endif
public:
@@ -177,45 +187,77 @@ class SubGridLocalGeometry
/** \brief Return the element type identifier
*/
GeometryType type () const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.type();}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.type();}, localGeometry_);
#endif
}
/** \brief Return true if the geometry mapping is affine and false otherwise */
bool affine() const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.affine();}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.affine();}, localGeometry_);
#endif
}
//! return the number of corners of this element. Corners are numbered 0...n-1
int corners () const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.corners();}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.corners();}, localGeometry_);
#endif
}
//! access to coordinates of corners. Index is the number of the corner
FieldVector<ctype, coorddim> corner (int i) const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.corner(i);}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.corner(i);}, localGeometry_);
#endif
}
/** \brief Maps a local coordinate within reference element to
* global coordinate in element */
FieldVector<ctype, coorddim> global (const FieldVector<ctype, mydim>& local) const{
FieldVector<ctype, coorddim> global (const FieldVector<ctype, mydim>& local) const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.global(local);}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.global(local);}, localGeometry_);
#endif
}
/** \brief Maps a global coordinate within the element to a
* local coordinate in its reference element */
FieldVector<ctype, mydim> local (const FieldVector<ctype, coorddim>& global) const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.local(global);}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.local(global);}, localGeometry_);
#endif
}
/**
*/
ctype integrationElement (const FieldVector<ctype, mydim>& local) const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.integrationElement(local);}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.integrationElement(local);}, localGeometry_);
#endif
}
/** \brief return volume of geometry */
ctype volume () const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.volume();}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.volume();}, localGeometry_);
#endif
}
/** \brief return center of geometry
@@ -228,17 +270,29 @@ class SubGridLocalGeometry
* find reasonably efficient ways to implement it properly.
*/
FieldVector<ctype, coorddim> center () const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return geom.center();}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return geom.center();}, localGeometry_);
#endif
}
//! The Jacobian matrix of the mapping from the reference element to this element
const JacobianTransposed jacobianTransposed (const FieldVector<ctype, mydim>& local) const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return static_cast<JacobianTransposed>(geom.jacobianTransposed(local));}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return static_cast<JacobianTransposed>(geom.jacobianTransposed(local));}, localGeometry_);
#endif
}
//! The inverse of the Jacobian matrix of the mapping from the reference element to this element
const JacobianInverseTransposed jacobianInverseTransposed (const FieldVector<ctype, mydim>& local) const {
#if DUNE_VERSION_LT(DUNE_COMMON,2,8)
return Std::visit([&](auto&& geom) {return static_cast<JacobianInverseTransposed>(geom.jacobianInverseTransposed(local));}, localGeometry_);
#else
return std::visit([&](auto&& geom) {return static_cast<JacobianInverseTransposed>(geom.jacobianInverseTransposed(local));}, localGeometry_);
#endif
}
private:
Loading