dune-localfunctions  2.9.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 // SPDX-FileCopyrightInfo: Copyright (C) DUNE Project contributors, see file LICENSE.md in module root
4 // SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
5 #ifndef DUNE_LOCALFUNCTIONS_COMMON_LOCALFINITEELEMENTVARIANTCACHE_HH
6 #define DUNE_LOCALFUNCTIONS_COMMON_LOCALFINITEELEMENTVARIANTCACHE_HH
7 
8 #include <vector>
9 #include <tuple>
10 #include <utility>
11 #include <type_traits>
12 
13 #include <dune/common/std/type_traits.hh>
14 #include <dune/common/exceptions.hh>
15 #include <dune/common/typelist.hh>
16 #include <dune/common/hybridutilities.hh>
17 
18 #include <dune/geometry/type.hh>
19 #include <dune/geometry/typeindex.hh>
20 
22 
23 
24 namespace Dune {
25 
26 namespace Impl {
27 
28  // This class provides the index method of LocalGeometryTypeIndex
29  // but throws a Dune::RangeError if the dimension does not match.
30  // This can be helpful to catch errors in a LocalFiniteElementVariantCache
31  // instance based on dimension specific GeometryType indices.
32  template<std::size_t dim>
33  struct FixedDimLocalGeometryTypeIndex {
34  inline static std::size_t index(const GeometryType &gt)
35  {
36  if (gt.dim() != dim)
37  DUNE_THROW(Dune::RangeError, "Asking for dim=" << dim << " specific index of GeometryType with dimension " << gt.dim());
38  return LocalGeometryTypeIndex::index(gt);
39  }
40  };
41 
42 } // end namespace Impl
43 
66 template<class Base>
68 {
69 
70  template<class LFEImplTuple>
71  struct GenerateLFEVariant;
72 
73  template<class Index, class... LFEImpl>
74  struct GenerateLFEVariant<std::tuple<std::pair<Index, LFEImpl>...>>
75  {
76  using type = UniqueTypes_t<LocalFiniteElementVariant, decltype(std::declval<LFEImpl>()())...>;
77  };
78 
79  using Base::getImplementations;
80  using Base::index;
81  using Implementations = decltype(std::declval<Base>().getImplementations());
82 
83 public:
84 
92  using FiniteElementType = typename GenerateLFEVariant<Implementations>::type;
93 
98  template<class... Args>
100  Base(std::forward<Args>(args)...)
101  {
102  Dune::Hybrid::forEach(getImplementations(), [&,this](auto feImpl) {
103  auto implIndex = feImpl.first;
104  if (cache_.size() < implIndex+1)
105  cache_.resize(implIndex+1);
106  cache_[implIndex] = feImpl.second();
107  });
108  }
109 
112 
115 
120  template<class... Key>
121  const auto& get(const Key&... key) const
122  {
123  auto implIndex = index(key...);
124  if (implIndex >= cache_.size())
125  DUNE_THROW(Dune::RangeError,"There is no LocalFiniteElement of the requested type.");
126  if (not(cache_[implIndex]))
127  DUNE_THROW(Dune::RangeError,"There is no LocalFiniteElement of the requested type.");
128  return cache_[implIndex];
129  }
130 
131 private:
132  std::vector<FiniteElementType> cache_;
133 };
134 
135 
136 
137 } // namespace Dune
138 
139 
140 
141 
142 #endif // DUNE_LOCALFUNCTIONS_COMMON_LOCALFINITEELEMENTVARIANT_HH
Definition: bdfmcube.hh:18
A cache storing a compile time selection of local finite element implementations.
Definition: localfiniteelementvariantcache.hh:68
LocalFiniteElementVariantCache(Args &&... args)
Default constructor.
Definition: localfiniteelementvariantcache.hh:99
LocalFiniteElementVariantCache(LocalFiniteElementVariantCache &&other)=default
Move constructor.
LocalFiniteElementVariantCache(const LocalFiniteElementVariantCache &other)=default
Copy constructor.
const auto & get(const Key &... key) const
Get the LocalFiniteElement for the given key data.
Definition: localfiniteelementvariantcache.hh:121
typename GenerateLFEVariant< Implementations >::type FiniteElementType
Type of exported LocalFiniteElement's.
Definition: localfiniteelementvariantcache.hh:92