...
  View open merge request
Commits (8)
  • Carsten Gräser's avatar
    Add StaticGeometryType · 134dafdc
    Carsten Gräser authored
    This encodes all the information of a GeometryType as template
    parameters. Some remarks on the implementation:
    * This may be a solition for using GeometryType's a template paremters
      and is thus a candidate for dune-geometry
    * All methods in Dune::GeometryTypes:: have analogouges in
      Dune::Functions::StaticGeometryTypes:: and additional template
      aliases for parametrized types. Hence you can e.g. use Simplex<dim>.
    * There is LocalHybridGeometryTypeIndex which works like LocalGeometryTypeIndex
      but supports GeometryType and StaticGeometryType.
    * Some topologyIds are different: If you want to overload for specific
      StaticGeometryType's, you must be sure that the same geometry type
      leads to the same type. Hence topologyIds are normalize to have 0
      in the first bit.
    * One may consider using uniqueTopologyId = (topologyId >> 1) as
      template parameter instead to get rid of the nasty non-uniqueness.
    134dafdc
  • Carsten Gräser's avatar
    Add LagrangeFiniteElementCache · 6d25947b
    Carsten Gräser authored
    This works like PQkFiniteElementCache in dune-loicalfunctions
    but also allows to obtain raw fe-implemenations without polymorphic
    wrapper if the geometry type is provided statically as StaticGeometryType.
    This allows to seamlessly switch between the polymorphic and non-polymorphic
    version by using GeometryType or StaticGeometryType.
    
    The StaticGeometryType imlementations are stored in a TupleVector.
    6d25947b
  • Carsten Gräser's avatar
    Add different GeometryTypeProviders · 60664304
    Carsten Gräser authored
    These allow to map an entity to its geometry type. Depending on the
    selected provider and grid view, this is encoded as GeometryType or
    StaticGeometryType. There are several implementations:
    * MixedGeometryTypeProvider will always act like a mixed element
      grid and simply forward the dyanamic GeometryType. You'd normally
      not want to use this one.
    * AutoGeometryTypeProvider will use a StaticGeometryType if the grid
      satisfies hasSingleGeometryType(). Otherwise GeometryType is used.
      This is the default choise which should work always.
    * SimplexGeometryTypeProvider and CubeGeometryTypeProvider will
      always provide the corresponding StaticGeometryType. This can be used
      to avoid the costly polymorphic interface if your grid supports mixed
      elements, but you know, that there are only simplicies or cubes.
    60664304
  • Carsten Gräser's avatar
    Add support for GeometryTypeProviders to LagrangeBasis · b68cd2cb
    Carsten Gräser authored
    * For single element type grids the polymorphic interface is by
      default no longer used.
    * For mixed element type grids the polymorphic interface is by
      default used as before.
    * When explicitly switching fixing a single element type
      (in case you have additional knowledge) you can now also
      avoid the costly polymorphic interface.
    b68cd2cb
  • Carsten Gräser's avatar
  • Carsten Gräser's avatar
    Normalize topologyId of grid capability · fc40bb56
    Carsten Gräser authored
    Otherwise overloading for the StaticGeometryType
    may fail if the same geometry type is encoded
    with 0 and 1 in the first digit.
    fc40bb56
  • Carsten Gräser's avatar
    047f538c
  • Carsten Gräser's avatar
    7381f53e
......@@ -10,6 +10,7 @@ install(FILES
functionconcepts.hh
indexaccess.hh
interfaces.hh
lagrangelfecache.hh
localfunction.hh
localfunction_imp.hh
optional.hh
......@@ -17,6 +18,7 @@ install(FILES
reserveddeque.hh
signature.hh
staticforloop.hh
staticgeometrytype.hh
treedata.hh
type_traits.hh
typeerasure.hh
......
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_FUNCTIONS_COMMON_LAGRANGELFECACHE_HH
#define DUNE_FUNCTIONS_COMMON_LAGRANGELFECACHE_HH
#include <map>
#include <dune/common/tuplevector.hh>
#include <dune/common/overloadset.hh>
#include <dune/geometry/type.hh>
#include <dune/localfunctions/common/virtualinterface.hh>
#include <dune/localfunctions/common/virtualwrappers.hh>
#include <dune/localfunctions/lagrange/p0.hh>
#include <dune/localfunctions/lagrange/pk.hh>
#include <dune/localfunctions/lagrange/qk.hh>
#include <dune/localfunctions/lagrange/prismp1.hh>
#include <dune/localfunctions/lagrange/prismp2.hh>
#include <dune/localfunctions/lagrange/pyramidp1.hh>
#include <dune/localfunctions/lagrange/pyramidp2.hh>
#include <dune/localfunctions/lagrange/pqkfactory.hh>
#include <dune/functions/common/staticgeometrytype.hh>
namespace Dune {
namespace Functions {
namespace Imp {
// This class is used as fallback if there's no matching Lagrange
// finite element implementation. By having the required static
// geometry type and order as template parameter, the compiler
// will print out this information in the error message in case
// this dummy is used. This is helpful when debugging the reason.
template<class SGT, std::size_t order>
struct DummyLocalFiniteElement {};
// Create Lagrange finite element for given geometry type, ctypes and order
template<class D, class R, std::size_t order, std::size_t dim, unsigned int topologyId, bool none>
constexpr auto lagrangeFiniteElement(const StaticGeometryType<topologyId, dim, none>& type)
{
// Using orderedOverload allows to avoid explicit SFINAE
// conditions, because the first match is used which
// avoids ambiguity by construction.
return Dune::orderedOverload(
[](Dune::index_constant<0> o, auto t) {
return P0LocalFiniteElement<D,R,dim>(t);
},
[](auto o, const StaticGeometryTypes::Simplex<dim>& t) {
return PkLocalFiniteElement<D,R,dim,o>();
},
[](auto o, const StaticGeometryTypes::Cube<dim>& t) {
return QkLocalFiniteElement<D,R,dim,o>();
},
[](Dune::index_constant<1> o, StaticGeometryTypes::Prism t) {
return PrismP1LocalFiniteElement<D,R>();
},
[](Dune::index_constant<2> o, StaticGeometryTypes::Prism t) {
return PrismP2LocalFiniteElement<D,R>();
},
[](Dune::index_constant<1> o, StaticGeometryTypes::Pyramid t) {
return PyramidP1LocalFiniteElement<D,R>();
},
[](Dune::index_constant<2> o, StaticGeometryTypes::Pyramid t) {
return PyramidP2LocalFiniteElement<D,R>();
},
[](auto o, auto t) {
return Imp::DummyLocalFiniteElement<decltype(t), o>();
}
)(Dune::index_constant<order>(), type);
}
// Don't call this method. It is only used to determine the return
// value of LagrangeFiniteElementCache::get() because clang
// refuses a simple decltype(...get(...)) due to an incomplete type.
template<class D, class R, std::size_t order, std::size_t dim>
constexpr auto lagrangeFiniteElement(const GeometryType& type)
-> const typename Dune::PQkLocalFiniteElementCache<D, R, dim, order>::FiniteElementType&;
} // namespace Imp
/** \brief A cache that stores all available Pk/Qk like local finite elements for the given dimension and order
*
* In contrast to the PQkLocalFiniteElementCache in dune-localfunctions this
* class also provides access using StaticGeometryType which will directly
* provide the finite element implementation without wrapping it into the
* polymorphic interface.
*
* An interface for dealing with different vertex orders is currently missing.
*
* \tparam D Type used for domain coordinates
* \tparam R Type used for shape function values
* \tparam dim Element dimension
* \tparam order Element order
*/
template<class D, class R, std::size_t dim, std::size_t order>
class LagrangeFiniteElementCache
{
static constexpr auto createStaticCache()
{
auto gtIndices = std::make_index_sequence<LocalHybridGeometryTypeIndex::size(dim)>();
return unpackIntegerSequence([&](auto... index) {
return Dune::makeTupleVector(
Imp::lagrangeFiniteElement<D, R, order>(LocalHybridGeometryTypeIndex::type<dim, index>())...);
}, gtIndices);
}
using StaticCache = decltype(createStaticCache());
using DynamicCache = Dune::PQkLocalFiniteElementCache<D, R, dim, order>;
public:
LagrangeFiniteElementCache() :
dynamicCache_(),
staticCache_(createStaticCache())
{}
LagrangeFiniteElementCache(const LagrangeFiniteElementCache& other) = default;
template<unsigned int topologyId>
const auto& get(const StaticGeometryType<topologyId, dim, false>& gt) const
{
return staticCache_[LocalHybridGeometryTypeIndex::index(gt)];
}
const auto& get(const GeometryType& gt) const
{
return dynamicCache_.get(gt);
}
// The next typedef is a little complicated since clang does not
// like the following because it uses members of an incomplete type.
//
// template<class GT>
// using FiniteElement = std::decay_t<decltype(std::declval<LagrangeFiniteElementCache>().get(std::declval<GT>()))>;
template<class GT>
using FiniteElement = std::decay_t<decltype(Imp::lagrangeFiniteElement<D,R,order,dim>(std::declval<GT>()))>;
private:
DynamicCache dynamicCache_;
StaticCache staticCache_;
};
} // namespace Functions
} // namespace Dune
#endif // DUNE_FUNCTIONS_COMMON_LAGRANGELFECACHE_HH
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_FUNCTIONS_COMMON_STATICGEOMETRYTYPE_HH
#define DUNE_FUNCTIONS_COMMON_STATICGEOMETRYTYPE_HH
#include <dune/common/std/type_traits.hh>
#include <dune/common/indices.hh>
#include <dune/common/keywords.hh>
#include <dune/geometry/type.hh>
namespace Dune {
namespace Functions {
template<unsigned int topologyId_T, std::size_t dim_T, bool none_T>
class StaticGeometryType
{
using NormalizedType = StaticGeometryType<((topologyId_T >> 1) << 1), dim_T, none_T>;
public:
/** \brief Return the topology id of the type */
static constexpr auto id()
{
return Dune::index_constant<topologyId_T>();
}
/** \brief Return dimension of the type */
static constexpr auto dim()
{
return Dune::index_constant<dim_T>();
}
/** \brief Return true if entity is a singular of any dimension */
static constexpr auto isNone()
{
return Dune::Std::bool_constant<none_T>();
}
/** \brief Return true if entity is a vertex */
static constexpr auto isVertex()
{
return Dune::Std::bool_constant<dim()==0>();
}
/** \brief Return true if entity is a line segment */
static constexpr auto isLine()
{
return Dune::Std::bool_constant<dim()==1>();
}
/** \brief Return true if entity is a triangle */
static constexpr auto isTriangle()
{
return std::is_same<NormalizedType, StaticGeometryType<0, 2, false>>();
}
/** \brief Return true if entity is a quadrilateral */
static constexpr auto isQuadrilateral()
{
return std::is_same<NormalizedType, StaticGeometryType<0b0010, 2, false>>();
}
/** \brief Return true if entity is a tetrahedron */
static constexpr auto isTetrahedron()
{
return std::is_same<NormalizedType, StaticGeometryType<0, 3, false>>();
}
/** \brief Return true if entity is a pyramid */
static constexpr auto isPyramid()
{
return std::is_same<NormalizedType, StaticGeometryType<0b0010, 3, false>>();
}
/** \brief Return true if entity is a prism */
static constexpr auto isPrism()
{
return std::is_same<NormalizedType, StaticGeometryType<0b0100, 3, false>>();
}
/** \brief Return true if entity is a hexahedron */
static constexpr auto isHexahedron()
{
return std::is_same<NormalizedType, StaticGeometryType<0b0110, 3, false>>();
}
/** \brief Return true if entity is a simplex of any dimension */
static constexpr auto isSimplex()
{
return Dune::Std::bool_constant<(!isNone()) && ((id() | 1) == 1)>();
}
/** \brief Return true if entity is a cube of any dimension */
static constexpr auto isCube()
{
return Dune::Std::bool_constant<(!isNone()) && ((id() ^ ((1 << dim())-1)) >> 1 == 0)>();
}
static constexpr Dune::GeometryType dynamicType()
{
return Dune::GeometryType(id(), dim(), isNone());
}
constexpr operator Dune::GeometryType () const
{
return dynamicType();
}
};
namespace StaticGeometryTypes {
// Simplex *******************************************************************
template<std::size_t dim>
using Simplex = StaticGeometryType<0,dim,false>;
template<std::size_t dim>
inline constexpr auto simplex()
{
return Simplex<dim>();
}
template<std::size_t dim>
inline constexpr auto simplex(Dune::index_constant<dim>)
{
return Simplex<dim>();
}
// Cube **********************************************************************
template<std::size_t dim>
using Cube = StaticGeometryType<((dim>1) ? ((1 << dim) - 2) : 0),dim,false>;
template<std::size_t dim>
inline constexpr auto cube()
{
return Cube<dim>();
}
template<std::size_t dim>
inline constexpr auto cube(Dune::index_constant<dim>)
{
return Cube<dim>();
}
// None **********************************************************************
template<std::size_t dim>
using None = StaticGeometryType<0, dim, true>;
template<std::size_t dim>
inline constexpr auto none()
{
return None<dim>();
}
template<std::size_t dim>
inline constexpr auto none(Dune::index_constant<dim>)
{
return None<dim>();
}
#ifndef __cpp_inline_variables
namespace {
#endif
using Vertex = StaticGeometryType<0,0,false>;
DUNE_INLINE_VARIABLE constexpr auto vertex = Vertex();
using Line = StaticGeometryType<0,1,false>;
DUNE_INLINE_VARIABLE constexpr auto line = Line();
using Triangle = StaticGeometryTypes::Simplex<2>;
DUNE_INLINE_VARIABLE constexpr auto triangle = Triangle();
using Quadrilateral = StaticGeometryTypes::Cube<2>;
DUNE_INLINE_VARIABLE constexpr auto quadrilateral = Quadrilateral();
using Tetrahedron = StaticGeometryTypes::Simplex<3>;
DUNE_INLINE_VARIABLE constexpr auto tetrahedron = Tetrahedron();
using Pyramid = StaticGeometryType<0b0010,3,false>;
DUNE_INLINE_VARIABLE constexpr auto pyramid = Pyramid();
using Prism = StaticGeometryType<0b0100,3,false>;
DUNE_INLINE_VARIABLE constexpr auto prism = Prism();
using Hexahedron = StaticGeometryTypes::Cube<3>;
DUNE_INLINE_VARIABLE constexpr auto hexahedron = Hexahedron();
#ifndef __cpp_inline_variables
}
#endif
} // namespace StaticGeometryTypes
//! Compute per-dimension indices for geometry types
class LocalHybridGeometryTypeIndex
{
inline static constexpr std::size_t regular_size(std::size_t dim)
{
return (1 << dim) - ((1 << dim) >> 1);
}
public:
/**
* \brief Compute total number of geometry types for the given dimension
*
* This includes irregular geometry types such as "None".
*/
inline static constexpr std::size_t size(std::size_t dim)
{
return regular_size(dim) + 1;
}
/**
* \brief Compute the index for the given GeometryType
*/
inline static constexpr std::size_t index(const GeometryType &gt)
{
return gt.isNone() ? regular_size(gt.dim()) : (gt.id() >> 1);
}
/**
* \brief Compute the index for the given GeometryType
*/
template<unsigned int topologyId, std::size_t dim, bool none>
inline static constexpr auto index(const StaticGeometryType<topologyId, dim, none>& gt)
{
return Dune::index_constant<none ? regular_size(dim) : (topologyId >> 1)>();
}
//! compute the geometry type for the given local index and dimension
inline static constexpr GeometryType type(std::size_t dim, std::size_t index) {
return (index == regular_size(dim)) ?
GeometryTypes::none(dim) :
GeometryType(static_cast< unsigned int >(index << 1), dim);
}
//! compute the geometry type for the given local index and dimension
template<std::size_t dim, std::size_t index>
inline static constexpr auto type()
{
return std::conditional_t<(index == regular_size(dim)),
StaticGeometryTypes::None<dim>,
StaticGeometryType<(index << 1), dim, false>>();
}
template<std::size_t dim, std::size_t index>
inline static constexpr auto type(Dune::index_constant<dim>, Dune::index_constant<index>)
{
return type<dim, index>();
}
};
} // namespace Functions
} // namespace Dune
#endif // DUNE_FUNCTIONS_COMMON_STATICGEOMETRYTYPE_HH
......@@ -13,6 +13,7 @@ install(FILES
defaultnodetorangemap.hh
flatmultiindex.hh
flatvectorview.hh
geometrytypeprovider.hh
hierarchicnodetorangemap.hh
hierarchicvectorwrapper.hh
interpolate.hh
......
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_GEOMETRYTYPEPROVIDER_HH
#define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_GEOMETRYTYPEPROVIDER_HH
#include <dune/common/std/type_traits.hh>
#include <dune/common/indices.hh>
#include <dune/common/keywords.hh>
#include <dune/geometry/type.hh>
#include <dune/functions/common/staticgeometrytype.hh>
namespace Dune {
namespace Functions {
struct MixedGeometryTypeProvider
{
template<class GridView, class Entity>
static constexpr auto type(const Entity& entity)
{
return entity.type();
}
template<class GridView>
using Type = GeometryType;
};
struct SimplexGeometryTypeProvider
{
template<class GridView, class Entity>
static constexpr auto type(const Entity& entity)
{
return StaticGeometryTypes::Simplex<GridView::dimension>();
}
template<class GridView>
using Type = StaticGeometryTypes::Simplex<GridView::dimension>;
};
struct CubeGeometryTypeProvider
{
template<class GridView, class Entity>
static constexpr auto type(const Entity& entity)
{
return StaticGeometryTypes::Cube<GridView::dimension>();
}
template<class GridView>
using Type = StaticGeometryTypes::Cube<GridView::dimension>;
};
struct AutoGeometryTypeProvider
{
template<class GridView, class Entity,
std::enable_if_t<Dune::Capabilities::hasSingleGeometryType<typename GridView::Grid>::v, int> = 0>
static constexpr auto type(const Entity& entity)
{
constexpr unsigned int topologyId = Dune::Capabilities::hasSingleGeometryType<typename GridView::Grid>::topologyId;
constexpr unsigned int normalizedTopologyId = (topologyId >> 1) << 1;
return StaticGeometryType<normalizedTopologyId, GridView::dimension, false>{};
}
template<class GridView, class Entity,
std::enable_if_t<not Dune::Capabilities::hasSingleGeometryType<typename GridView::Grid>::v, int> = 0>
static constexpr auto type(const Entity& entity)
{
return entity.type();
}
template<class GridView>
using Type = decltype(type<GridView>(std::declval<typename GridView::template Codim<0>::Entity>()));
};
} // namespace Functions
} // namespace Dune
#endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_GEOMETRYTYPEPROVIDER_HH