dune-localfunctions  2.8.0
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_LOCALFUNCTIONS_COMMON_VIRTUALINTERFACE_HH
4 #define DUNE_LOCALFUNCTIONS_COMMON_VIRTUALINTERFACE_HH
5 
6 #include <type_traits>
7 #include <array>
8 #include <vector>
9 #include <functional>
10 
11 #include <dune/geometry/type.hh>
12 
17 
18 namespace Dune
19 {
20 
21  // forward declaration needed by the helper traits
22  template<class DomainType, class RangeType>
23  class LocalInterpolationVirtualInterface;
24 
25  // -----------------------------------------------------------------
26  // Helper traits classes
27  // -----------------------------------------------------------------
28 
41  template<class FE>
42  class
43  [[deprecated("Dune::LocalFiniteElementFunctionBase is deprecated after Dune 2.7. You can now pass functions providing operator() to interpolate.")]]
45  {
46  typedef typename FE::Traits::LocalBasisType::Traits::DomainType Domain;
47  typedef typename FE::Traits::LocalBasisType::Traits::RangeType Range;
48 
49  // Hack: Keep a copy of Dune::Function here. This allows to avoid depending
50  // on the deprecated dune-common header while still keeping the LocalFiniteElementFunctionBase
51  // mechanism working during its deprecation period.
52  class FunctionBaseDummy
53  {
54  public:
55 
56  using RangeType = Range;
57  using DomainType = Domain;
58 
59  struct Traits
60  {
61  using RangeType = Range;
62  using DomainType = Domain;
63  };
64 
65  void evaluate(const DomainType& x, RangeType& y) const;
66  };
67 
68  public:
69 
70  using VirtualFunctionBase = FunctionBaseDummy;
71  using FunctionBase = FunctionBaseDummy;
72 
78  using type = FunctionBaseDummy;
79  };
80 
81 
82 
83  // -----------------------------------------------------------------
84  // Basis
85  // -----------------------------------------------------------------
86 
93  template<class T>
95  {
96  public:
97  using Traits = T;
98 
99 
101 
103  virtual unsigned int size () const = 0;
104 
106  virtual unsigned int order () const = 0;
107 
113  virtual void evaluateFunction (const typename Traits::DomainType& in,
114  std::vector<typename Traits::RangeType>& out) const = 0;
115 
124  virtual void evaluateJacobian(const typename Traits::DomainType& in, // position
125  std::vector<typename Traits::JacobianType>& out) const = 0;
126 
132  virtual void partial(const std::array<unsigned int,Traits::dimDomain>& order,
133  const typename Traits::DomainType& in,
134  std::vector<typename Traits::RangeType>& out) const = 0;
135  };
136 
137 
138 
139  // -----------------------------------------------------------------
140  // Interpolation
141  // -----------------------------------------------------------------
142 
155  template<class DomainType, class RangeType>
157  {
158  public:
159 
161  using FunctionType = std::function<RangeType(DomainType)>;
162 
164  typedef typename RangeType::field_type CoefficientType;
165 
167 
175  virtual void interpolate (const FunctionType& f, std::vector<CoefficientType>& out) const = 0;
176  };
177 
185  template<class DomainType, class RangeType>
187  : public LocalInterpolationVirtualInterfaceBase<DomainType, RangeType>
188  {
189  public:
190 
192  using FunctionType = std::function<RangeType(DomainType)>;
193 
195  typedef typename RangeType::field_type CoefficientType;
196 
197 
199 
200  // This method is only noted again for to make the documentation complete.
201 
209  virtual void interpolate (const FunctionType& f, std::vector<CoefficientType>& out) const = 0;
210 
216  template<class F,
218  void interpolate (const F& ff, std::vector<CoefficientType>& out) const
219  {
220  const auto& f = Impl::makeFunctionWithCallOperator<DomainType>(ff);
221 
223  asBase.interpolate(FunctionType(std::cref(f)),out);
224  }
225 
231  template<class F, class C>
232  void interpolate (const F& ff, std::vector<C>& out) const
233  {
234  const auto& f = Impl::makeFunctionWithCallOperator<DomainType>(ff);
235 
236  std::vector<CoefficientType> outDummy;
238  asBase.interpolate(FunctionType(std::cref(f)),outDummy);
239  out.resize(outDummy.size());
240  for(typename std::vector<CoefficientType>::size_type i=0; i<outDummy.size(); ++i)
241  out[i] = outDummy[i];
242  }
243  };
244 
245 
246 
247  // -----------------------------------------------------------------
248  // Coefficients
249  // -----------------------------------------------------------------
250 
257  {
258  public:
259 
261 
263  virtual std::size_t size () const = 0;
264 
266  const virtual LocalKey& localKey (std::size_t i) const = 0;
267 
268  };
269 
270 
271 
272  // -----------------------------------------------------------------
273  // Finite Element
274  // -----------------------------------------------------------------
275 
276 
282  template<class T>
284  {
285  using LocalBasisTraits = T;
286  public:
287  typedef LocalFiniteElementTraits<
293 
295 
297  virtual const typename Traits::LocalBasisType& localBasis () const = 0;
298 
300  virtual const typename Traits::LocalCoefficientsType& localCoefficients () const = 0;
301 
303  virtual const typename Traits::LocalInterpolationType& localInterpolation () const = 0;
304 
306  virtual unsigned int size () const = 0;
307 
309  virtual const GeometryType type () const = 0;
310 
312  };
313 }
314 #endif
Definition: bdfmcube.hh:16
@ value
Definition: tensor.hh:166
D DomainType
domain type
Definition: common/localbasis.hh:43
R RangeType
range type
Definition: common/localbasis.hh:55
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:188
virtual ~LocalInterpolationVirtualInterface()
Definition: virtualinterface.hh:198
std::function< RangeType(DomainType)> FunctionType
type of function to interpolate
Definition: virtualinterface.hh:192
RangeType::field_type CoefficientType
type of the coefficient vector in the interpolate method
Definition: virtualinterface.hh:195
void interpolate(const F &ff, std::vector< CoefficientType > &out) const
determine coefficients interpolating a given function
Definition: virtualinterface.hh:218
virtual void interpolate(const FunctionType &f, std::vector< CoefficientType > &out) const =0
determine coefficients interpolating a given function
void interpolate(const F &ff, std::vector< C > &out) const
determine coefficients interpolating a given function
Definition: virtualinterface.hh:232
Return a proper base class for functions to use with LocalInterpolation.
Definition: virtualinterface.hh:45
FunctionBaseDummy FunctionBase
Definition: virtualinterface.hh:71
FunctionBaseDummy VirtualFunctionBase
Definition: virtualinterface.hh:70
FunctionBaseDummy type
Base class type for functions to use with LocalInterpolation.
Definition: virtualinterface.hh:78
Domain DomainType
Definition: virtualinterface.hh:62
virtual base class for a local basis
Definition: virtualinterface.hh:95
virtual unsigned int order() const =0
Polynomial order of the 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 unsigned int size() const =0
Number of shape functions.
T Traits
Definition: virtualinterface.hh:97
virtual void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const =0
Evaluate all basis function at given position.
virtual ~LocalBasisVirtualInterface()
Definition: virtualinterface.hh:100
virtual void partial(const std::array< unsigned int, Traits::dimDomain > &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 base class for a local interpolation
Definition: virtualinterface.hh:157
std::function< RangeType(DomainType)> FunctionType
type of function to interpolate
Definition: virtualinterface.hh:161
RangeType::field_type CoefficientType
type of the coefficient vector in the interpolate method
Definition: virtualinterface.hh:164
virtual ~LocalInterpolationVirtualInterfaceBase()
Definition: virtualinterface.hh:166
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:257
virtual ~LocalCoefficientsVirtualInterface()
Definition: virtualinterface.hh:260
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:284
virtual const Traits::LocalBasisType & localBasis() const =0
virtual unsigned int size() const =0
virtual const GeometryType type() const =0
virtual const Traits::LocalCoefficientsType & localCoefficients() const =0
virtual LocalFiniteElementVirtualInterface< T > * clone() const =0
LocalFiniteElementTraits< LocalBasisVirtualInterface< LocalBasisTraits >, LocalCoefficientsVirtualInterface, LocalInterpolationVirtualInterface< typename LocalBasisTraits::DomainType, typename LocalBasisTraits::RangeType > > Traits
Definition: virtualinterface.hh:292
virtual ~LocalFiniteElementVirtualInterface()
Definition: virtualinterface.hh:294
virtual const Traits::LocalInterpolationType & localInterpolation() const =0