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>
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,
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;
179 virtual void evaluate (
180 const typename std::template array<int,Traits::diffOrder>& directions,
182 std::vector<typename Traits::RangeType>& out)
const = 0;
197 public virtual LocalBasisVirtualInterfaceBase<T>
199 typedef LocalBasisVirtualInterfaceBase<T> BaseInterface;
206 const typename std::template array<int,k>& directions,
208 std::vector<typename Traits::RangeType>& out)
const 210 typedef LocalBasisVirtualInterfaceBase<typename FixedOrderLocalBasisTraits<T,k>::Traits > OrderKBaseInterface;
211 const OrderKBaseInterface& asBase = *
this;
212 asBase.evaluate(directions, in, out);
215 using BaseInterface::size;
216 using BaseInterface::order;
217 using BaseInterface::evaluateFunction;
218 using BaseInterface::evaluateJacobian;
221 #ifndef __INTEL_COMPILER 222 using BaseInterface::evaluate;
245 template<
class DomainType,
class RangeType>
265 virtual void interpolate (
const FunctionType& f, std::vector<CoefficientType>& out)
const = 0;
275 template<
class DomainType,
class RangeType>
299 virtual void interpolate (
const FunctionType& f, std::vector<CoefficientType>& out)
const = 0;
304 void interpolate (
const F& f, std::vector<CoefficientType>& out)
const 307 asBase.
interpolate(VirtualFunctionWrapper<F>(f),out);
310 template<
class F,
class C>
313 std::vector<CoefficientType> outDummy;
315 asBase.
interpolate(VirtualFunctionWrapper<F>(f),outDummy);
316 out.resize(outDummy.size());
317 for(
typename std::vector<CoefficientType>::size_type i=0; i<outDummy.size(); ++i)
318 out[i] = outDummy[i];
323 template <
typename F>
324 struct VirtualFunctionWrapper
325 :
public FunctionType
328 VirtualFunctionWrapper(
const F &f)
332 virtual ~VirtualFunctionWrapper() {}
334 virtual void evaluate(
const DomainType& x, RangeType& y)
const 361 virtual std::size_t size ()
const = 0;
364 const virtual LocalKey& localKey (std::size_t i)
const = 0;
391 typename T::DomainType,
397 using BaseInterface::localBasis;
398 using BaseInterface::localCoefficients;
399 using BaseInterface::localInterpolation;
400 using BaseInterface::type;
412 template<
class DF,
int n,
class D,
class RF,
int m,
class R,
class J>
437 virtual unsigned int size ()
const = 0;
440 virtual const GeometryType type ()
const = 0;
Definition: tensor.hh:165
void interpolate(const F &f, std::vector< CoefficientType > &out) const
determine coefficients interpolating a given function
Definition: virtualinterface.hh:304
traits helper struct
Definition: localfiniteelementtraits.hh:10
LocalFiniteElementTraits< LocalBasisVirtualInterface< T >, LocalCoefficientsVirtualInterface, LocalInterpolationVirtualInterface< typename T::DomainType, typename T::RangeType > > Traits
Definition: virtualinterface.hh:392
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
Dune::VirtualFunction< DomainType, RangeType > FunctionType
type of virtual function to interpolate
Definition: virtualinterface.hh:282
VirtualFunction< DomainType, RangeType > VirtualFunctionBase
Definition: virtualinterface.hh:87
virtual ~LocalInterpolationVirtualInterface()
Definition: virtualinterface.hh:288
Return a proper base class for functions to use with LocalInterpolation.
Definition: virtualinterface.hh:77
conditional< IsBaseOf< Interface, Implementation >::value, VirtualFunctionBase, FunctionBase >::type type
Base class type for functions to use with LocalInterpolation.
Definition: virtualinterface.hh:95
R RangeType
range type
Definition: localbasis.hh:61
virtual base class for local finite elements with functions
Definition: virtualinterface.hh:381
RangeType::field_type CoefficientType
type of the coefficient vector in the interpolate method
Definition: virtualinterface.hh:254
Describe position of one degree of freedom.
Definition: localkey.hh:21
virtual ~LocalBasisVirtualInterfaceBase()
Definition: virtualinterface.hh:151
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
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
virtual base class for local coefficients
Definition: virtualinterface.hh:354
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
LI LocalInterpolationType
Definition: localfiniteelementtraits.hh:22
virtual base class for a local basis
Definition: virtualinterface.hh:24
Construct LocalBasisTraits with fixed diff order.
Definition: virtualinterface.hh:57
LC LocalCoefficientsType
Definition: localfiniteelementtraits.hh:18
Dune::VirtualFunction< DomainType, RangeType > FunctionType
type of virtual function to interpolate
Definition: virtualinterface.hh:251
RangeType::field_type CoefficientType
type of the coefficient vector in the interpolate method
Definition: virtualinterface.hh:285
virtual void interpolate(const FunctionType &f, std::vector< CoefficientType > &out) const =0
determine coefficients interpolating a given function
virtual ~LocalCoefficientsVirtualInterface()
Definition: virtualinterface.hh:358
void interpolate(const F &f, std::vector< C > &out) const
Definition: virtualinterface.hh:311
LB LocalBasisType
Definition: localfiniteelementtraits.hh:14
virtual base class for a local interpolation
Definition: virtualinterface.hh:246
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:205
D DomainType
domain type
Definition: localbasis.hh:49
virtual ~LocalInterpolationVirtualInterfaceBase()
Definition: virtualinterface.hh:256
LocalBasisTraits< DF, n, D, RF, m, R, J, 0 > Traits
Definition: virtualinterface.hh:149
virtual ~LocalFiniteElementVirtualInterface()
Definition: virtualinterface.hh:425
Construct LocalBasisTraits with one diff order lower.
Definition: virtualinterface.hh:36
T Traits
Definition: virtualinterface.hh:201
virtual base class for a local interpolation
Definition: virtualinterface.hh:21
Function< const DomainType &, RangeType & > FunctionBase
Definition: virtualinterface.hh:88
LocalFiniteElementTraits< LocalBasisVirtualInterface< T >, LocalCoefficientsVirtualInterface, LocalInterpolationVirtualInterface< typename T::DomainType, typename T::RangeType > > Traits
Definition: virtualinterface.hh:423