dune-localfunctions  2.7.0
localfiniteelementvariantcache.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_LOCALFUNCTIONS_COMMON_LOCALFINITEELEMENTVARIANTCACHE_HH
4 #define DUNE_LOCALFUNCTIONS_COMMON_LOCALFINITEELEMENTVARIANTCACHE_HH
5 
6 #include <vector>
7 #include <tuple>
8 #include <utility>
9 #include <type_traits>
10 
11 #include <dune/common/std/type_traits.hh>
12 #include <dune/common/exceptions.hh>
13 #include <dune/common/typelist.hh>
14 #include <dune/common/hybridutilities.hh>
15 
16 #include <dune/geometry/type.hh>
17 #include <dune/geometry/typeindex.hh>
18 
20 
21 
22 namespace Dune {
23 
24 namespace Impl {
25 
26  // This class provides the index method of LocalGeometryTypeIndex
27  // but throws a Dune::RangeError if the dimension does not match.
28  // This can be helpful to catch errors in a LocalFiniteElementVariantCache
29  // instance based on dimension specific GeometryType indices.
30  template<std::size_t dim>
31  struct FixedDimLocalGeometryTypeIndex {
32  inline static std::size_t index(const GeometryType &gt)
33  {
34  if (gt.dim() != dim)
35  DUNE_THROW(Dune::RangeError, "Asking for dim=" << dim << " specific index of GeometryType with dimension " << gt.dim());
36  return LocalGeometryTypeIndex::index(gt);
37  }
38  };
39 
40 } // end namespace Impl
41 
64 template<class Base>
66 {
67 
68  template<class LFEImplTuple>
69  struct GenerateLFEVariant;
70 
71  template<class Index, class... LFEImpl>
72  struct GenerateLFEVariant<std::tuple<std::pair<Index, LFEImpl>...>>
73  {
74  using type = UniqueTypes_t<LocalFiniteElementVariant, decltype(std::declval<LFEImpl>()())...>;
75  };
76 
77  using Base::getImplementations;
78  using Base::index;
79  using Implementations = decltype(std::declval<Base>().getImplementations());
80 
81 public:
82 
90  using FiniteElementType = typename GenerateLFEVariant<Implementations>::type;
91 
96  template<class... Args>
98  Base(std::forward<Args>(args)...)
99  {
100  Dune::Hybrid::forEach(getImplementations(), [&,this](auto feImpl) {
101  auto implIndex = feImpl.first;
102  if (cache_.size() < implIndex+1)
103  cache_.resize(implIndex+1);
104  cache_[implIndex] = feImpl.second();
105  });
106  }
107 
110 
113 
118  template<class... Key>
119  const auto& get(const Key&... key) const
120  {
121  auto implIndex = index(key...);
122  if (implIndex >= cache_.size())
123  DUNE_THROW(Dune::RangeError,"There is no LocalFiniteElement of the requested type.");
124  if (not(cache_[implIndex]))
125  DUNE_THROW(Dune::RangeError,"There is no LocalFiniteElement of the requested type.");
126  return cache_[implIndex];
127  }
128 
129 private:
130  std::vector<FiniteElementType> cache_;
131 };
132 
133 
134 
135 } // namespace Dune
136 
137 
138 
139 
140 #endif // DUNE_LOCALFUNCTIONS_COMMON_LOCALFINITEELEMENTVARIANT_HH
Dune::LocalFiniteElementVariantCache
A cache storing a compile time selection of local finite element implementations.
Definition: localfiniteelementvariantcache.hh:65
Dune::LocalFiniteElementVariantCache::FiniteElementType
typename GenerateLFEVariant< Implementations >::type FiniteElementType
Type of exported LocalFiniteElement's.
Definition: localfiniteelementvariantcache.hh:90
Dune
Definition: bdfmcube.hh:15
Dune::LocalFiniteElementVariantCache::LocalFiniteElementVariantCache
LocalFiniteElementVariantCache(Args &&... args)
Default constructor.
Definition: localfiniteelementvariantcache.hh:97
Dune::LocalFiniteElementVariantCache::get
const auto & get(const Key &... key) const
Get the LocalFiniteElement for the given key data.
Definition: localfiniteelementvariantcache.hh:119
localfiniteelementvariant.hh