3#ifndef DUNE_VIRTUALINTERFACE_HH
4#define DUNE_VIRTUALINTERFACE_HH
8#include <dune/common/function.hh>
10#include <dune/geometry/type.hh>
20 template<
class DomainType,
class RangeType>
21 class LocalInterpolationVirtualInterface;
24 class LocalBasisVirtualInterface;
40 typename T::DomainFieldType,
42 typename T::DomainType,
43 typename T::RangeFieldType,
45 typename T::RangeType,
46 typename T::JacobianType,
56 template<
class T,
int order>
61 typename T::DomainFieldType,
63 typename T::DomainType,
64 typename T::RangeFieldType,
66 typename T::RangeType,
67 typename T::JacobianType,
79 typedef typename FE::Traits::LocalBasisType::Traits::DomainType DomainType;
80 typedef typename FE::Traits::LocalBasisType::Traits::RangeType RangeType;
83 typedef typename FE::Traits::LocalInterpolationType Implementation;
122 class LocalBasisVirtualInterfaceBase :
130 virtual void evaluate (
131 const typename std::template array<int,Traits::diffOrder>& directions,
132 const typename Traits::DomainType& in,
133 std::vector<typename Traits::RangeType>& out)
const = 0;
135 using BaseInterface::evaluate;
145 template<
class DF,
int n,
class D,
class RF,
int m,
class R,
class J>
154 virtual unsigned int size ()
const = 0;
157 virtual unsigned int order ()
const = 0;
165 std::vector<typename Traits::RangeType>& out)
const = 0;
176 std::vector<typename Traits::JacobianType>& out)
const = 0;
183 virtual void partial(
const std::array<unsigned int,n>& order,
185 std::vector<typename Traits::RangeType>& out)
const = 0;
189 const typename std::template array<int,Traits::diffOrder>& directions,
191 std::vector<typename Traits::RangeType>& out)
const = 0;
206 public virtual LocalBasisVirtualInterfaceBase<T>
208 typedef LocalBasisVirtualInterfaceBase<T> BaseInterface;
215 const typename std::template array<int,k>& directions,
216 const typename Traits::DomainType& in,
217 std::vector<typename Traits::RangeType>& out)
const
219 typedef LocalBasisVirtualInterfaceBase<typename FixedOrderLocalBasisTraits<T,k>::Traits > OrderKBaseInterface;
220 const OrderKBaseInterface& asBase = *
this;
221 asBase.evaluate(directions, in, out);
224 using BaseInterface::size;
225 using BaseInterface::order;
226 using BaseInterface::partial;
227 using BaseInterface::evaluateFunction;
228 using BaseInterface::evaluateJacobian;
231#ifndef __INTEL_COMPILER
232 using BaseInterface::evaluate;
255 template<
class DomainType,
class RangeType>
285 template<
class DomainType,
class RangeType>
314 void interpolate (
const F& f, std::vector<CoefficientType>& out)
const
317 asBase.
interpolate(VirtualFunctionWrapper<F>(f),out);
320 template<
class F,
class C>
323 std::vector<CoefficientType> outDummy;
325 asBase.
interpolate(VirtualFunctionWrapper<F>(f),outDummy);
326 out.resize(outDummy.size());
327 for(
typename std::vector<CoefficientType>::size_type i=0; i<outDummy.size(); ++i)
328 out[i] = outDummy[i];
333 template <
typename F>
334 struct VirtualFunctionWrapper
338 VirtualFunctionWrapper(
const F &f)
342 virtual ~VirtualFunctionWrapper() {}
344 virtual void evaluate(
const DomainType& x, RangeType& y)
const
371 virtual std::size_t
size ()
const = 0;
401 typename T::DomainType,
408 using BaseInterface::localCoefficients;
409 using BaseInterface::localInterpolation;
410 using BaseInterface::type;
422 template<
class DF,
int n,
class D,
class RF,
int m,
class R,
class J>
447 virtual unsigned int size ()
const = 0;
450 virtual const GeometryType
type ()
const = 0;
Definition brezzidouglasmarini1cube2d.hh:14
Type traits for LocalBasisVirtualInterface.
Definition localbasis.hh:38
R RangeType
range type
Definition localbasis.hh:61
D DomainType
domain type
Definition localbasis.hh:49
traits helper struct
Definition localfiniteelementtraits.hh:11
LB LocalBasisType
Definition localfiniteelementtraits.hh:14
LC LocalCoefficientsType
Definition localfiniteelementtraits.hh:18
LI LocalInterpolationType
Definition localfiniteelementtraits.hh:22
Describe position of one degree of freedom.
Definition localkey.hh:21
virtual base class for a local interpolation
Definition virtualinterface.hh:288
virtual ~LocalInterpolationVirtualInterface()
Definition virtualinterface.hh:298
void interpolate(const F &f, std::vector< C > &out) const
Definition virtualinterface.hh:321
RangeType::field_type CoefficientType
type of the coefficient vector in the interpolate method
Definition virtualinterface.hh:295
Dune::VirtualFunction< DomainType, RangeType > FunctionType
type of virtual function to interpolate
Definition virtualinterface.hh:292
virtual void interpolate(const FunctionType &f, std::vector< CoefficientType > &out) const =0
determine coefficients interpolating a given function
void interpolate(const F &f, std::vector< CoefficientType > &out) const
determine coefficients interpolating a given function
Definition virtualinterface.hh:314
virtual base class for a local basis
Definition virtualinterface.hh:207
void evaluate(const typename std::template array< int, k > &directions, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Definition virtualinterface.hh:214
T Traits
Definition virtualinterface.hh:210
Construct LocalBasisTraits with one diff order lower.
Definition virtualinterface.hh:37
LocalBasisTraits< typename T::DomainFieldType, T::dimDomain, typename T::DomainType, typename T::RangeFieldType, T::dimRange, typename T::RangeType, typename T::JacobianType, T::diffOrder-1 > Traits
The LocalBasisTraits with one order lower.
Definition virtualinterface.hh:47
Construct LocalBasisTraits with fixed diff order.
Definition virtualinterface.hh:58
LocalBasisTraits< typename T::DomainFieldType, T::dimDomain, typename T::DomainType, typename T::RangeFieldType, T::dimRange, typename T::RangeType, typename T::JacobianType, order > Traits
The LocalBasisTraits specified order.
Definition virtualinterface.hh:68
Return a proper base class for functions to use with LocalInterpolation.
Definition virtualinterface.hh:78
std::conditional< std::is_base_of< Interface, Implementation >::value, VirtualFunctionBase, FunctionBase >::type type
Base class type for functions to use with LocalInterpolation.
Definition virtualinterface.hh:95
Function< const DomainType &, RangeType & > FunctionBase
Definition virtualinterface.hh:88
VirtualFunction< DomainType, RangeType > VirtualFunctionBase
Definition virtualinterface.hh:87
LocalBasisTraits< DF, n, D, RF, m, R, J, 0 > Traits
Definition virtualinterface.hh:149
virtual ~LocalBasisVirtualInterfaceBase()
Definition virtualinterface.hh:151
virtual unsigned int size() const =0
Number of shape functions.
virtual void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const =0
Evaluate jacobian of all shape functions at given position.
virtual void partial(const std::array< unsigned int, n > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const =0
Evaluate partial derivatives of any order of all shape functions.
virtual unsigned int order() const =0
Polynomial order of the shape functions.
virtual void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const =0
Evaluate all basis function at given position.
virtual void evaluate(const typename std::template array< int, Traits::diffOrder > &directions, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const =0
virtual base class for a local interpolation
Definition virtualinterface.hh:257
Dune::VirtualFunction< DomainType, RangeType > FunctionType
type of virtual function to interpolate
Definition virtualinterface.hh:261
RangeType::field_type CoefficientType
type of the coefficient vector in the interpolate method
Definition virtualinterface.hh:264
virtual ~LocalInterpolationVirtualInterfaceBase()
Definition virtualinterface.hh:266
virtual void interpolate(const FunctionType &f, std::vector< CoefficientType > &out) const =0
determine coefficients interpolating a given function
virtual base class for local coefficients
Definition virtualinterface.hh:365
virtual ~LocalCoefficientsVirtualInterface()
Definition virtualinterface.hh:368
virtual std::size_t size() const =0
number of coefficients
virtual const LocalKey & localKey(std::size_t i) const =0
get i'th index
virtual base class for local finite elements with functions
Definition virtualinterface.hh:393
virtual LocalFiniteElementVirtualInterface< T > * clone() const =0
virtual const Traits::LocalBasisType & localBasis() const =0
LocalFiniteElementTraits< LocalBasisVirtualInterface< T >, LocalCoefficientsVirtualInterface, LocalInterpolationVirtualInterface< typename T::DomainType, typename T::RangeType > > Traits
Definition virtualinterface.hh:402
virtual LocalFiniteElementVirtualInterface< T > * clone() const =0
virtual const GeometryType type() const =0
virtual unsigned int size() const =0
virtual const Traits::LocalCoefficientsType & localCoefficients() const =0
virtual const Traits::LocalBasisType & localBasis() const =0
virtual ~LocalFiniteElementVirtualInterface()
Definition virtualinterface.hh:435
virtual const Traits::LocalInterpolationType & localInterpolation() const =0
LocalFiniteElementTraits< LocalBasisVirtualInterface< T >, LocalCoefficientsVirtualInterface, LocalInterpolationVirtualInterface< typename T::DomainType, typename T::RangeType > > Traits
Definition virtualinterface.hh:433