dune-localfunctions 3.0-git
virtualinterface.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#ifndef DUNE_VIRTUALINTERFACE_HH
4#define DUNE_VIRTUALINTERFACE_HH
5
6#include <array>
7
8#include <dune/common/function.hh>
9
10#include <dune/geometry/type.hh>
11
15
16namespace Dune
17{
18
19 // forward declaration needed by the helper traits
20 template<class DomainType, class RangeType>
21 class LocalInterpolationVirtualInterface;
22
23 template<class T>
24 class LocalBasisVirtualInterface;
25
26 // -----------------------------------------------------------------
27 // Helper traits classes
28 // -----------------------------------------------------------------
29
35 template<class T>
37 {
39 typedef LocalBasisTraits<
40 typename T::DomainFieldType,
41 T::dimDomain,
42 typename T::DomainType,
43 typename T::RangeFieldType,
44 T::dimRange,
45 typename T::RangeType,
46 typename T::JacobianType,
47 T::diffOrder-1> Traits;
48 };
49
56 template<class T, int order>
58 {
60 typedef LocalBasisTraits<
61 typename T::DomainFieldType,
62 T::dimDomain,
63 typename T::DomainType,
64 typename T::RangeFieldType,
65 T::dimRange,
66 typename T::RangeType,
67 typename T::JacobianType,
68 order> Traits;
69 };
70
76 template<class FE>
78 {
79 typedef typename FE::Traits::LocalBasisType::Traits::DomainType DomainType;
80 typedef typename FE::Traits::LocalBasisType::Traits::RangeType RangeType;
81
83 typedef typename FE::Traits::LocalInterpolationType Implementation;
84
85 public:
86
87 typedef VirtualFunction<DomainType, RangeType> VirtualFunctionBase;
88 typedef Function<const DomainType&, RangeType&> FunctionBase;
89
95 typedef typename std::conditional<std::is_base_of<Interface, Implementation>::value, VirtualFunctionBase, FunctionBase>::type type;
96 };
97
98
99
100 // -----------------------------------------------------------------
101 // Basis
102 // -----------------------------------------------------------------
103
104 // current versions of doxygen (<= 1.6.2) enter an infinite loop when parsing
105 // the following class
106#ifndef DOXYGEN
121 template<class T>
122 class LocalBasisVirtualInterfaceBase :
123 public virtual LocalBasisVirtualInterface<typename LowerOrderLocalBasisTraits<T>::Traits>
124 {
126 public:
127 typedef T Traits;
128
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;
134
135 using BaseInterface::evaluate;
136 };
137#endif // DOXYGEN
138
145 template<class DF, int n, class D, class RF, int m, class R, class J>
146 class LocalBasisVirtualInterfaceBase<LocalBasisTraits<DF,n,D,RF,m,R,J,0> >
147 {
148 public:
150
152
154 virtual unsigned int size () const = 0;
155
157 virtual unsigned int order () const = 0;
158
164 virtual void evaluateFunction (const typename Traits::DomainType& in,
165 std::vector<typename Traits::RangeType>& out) const = 0;
166
175 virtual void evaluateJacobian(const typename Traits::DomainType& in, // position
176 std::vector<typename Traits::JacobianType>& out) const = 0;
177
183 virtual void partial(const std::array<unsigned int,n>& order,
184 const typename Traits::DomainType& in,
185 std::vector<typename Traits::RangeType>& out) const = 0;
186
188 virtual void evaluate (
189 const typename std::template array<int,Traits::diffOrder>& directions,
190 const typename Traits::DomainType& in,
191 std::vector<typename Traits::RangeType>& out) const = 0;
192
193 };
194
204 template<class T>
206 public virtual LocalBasisVirtualInterfaceBase<T>
207 {
208 typedef LocalBasisVirtualInterfaceBase<T> BaseInterface;
209 public:
210 typedef T Traits;
211
213 template <int k>
214 void evaluate (
215 const typename std::template array<int,k>& directions,
216 const typename Traits::DomainType& in,
217 std::vector<typename Traits::RangeType>& out) const
218 {
219 typedef LocalBasisVirtualInterfaceBase<typename FixedOrderLocalBasisTraits<T,k>::Traits > OrderKBaseInterface;
220 const OrderKBaseInterface& asBase = *this;
221 asBase.evaluate(directions, in, out);
222 }
223
224 using BaseInterface::size;
225 using BaseInterface::order;
226 using BaseInterface::partial;
227 using BaseInterface::evaluateFunction;
228 using BaseInterface::evaluateJacobian;
229 /* Unfortunately, the intel compiler cannot use the different evaluate
230 * methods with varying argument lists. :-( */
231#ifndef __INTEL_COMPILER
232 using BaseInterface::evaluate;
233#endif
234 };
235
236
237
238
239 // -----------------------------------------------------------------
240 // Interpolation
241 // -----------------------------------------------------------------
242
255 template<class DomainType, class RangeType>
257 {
258 public:
259
261 typedef Dune::VirtualFunction<DomainType, RangeType> FunctionType;
262
264 typedef typename RangeType::field_type CoefficientType;
265
267
275 virtual void interpolate (const FunctionType& f, std::vector<CoefficientType>& out) const = 0;
276 };
277
285 template<class DomainType, class RangeType>
287 : public LocalInterpolationVirtualInterfaceBase<DomainType, RangeType>
288 {
289 public:
290
292 typedef Dune::VirtualFunction<DomainType, RangeType> FunctionType;
293
295 typedef typename RangeType::field_type CoefficientType;
296
297
299
300 // This method is only notet again for to make the documentation complete.
301
309 virtual void interpolate (const FunctionType& f, std::vector<CoefficientType>& out) const = 0;
310
313 template<class F>
314 void interpolate (const F& f, std::vector<CoefficientType>& out) const
315 {
317 asBase.interpolate(VirtualFunctionWrapper<F>(f),out);
318 }
319
320 template<class F, class C>
321 void interpolate (const F& f, std::vector<C>& out) const
322 {
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];
329 }
330
331 private:
332
333 template <typename F>
334 struct VirtualFunctionWrapper
335 : public FunctionType
336 {
337 public:
338 VirtualFunctionWrapper(const F &f)
339 : f_(f)
340 {}
341
342 virtual ~VirtualFunctionWrapper() {}
343
344 virtual void evaluate(const DomainType& x, RangeType& y) const
345 {
346 f_.evaluate(x,y);
347 }
348
349 const F &f_;
350 };
351 };
352
353
354
355 // -----------------------------------------------------------------
356 // Coefficients
357 // -----------------------------------------------------------------
358
365 {
366 public:
367
369
371 virtual std::size_t size () const = 0;
372
374 const virtual LocalKey& localKey (std::size_t i) const = 0;
375
376 };
377
378
379
380 // -----------------------------------------------------------------
381 // Finite Element
382 // -----------------------------------------------------------------
383
390 template<class T>
392 : public virtual LocalFiniteElementVirtualInterface<typename LowerOrderLocalBasisTraits<T>::Traits >
393 {
395
396 public:
401 typename T::DomainType,
402 typename T::RangeType> > Traits;
403
405 virtual const typename Traits::LocalBasisType& localBasis () const = 0;
406
408 using BaseInterface::localCoefficients;
409 using BaseInterface::localInterpolation;
410 using BaseInterface::type;
411
413 };
414
415
422 template<class DF, int n, class D, class RF, int m, class R, class J>
424 {
426
427 public:
432 typename T::DomainType,
433 typename T::RangeType> > Traits;
434
436
438 virtual const typename Traits::LocalBasisType& localBasis () const = 0;
439
441 virtual const typename Traits::LocalCoefficientsType& localCoefficients () const = 0;
442
444 virtual const typename Traits::LocalInterpolationType& localInterpolation () const = 0;
445
447 virtual unsigned int size () const = 0;
448
450 virtual const GeometryType type () const = 0;
451
453 };
454
455}
456#endif
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 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 Traits::LocalCoefficientsType & localCoefficients() const =0
virtual const Traits::LocalInterpolationType & localInterpolation() const =0
LocalFiniteElementTraits< LocalBasisVirtualInterface< T >, LocalCoefficientsVirtualInterface, LocalInterpolationVirtualInterface< typename T::DomainType, typename T::RangeType > > Traits
Definition virtualinterface.hh:433