dune-localfunctions  2.2.1
virtualinterface.hh
Go to the documentation of this file.
1 // -*- tab-width: 8; indent-tabs-mode: nil -*-
2 // vi: set ts=8 sw=2 et sts=2:
3 #ifndef DUNE_VIRTUALINTERFACE_HH
4 #define DUNE_VIRTUALINTERFACE_HH
5 
6 #include<dune/common/function.hh>
7 
8 #include<dune/geometry/type.hh>
9 
13 
14 namespace Dune
15 {
16 
17  // forward declaration needed by the helper traits
18  template<class DomainType, class RangeType>
20 
21  template<class T>
23 
24  // -----------------------------------------------------------------
25  // Helper traits classes
26  // -----------------------------------------------------------------
27 
33  template<class T>
35  {
37  typedef LocalBasisTraits<
38  typename T::DomainFieldType,
39  T::dimDomain,
40  typename T::DomainType,
41  typename T::RangeFieldType,
42  T::dimRange,
43  typename T::RangeType,
44  typename T::JacobianType,
45  T::diffOrder-1> Traits;
46  };
47 
54  template<class T, int order>
56  {
58  typedef LocalBasisTraits<
59  typename T::DomainFieldType,
60  T::dimDomain,
61  typename T::DomainType,
62  typename T::RangeFieldType,
63  T::dimRange,
64  typename T::RangeType,
65  typename T::JacobianType,
66  order> Traits;
67  };
68 
74  template<class FE>
76  {
77  typedef typename FE::Traits::LocalBasisType::Traits::DomainType DomainType;
78  typedef typename FE::Traits::LocalBasisType::Traits::RangeType RangeType;
79 
81  typedef typename FE::Traits::LocalInterpolationType Implementation;
82 
83  public:
84 
85  typedef VirtualFunction<DomainType, RangeType> VirtualFunctionBase;
86  typedef Function<const DomainType&, RangeType&> FunctionBase;
87 
94  };
95 
96 
97 
98  // -----------------------------------------------------------------
99  // Basis
100  // -----------------------------------------------------------------
101 
102 // current versions of doxygen (<= 1.6.2) enter an infinite loop when parsing
103 // the following class
104 #ifndef DOXYGEN
105 
119  template<class T>
120  class LocalBasisVirtualInterfaceBase :
121  public virtual LocalBasisVirtualInterface<typename LowerOrderLocalBasisTraits<T>::Traits>
122  {
124  public:
125  typedef T Traits;
126 
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;
132 
133  using BaseInterface::evaluate;
134  };
135 #endif // DOXYGEN
136 
143  template<class DF, int n, class D, class RF, int m, class R, class J>
144  class LocalBasisVirtualInterfaceBase<LocalBasisTraits<DF,n,D,RF,m,R,J,0> >
145  {
146  public:
148 
150 
152  virtual unsigned int size () const = 0;
153 
155  virtual unsigned int order () const = 0;
156 
162  virtual void evaluateFunction (const typename Traits::DomainType& in,
163  std::vector<typename Traits::RangeType>& out) const = 0;
164 
173  virtual void evaluateJacobian(const typename Traits::DomainType& in, // position
174  std::vector<typename Traits::JacobianType>& out) const = 0;
175 
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;
181 
182  };
183 
193  template<class T>
194  class LocalBasisVirtualInterface :
195  public virtual LocalBasisVirtualInterfaceBase<T>
196  {
197  typedef LocalBasisVirtualInterfaceBase<T> BaseInterface;
198  public:
199  typedef T Traits;
200 
202  template <int k>
203  void evaluate (
204  const typename Dune::template array<int,k>& directions,
205  const typename Traits::DomainType& in,
206  std::vector<typename Traits::RangeType>& out) const
207  {
208  typedef LocalBasisVirtualInterfaceBase<typename FixedOrderLocalBasisTraits<T,k>::Traits > OrderKBaseInterface;
209  const OrderKBaseInterface& asBase = *this;
210  asBase.evaluate(directions, in, out);
211  }
212 
213  using BaseInterface::size;
214  using BaseInterface::order;
215  using BaseInterface::evaluateFunction;
216  using BaseInterface::evaluateJacobian;
217  /* Unfortunately, the intel compiler cannot use the different evaluate
218  * methods with varying argument lists. :-( */
219 #ifndef __INTEL_COMPILER
220  using BaseInterface::evaluate;
221 #endif
222  };
223 
224 
225 
226 
227  // -----------------------------------------------------------------
228  // Interpolation
229  // -----------------------------------------------------------------
230 
243  template<class DomainType, class RangeType>
245  {
246  public:
247 
249  typedef Dune::VirtualFunction<DomainType, RangeType> FunctionType;
250 
252  typedef typename RangeType::field_type CoefficientType;
253 
255 
263  virtual void interpolate (const FunctionType& f, std::vector<CoefficientType>& out) const = 0;
264  };
265 
273  template<class DomainType, class RangeType>
274  class LocalInterpolationVirtualInterface
275  : public LocalInterpolationVirtualInterfaceBase<DomainType, RangeType>
276  {
277  public:
278 
280  typedef Dune::VirtualFunction<DomainType, RangeType> FunctionType;
281 
283  typedef typename RangeType::field_type CoefficientType;
284 
285 
287 
288  // This method is only notet again for to make the documentation complete.
289 
297  virtual void interpolate (const FunctionType& f, std::vector<CoefficientType>& out) const = 0;
298 
301  template<class F>
302  void interpolate (const F& f, std::vector<CoefficientType>& out) const
303  {
305  asBase.interpolate(VirtualFunctionWrapper<F>(f),out);
306  }
307 
308  template<class F, class C>
309  void interpolate (const F& f, std::vector<C>& out) const
310  {
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];
317  }
318 
319  private:
320 
321  template <typename F>
322  struct VirtualFunctionWrapper
323  : public FunctionType
324  {
325  public:
326  VirtualFunctionWrapper(const F &f)
327  : f_(f)
328  {}
329 
330  virtual ~VirtualFunctionWrapper() {}
331 
332  virtual void evaluate(const DomainType& x, RangeType& y) const
333  {
334  f_.evaluate(x,y);
335  }
336 
337  const F &f_;
338  };
339  };
340 
341 
342 
343  // -----------------------------------------------------------------
344  // Coefficients
345  // -----------------------------------------------------------------
346 
353  {
354  public:
355 
357 
359  virtual std::size_t size () const = 0;
360 
362  const virtual LocalKey& localKey (std::size_t i) const = 0;
363 
364  };
365 
366 
367 
368  // -----------------------------------------------------------------
369  // Finite Element
370  // -----------------------------------------------------------------
371 
378  template<class T>
380  : public virtual LocalFiniteElementVirtualInterface<typename LowerOrderLocalBasisTraits<T>::Traits >
381  {
383 
384  public:
385  typedef LocalFiniteElementTraits<
389  typename T::DomainType,
390  typename T::RangeType> > Traits;
391 
393  virtual const typename Traits::LocalBasisType& localBasis () const = 0;
394 
396  using BaseInterface::localCoefficients;
397  using BaseInterface::localInterpolation;
398  using BaseInterface::type;
399 
400  virtual LocalFiniteElementVirtualInterface<T>* clone() const = 0;
401  };
402 
403 
410  template<class DF, int n, class D, class RF, int m, class R, class J>
412  {
414 
415  public:
416  typedef LocalFiniteElementTraits<
420  typename T::DomainType,
421  typename T::RangeType> > Traits;
422 
424 
426  virtual const typename Traits::LocalBasisType& localBasis () const = 0;
427 
429  virtual const typename Traits::LocalCoefficientsType& localCoefficients () const = 0;
430 
432  virtual const typename Traits::LocalInterpolationType& localInterpolation () const = 0;
433 
435  virtual const GeometryType type () const = 0;
436 
437  virtual LocalFiniteElementVirtualInterface<T>* clone() const = 0;
438  };
439 
440 }
441 #endif