3 #ifndef DUNE_LOCALFUNCTIONS_COMMON_LOCALFINITEELEMENTVARIANT_HH
4 #define DUNE_LOCALFUNCTIONS_COMMON_LOCALFINITEELEMENTVARIANT_HH
10 #include <dune/common/typeutilities.hh>
11 #include <dune/common/std/type_traits.hh>
12 #include <dune/common/overloadset.hh>
14 #include <dune/geometry/type.hh>
30 template<
class Visitor,
class Variant>
31 void visitIf(Visitor&& visitor, Variant&& variant)
33 auto visitorWithFallback = overload([&](std::monostate& impl) {}, [&](
const std::monostate& impl) {}, visitor);
34 std::visit(visitorWithFallback, variant);
37 template<
class... Implementations>
38 class LocalBasisVariant
41 template<
class I0,
class... II>
45 using FirstImpTraits =
typename FirstType<Implementations...>::type::Traits;
53 typename FirstImpTraits::DomainFieldType,
54 FirstImpTraits::dimDomain,
55 typename FirstImpTraits::DomainType,
56 typename FirstImpTraits::RangeFieldType,
57 FirstImpTraits::dimRange,
58 typename FirstImpTraits::RangeType,
59 typename FirstImpTraits::JacobianType>;
61 template<
class Implementation>
62 LocalBasisVariant(
const Implementation& impl) :
68 LocalBasisVariant() =
default;
69 LocalBasisVariant(
const LocalBasisVariant& other) =
default;
70 LocalBasisVariant(LocalBasisVariant&& other) =
default;
71 LocalBasisVariant& operator=(
const LocalBasisVariant& other) =
default;
72 LocalBasisVariant& operator=(LocalBasisVariant&& other) =
default;
77 unsigned int size()
const
85 unsigned int order()
const
93 inline void evaluateFunction(
94 const typename Traits::DomainType& x,
95 std::vector<typename Traits::RangeType>& out)
const
97 Impl::visitIf([&](
const auto* impl) { impl->evaluateFunction(x, out); }, impl_);
103 inline void evaluateJacobian(
104 const typename Traits::DomainType& x,
105 std::vector<typename Traits::JacobianType>& out)
const
107 Impl::visitIf([&](
const auto* impl) { impl->evaluateJacobian(x, out); }, impl_);
118 const std::array<unsigned int,Traits::dimDomain>& order,
119 const typename Traits::DomainType& x,
120 std::vector<typename Traits::RangeType>& out)
const
122 Impl::visitIf([&](
const auto* impl) { impl->partial(order, x, out); }, impl_);
126 std::variant<std::monostate,
const Implementations*...> impl_;
132 template<
class... Implementations>
133 class LocalCoefficientsVariant
137 template<
class Implementation>
138 LocalCoefficientsVariant(
const Implementation& impl) :
143 LocalCoefficientsVariant() =
default;
144 LocalCoefficientsVariant(
const LocalCoefficientsVariant& other) =
default;
145 LocalCoefficientsVariant(LocalCoefficientsVariant&& other) =
default;
146 LocalCoefficientsVariant& operator=(
const LocalCoefficientsVariant& other) =
default;
147 LocalCoefficientsVariant& operator=(LocalCoefficientsVariant&& other) =
default;
152 unsigned int size()
const
164 return std::visit(overload(
165 [&](
const std::monostate& impl) -> decltype(
auto) {
return (dummyLocalKey);},
166 [&](
const auto* impl) -> decltype(
auto) {
return impl->localKey(i); }), impl_);
170 std::variant<std::monostate,
const Implementations*...> impl_;
175 template<
class... Implementations>
176 class LocalInterpolationVariant
180 template<
class Implementation>
181 LocalInterpolationVariant(
const Implementation& impl) :
185 LocalInterpolationVariant() =
default;
186 LocalInterpolationVariant(
const LocalInterpolationVariant& other) =
default;
187 LocalInterpolationVariant(LocalInterpolationVariant&& other) =
default;
188 LocalInterpolationVariant& operator=(
const LocalInterpolationVariant& other) =
default;
189 LocalInterpolationVariant& operator=(LocalInterpolationVariant&& other) =
default;
191 template<
typename F,
typename C>
192 void interpolate (
const F& ff, std::vector<C>& out)
const
194 Impl::visitIf([&](
const auto* impl) { impl->interpolate(ff, out); }, impl_);
198 std::variant<std::monostate,
const Implementations*...> impl_;
232 template<
class... Implementations>
239 using LocalBasis = Impl::LocalBasisVariant<
typename Implementations::Traits::LocalBasisType...>;
240 using LocalCoefficients = Impl::LocalCoefficientsVariant<
typename Implementations::Traits::LocalCoefficientsType...>;
241 using LocalInterpolation = Impl::LocalInterpolationVariant<
typename Implementations::Traits::LocalInterpolationType...>;
247 [&](std::monostate&) {
248 localBasis_ = LocalBasis();
249 localCoefficients_ = LocalCoefficients();
250 localInterpolation_ = LocalInterpolation();
252 geometryType_ = GeometryType{};
253 }, [&](
auto&& impl) {
254 localBasis_ = LocalBasis(impl.localBasis());
255 localCoefficients_ = LocalCoefficients(impl.localCoefficients());
256 localInterpolation_ = LocalInterpolation(impl.localInterpolation());
258 geometryType_ = impl.type();
286 template<
class Implementation,
287 std::enable_if_t<std::disjunction<std::is_same<std::decay_t<Implementation>, Implementations>...>
::value,
int> = 0>
289 impl_(std::forward<Implementation>(impl))
307 impl_(std::move(other.impl_))
327 impl_ = std::move(other.impl_);
335 template<
class Implementation,
336 std::enable_if_t<std::disjunction<std::is_same<std::decay_t<Implementation>, Implementations>...>
::value,
int> = 0>
339 impl_ = std::forward<Implementation>(impl);
358 return localCoefficients_;
366 return localInterpolation_;
380 constexpr GeometryType
type()
const
382 return geometryType_;
406 operator bool ()
const
408 return not(std::holds_alternative<std::monostate>(
variant()));
412 std::variant<std::monostate, Implementations...> impl_;
414 GeometryType geometryType_;
415 LocalBasis localBasis_;
416 LocalCoefficients localCoefficients_;
417 LocalInterpolation localInterpolation_;
Definition: bdfmcube.hh:16
@ value
Definition: tensor.hh:166
Type traits for LocalBasisVirtualInterface.
Definition: common/localbasis.hh:32
traits helper struct
Definition: localfiniteelementtraits.hh:11
Type erasure class for wrapping LocalFiniteElement classes.
Definition: localfiniteelementvariant.hh:234
typename Dune::LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Export LocalFiniteElementTraits.
Definition: localfiniteelementvariant.hh:267
LocalFiniteElementVariant & operator=(const LocalFiniteElementVariant &other)
Copy assignment.
Definition: localfiniteelementvariant.hh:315
const Traits::LocalInterpolationType & localInterpolation() const
Provide access to LocalInterpolation implementation of this LocalFiniteElement.
Definition: localfiniteelementvariant.hh:364
unsigned int size() const
Number of shape functions.
Definition: localfiniteelementvariant.hh:372
constexpr GeometryType type() const
Number of shape functions.
Definition: localfiniteelementvariant.hh:380
LocalFiniteElementVariant(LocalFiniteElementVariant &&other)
Move constructor.
Definition: localfiniteelementvariant.hh:306
LocalFiniteElementVariant(const std::monostate &monostate)
Construct empty LocalFiniteElementVariant.
Definition: localfiniteelementvariant.hh:277
const Traits::LocalCoefficientsType & localCoefficients() const
Provide access to LocalCoefficients implementation of this LocalFiniteElement.
Definition: localfiniteelementvariant.hh:356
LocalFiniteElementVariant(const LocalFiniteElementVariant &other)
Copy constructor.
Definition: localfiniteelementvariant.hh:297
LocalFiniteElementVariant(Implementation &&impl)
Construct LocalFiniteElementVariant.
Definition: localfiniteelementvariant.hh:288
const auto & variant() const
Provide access to underlying std::variant.
Definition: localfiniteelementvariant.hh:396
LocalFiniteElementVariant & operator=(LocalFiniteElementVariant &&other)
Move assignment.
Definition: localfiniteelementvariant.hh:325
const Traits::LocalBasisType & localBasis() const
Provide access to LocalBasis implementation of this LocalFiniteElement.
Definition: localfiniteelementvariant.hh:348
LocalFiniteElementVariant & operator=(Implementation &&impl)
Assignment from implementation.
Definition: localfiniteelementvariant.hh:337
LocalFiniteElementVariant()=default
Construct empty LocalFiniteElementVariant.
Describe position of one degree of freedom.
Definition: localkey.hh:21