dune-fem  2.4.1-rc
converter.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_FUNCTION_LOCALFUNCTION_CONVERTER_HH
2 #define DUNE_FEM_FUNCTION_LOCALFUNCTION_CONVERTER_HH
3 
4 #include <functional>
5 #include <type_traits>
6 #include <utility>
7 
11 
12 
13 namespace Dune
14 {
15 
16  namespace Fem
17  {
18 
19  // LocalFunctionConverter
20  // --------------------
21 
71  template< class HostLocalFunction, class Converter, template< class > class Storage = __InstationaryFunction::HoldCopy >
73  : private Storage< HostLocalFunction >
74  {
76  typedef Storage< HostLocalFunction > BaseType;
77 
78  typedef typename HostLocalFunction::RangeType HostRangeType;
79  typedef typename HostLocalFunction::JacobianRangeType HostJacobianRangeType;
80  typedef typename HostLocalFunction::HessianRangeType HostHessianRangeType;
81 
82  public:
83  // obtain new dimRange from Converter
84  static const int dimRange = decltype( std::declval< Converter >() ( std::declval< HostRangeType >() ) ) ::dimension;
85 
86  // define new FunctionSpace
88 
89  // types from HostLocalFunction
90  typedef typename HostLocalFunction::EntityType EntityType;
91 
92  // types from FunctionSpace
93  typedef typename FunctionSpaceType::DomainType DomainType;
94  typedef typename FunctionSpaceType::RangeType RangeType;
95  typedef typename FunctionSpaceType::JacobianRangeType JacobianRangeType;
96  typedef typename FunctionSpaceType::HessianRangeType HessianRangeType;
97  typedef typename FunctionSpaceType::DomainFieldType DomainFieldType;
98  typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
99 
100  static const int dimDomain = FunctionSpaceType::dimDomain;
101 
102  LocalFunctionConverter ( const HostLocalFunction &hostLocalFunction, const Converter &converter = Converter() )
103  : BaseType( hostLocalFunction ), converter_( converter )
104  {}
105 
106  LocalFunctionConverter ( HostLocalFunction &&hostLocalFunction, const Converter &converter = Converter() )
107  : BaseType( std::move( hostLocalFunction ) ), converter_( converter )
108  {}
109 
110  template< class Point >
111  void evaluate ( const Point &p, RangeType &ret ) const
112  {
113  HostRangeType hRet;
114  this->get().evaluate( p, hRet );
115  ret = converter_( hRet );
116  }
117 
118  template< class Point >
119  void jacobian ( const Point &p, JacobianRangeType &jac ) const
120  {
121  HostJacobianRangeType hJac;
122  this->get().jacobian( p, hJac );
123  jac = converter_( hJac );
124  }
125 
126  template< class Point >
127  void hessian ( const Point &p, HessianRangeType &hes ) const
128  {
129  HostHessianRangeType hHes;
130  this->get().hessian( p, hHes );
131  hes = converter_( hHes );
132  }
133 
134  template< class Quadrature, class Vector >
135  void evaluateQuadrature ( const Quadrature &quad, Vector &vector ) const
136  {
137  evaluateQuadratureImp( quad, vector, vector[ 0 ] );
138  }
139 
140  int order () const { return this->get().order(); }
141 
142  const EntityType &entity () const { return this->get().entity(); }
143 
144  void init ( const EntityType &entity ) { this->get().init( entity ); }
145 
146  protected:
147  template< class QuadratureType, class VectorType >
148  void evaluateQuadratureImp ( const QuadratureType &quadrature, VectorType &values, const RangeType & ) const
149  {
150  const unsigned int nop = quadrature.nop();
151  for( unsigned int qp = 0; qp < nop; ++qp )
152  evaluate( quadrature[ qp ], values[ qp ] );
153  }
154 
155  template< class QuadratureType, class VectorType >
156  void evaluateQuadratureImp ( const QuadratureType &quadrature, VectorType &values, const JacobianRangeType & ) const
157  {
158  const unsigned int nop = quadrature.nop();
159  for( unsigned int qp = 0; qp < nop; ++qp )
160  jacobian( quadrature[ qp ], values[ qp ] );
161  }
162 
163  Converter converter_;
164  };
165 
166 
167 
168  // localFunctionConverter
169  // ----------------------
170 
171  template< class HostLocalFunction, class Converter >
173  localFunctionConverter ( HostLocalFunction hostLocalFunction, const Converter &converter = Converter() )
174  {
176  return LocalFunctionConverterType( std::move( hostLocalFunction ), converter );
177  }
178 
179  template< class HostLocalFunction, class Converter >
180  LocalFunctionConverter< typename std::remove_const< HostLocalFunction >::type, Converter, __InstationaryFunction::HoldReference >
181  localFunctionConverter ( std::reference_wrapper< HostLocalFunction > hostLocalFunction, const Converter &converter = Converter() )
182  {
183  typedef LocalFunctionConverter< typename std::remove_const< HostLocalFunction >::type, Converter, __InstationaryFunction::HoldReference > LocalFunctionConverterType;
184  return LocalFunctionConverterType( hostLocalFunction.get(), converter );
185  }
186 
187  } // namespace Fem
188 
189 } // namespace Dune
190 
191 #endif // #ifndef DUNE_FEM_FUNCTION_LOCALFUNCTION_CONVERTER_HH
LocalFunctionConverter(const HostLocalFunction &hostLocalFunction, const Converter &converter=Converter())
Definition: converter.hh:102
FunctionSpaceType::RangeFieldType RangeFieldType
Definition: converter.hh:98
void jacobian(const Point &p, JacobianRangeType &jac) const
Definition: converter.hh:119
FunctionSpaceType::RangeType RangeType
Definition: converter.hh:94
LocalFunctionConverter(HostLocalFunction &&hostLocalFunction, const Converter &converter=Converter())
Definition: converter.hh:106
static const int dimDomain
Definition: converter.hh:100
HostLocalFunction::EntityType EntityType
Definition: converter.hh:90
void evaluate(const Point &p, RangeType &ret) const
Definition: converter.hh:111
Converter converter_
Definition: converter.hh:163
void hessian(const Point &p, HessianRangeType &hes) const
Definition: converter.hh:127
ToNewDimRangeFunctionSpace< typename HostLocalFunction::FunctionSpaceType, dimRange >::Type FunctionSpaceType
Definition: converter.hh:87
static const int dimRange
Definition: converter.hh:84
int order() const
Definition: converter.hh:140
Definition: coordinate.hh:4
LocalFunctionConverter< HostLocalFunction, Converter, __InstationaryFunction::HoldCopy > localFunctionConverter(HostLocalFunction hostLocalFunction, const Converter &converter=Converter())
Definition: converter.hh:173
void evaluateQuadrature(const Quadrature &quad, Vector &vector) const
Definition: converter.hh:135
STL namespace.
FunctionSpaceType::JacobianRangeType JacobianRangeType
Definition: converter.hh:95
FunctionSpaceType::DomainFieldType DomainFieldType
Definition: converter.hh:97
void evaluateQuadratureImp(const QuadratureType &quadrature, VectorType &values, const JacobianRangeType &) const
Definition: converter.hh:156
implementation of a Dune::Fem::LocalFunction on a FunctionSpace V restircted/prolongated from an othe...
Definition: converter.hh:72
void move(ArrayInterface< T > &array, const unsigned int oldOffset, const unsigned int newOffset, const unsigned int length)
Definition: array_inline.hh:38
void evaluateQuadratureImp(const QuadratureType &quadrature, VectorType &values, const RangeType &) const
Definition: converter.hh:148
FunctionSpaceType::DomainType DomainType
Definition: converter.hh:93
const EntityType & entity() const
Definition: converter.hh:142
void init(const EntityType &entity)
Definition: converter.hh:144
FunctionSpaceType::HessianRangeType HessianRangeType
Definition: converter.hh:96
convert functions space to space with new dim range
Definition: functionspace.hh:246
actual interface class for quadratures
Definition: quadrature.hh:320