3 #ifndef DUNE_VIRTUALINTERFACE_HH
4 #define DUNE_VIRTUALINTERFACE_HH
6 #include <dune/common/function.hh>
8 #include <dune/geometry/type.hh>
18 template<
class DomainType,
class RangeType>
38 typename T::DomainFieldType,
40 typename T::DomainType,
41 typename T::RangeFieldType,
43 typename T::RangeType,
44 typename T::JacobianType,
54 template<
class T,
int order>
59 typename T::DomainFieldType,
61 typename T::DomainType,
62 typename T::RangeFieldType,
64 typename T::RangeType,
65 typename T::JacobianType,
77 typedef typename FE::Traits::LocalBasisType::Traits::DomainType DomainType;
78 typedef typename FE::Traits::LocalBasisType::Traits::RangeType RangeType;
81 typedef typename FE::Traits::LocalInterpolationType Implementation;
120 class LocalBasisVirtualInterfaceBase :
128 virtual void evaluate (
129 const typename Dune::template array<int,Traits::diffOrder>& directions,
130 const typename Traits::DomainType& in,
131 std::vector<typename Traits::RangeType>& out)
const = 0;
133 using BaseInterface::evaluate;
143 template<
class DF,
int n,
class D,
class RF,
int m,
class R,
class J>
152 virtual unsigned int size ()
const = 0;
155 virtual unsigned int order ()
const = 0;
162 virtual void evaluateFunction (
const typename Traits::DomainType& in,
163 std::vector<typename Traits::RangeType>& out)
const = 0;
173 virtual void evaluateJacobian(
const typename Traits::DomainType& in,
174 std::vector<typename Traits::JacobianType>& out)
const = 0;
177 virtual void evaluate (
178 const typename Dune::template array<int,Traits::diffOrder>& directions,
179 const typename Traits::DomainType& in,
180 std::vector<typename Traits::RangeType>& out)
const = 0;
194 class LocalBasisVirtualInterface :
195 public virtual LocalBasisVirtualInterfaceBase<T>
197 typedef LocalBasisVirtualInterfaceBase<T> BaseInterface;
204 const typename Dune::template array<int,k>& directions,
206 std::vector<typename Traits::RangeType>& out)
const
208 typedef LocalBasisVirtualInterfaceBase<typename FixedOrderLocalBasisTraits<T,k>::Traits > OrderKBaseInterface;
209 const OrderKBaseInterface& asBase = *
this;
210 asBase.evaluate(directions, in, out);
213 using BaseInterface::size;
214 using BaseInterface::order;
215 using BaseInterface::evaluateFunction;
216 using BaseInterface::evaluateJacobian;
219 #ifndef __INTEL_COMPILER
220 using BaseInterface::evaluate;
243 template<
class DomainType,
class RangeType>
273 template<
class DomainType,
class RangeType>
274 class LocalInterpolationVirtualInterface
275 :
public LocalInterpolationVirtualInterfaceBase<DomainType, RangeType>
302 void interpolate (
const F& f, std::vector<CoefficientType>& out)
const
305 asBase.
interpolate(VirtualFunctionWrapper<F>(f),out);
308 template<
class F,
class C>
311 std::vector<CoefficientType> outDummy;
313 asBase.
interpolate(VirtualFunctionWrapper<F>(f),outDummy);
314 out.resize(outDummy.size());
315 for(
typename std::vector<CoefficientType>::size_type i=0; i<outDummy.size(); ++i)
316 out[i] = outDummy[i];
321 template <
typename F>
322 struct VirtualFunctionWrapper
326 VirtualFunctionWrapper(
const F &f)
330 virtual ~VirtualFunctionWrapper() {}
332 virtual void evaluate(
const DomainType& x, RangeType& y)
const
359 virtual std::size_t
size ()
const = 0;
389 typename T::DomainType,
396 using BaseInterface::localCoefficients;
397 using BaseInterface::localInterpolation;
398 using BaseInterface::type;
410 template<
class DF,
int n,
class D,
class RF,
int m,
class R,
class J>
435 virtual const GeometryType type ()
const = 0;
void interpolate(const F &f, std::vector< CoefficientType > &out) const
determine coefficients interpolating a given function
Definition: virtualinterface.hh:302
traits helper struct
Definition: localfiniteelementtraits.hh:10
LocalFiniteElementTraits< LocalBasisVirtualInterface< T >, LocalCoefficientsVirtualInterface, LocalInterpolationVirtualInterface< typename T::DomainType, typename T::RangeType > > Traits
Definition: virtualinterface.hh:390
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:39
Dune::VirtualFunction< DomainType, RangeType > FunctionType
type of virtual function to interpolate
Definition: virtualinterface.hh:280
VirtualFunction< DomainType, RangeType > VirtualFunctionBase
Definition: virtualinterface.hh:85
virtual std::size_t size() const =0
number of coefficients
virtual ~LocalInterpolationVirtualInterface()
Definition: virtualinterface.hh:286
virtual const Traits::LocalBasisType & localBasis() const =0
Return a proper base class for functions to use with LocalInterpolation.
Definition: virtualinterface.hh:75
conditional< IsBaseOf< Interface, Implementation >::value, VirtualFunctionBase, FunctionBase >::type type
Base class type for functions to use with LocalInterpolation.
Definition: virtualinterface.hh:93
R RangeType
range type
Definition: localbasis.hh:63
virtual base class for local finite elements with functions
Definition: virtualinterface.hh:379
RangeType::field_type CoefficientType
type of the coefficient vector in the interpolate method
Definition: virtualinterface.hh:252
virtual const LocalKey & localKey(std::size_t i) const =0
get i'th index
Describe position of one degree of freedom.
Definition: localkey.hh:20
virtual ~LocalBasisVirtualInterfaceBase()
Definition: virtualinterface.hh:149
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:45
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:66
virtual base class for local coefficients
Definition: virtualinterface.hh:352
LI LocalInterpolationType
Definition: localfiniteelementtraits.hh:22
virtual base class for a local basis
Definition: virtualinterface.hh:22
void evaluate(const typename Dune::template array< int, k > &directions, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Definition: virtualinterface.hh:203
Construct LocalBasisTraits with fixed diff order.
Definition: virtualinterface.hh:55
LC LocalCoefficientsType
Definition: localfiniteelementtraits.hh:18
Dune::VirtualFunction< DomainType, RangeType > FunctionType
type of virtual function to interpolate
Definition: virtualinterface.hh:249
virtual void interpolate(const FunctionType &f, std::vector< CoefficientType > &out) const =0
determine coefficients interpolating a given function
virtual LocalFiniteElementVirtualInterface< T > * clone() const =0
Definition: tensor.hh:165
RangeType::field_type CoefficientType
type of the coefficient vector in the interpolate method
Definition: virtualinterface.hh:283
virtual void interpolate(const FunctionType &f, std::vector< CoefficientType > &out) const =0
determine coefficients interpolating a given function
virtual ~LocalCoefficientsVirtualInterface()
Definition: virtualinterface.hh:356
void interpolate(const F &f, std::vector< C > &out) const
Definition: virtualinterface.hh:309
LB LocalBasisType
Definition: localfiniteelementtraits.hh:14
virtual base class for a local interpolation
Definition: virtualinterface.hh:244
D DomainType
domain type
Definition: localbasis.hh:51
virtual ~LocalInterpolationVirtualInterfaceBase()
Definition: virtualinterface.hh:254
LocalBasisTraits< DF, n, D, RF, m, R, J, 0 > Traits
Definition: virtualinterface.hh:147
virtual ~LocalFiniteElementVirtualInterface()
Definition: virtualinterface.hh:423
Construct LocalBasisTraits with one diff order lower.
Definition: virtualinterface.hh:34
T Traits
Definition: virtualinterface.hh:199
virtual base class for a local interpolation
Definition: virtualinterface.hh:19
Function< const DomainType &, RangeType & > FunctionBase
Definition: virtualinterface.hh:86
LocalFiniteElementTraits< LocalBasisVirtualInterface< T >, LocalCoefficientsVirtualInterface, LocalInterpolationVirtualInterface< typename T::DomainType, typename T::RangeType > > Traits
Definition: virtualinterface.hh:421