dune-fem  2.8-git
space/basisfunctionset/default.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH
2 #define DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH
3 
4 // C++ includes
5 #include <cassert>
6 #include <cstddef>
7 #include <memory>
8 #include <type_traits>
9 #include <utility>
10 
11 // dune-geometry includes
12 #include <dune/geometry/referenceelements.hh>
13 #include <dune/geometry/type.hh>
14 
16 
17 // dune-fem includes
23 #include <dune/fem/version.hh>
24 
27 
28 namespace Dune
29 {
30 
31  namespace Fem
32  {
33 
34  // DefaultBasisFunctionSet
35  // -----------------------
36 
50  template< class Entity, class ShapeFunctionSet >
52  {
54 
55  public:
57  typedef Entity EntityType;
60 
61  // if underlying shape function set was created with storage CodegenStorage
62  // then this value should be true (see selectcaching.hh)
63  static constexpr bool codegenShapeFunctionSet = detail::IsCodegenShapeFunctionSet< ShapeFunctionSetType >::value;
64 
65  protected:
69 
71 
72  typedef typename EntityType::Geometry Geometry ;
73 
74  typedef typename Geometry::ctype ctype;
75  public:
76  // slight misuse of struct ToLocalFunctionSpace!!!
79 
81  typedef typename FunctionSpaceType::RangeType RangeType;
83  typedef typename FunctionSpaceType::DomainType DomainType;
85  typedef typename FunctionSpaceType::JacobianRangeType JacobianRangeType;
87  typedef typename FunctionSpaceType::HessianRangeType HessianRangeType;
88 
89  typedef typename FunctionSpaceType::ScalarFunctionSpaceType ScalarFunctionSpaceType;
90 
91  typedef typename ScalarFunctionSpaceType::RangeType ScalarRangeType;
92  typedef typename ScalarFunctionSpaceType::JacobianRangeType ScalarJacobianRangeType;
93 
95  typedef std::decay_t< decltype( Dune::ReferenceElements< ctype, Geometry::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) ) > ReferenceElementType;
96 
97  static const int pointSetId = detail::SelectPointSetId< ShapeFunctionSetType >::value;
98 
101 
104  : entity_( &entity ),
106  {
107  // Note that this should be geometry_ = entity.geometry()
108  // But Dune::Geometries are not assignable ...
109  geometry_.reset();
110  geometry_.emplace( entity.geometry() );
111  }
112 
114  : entity_( other.entity_ ),
116  {
117  // Note that this should be geometry_ = entity.geometry()
118  // But Dune::Geometries are not assignable ...
119  geometry_.reset();
120  if( other.geometry_ )
121  geometry_.emplace( other.geometry_.value() );
122  }
123 
125  {
126  entity_ = other.entity_;
128 
129  // Note that this should be geometry_ = entity.geometry()
130  // But Dune::Geometries are not assignable ...
131  geometry_.reset();
132  if( other.geometry_ )
133  geometry_.emplace( other.geometry_.value() );
134  return *this;
135  }
136 
137  // Basis Function Set Interface Methods
138  // ------------------------------------
139 
141  int order () const { return shapeFunctionSet().order(); }
142 
144  std::size_t size () const { return shapeFunctionSet().size(); }
145 
147  auto referenceElement () const
148  -> decltype( Dune::ReferenceElements< ctype, Geometry::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) )
149  {
150  return Dune::ReferenceElements< ctype, Geometry::coorddimension >::general( type() );
151  }
152 
156  template< class QuadratureType, class Vector, class DofVector >
157  void axpy ( const QuadratureType &quad, const Vector &values, DofVector &dofs ) const
158  {
159  axpyImpl( quad, values, dofs, values[ 0 ] );
160  }
161 
168  template< class QuadratureType, class VectorA, class VectorB, class DofVector >
169  void axpy ( const QuadratureType &quad, const VectorA &valuesA, const VectorB &valuesB, DofVector &dofs ) const
170  {
171  assert( valuesA.size() > 0 );
172  assert( valuesB.size() > 0 );
173 
174  axpyImpl( quad, valuesA, dofs, valuesA[ 0 ] );
175  axpyImpl( quad, valuesB, dofs, valuesB[ 0 ] );
176  }
177 
181  template< class Point, class DofVector >
182  void axpy ( const Point &x, const RangeType &valueFactor, DofVector &dofs ) const
183  {
185  shapeFunctionSet().evaluateEach( x, f );
186  }
187 
191  template< class Point, class DofVector >
192  void axpy ( const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs ) const
193  {
194  typedef typename Geometry::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
195  const Geometry &geo = geometry();
196  const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed( coordinate( x ) );
197  LocalJacobianRangeType tmpJacobianFactor;
198  for( int r = 0; r < FunctionSpaceType::dimRange; ++r )
199  gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
200 
202  shapeFunctionSet().jacobianEach( x, f );
203  }
204 
208  template< class Point, class DofVector >
209  void axpy ( const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor,
210  DofVector &dofs ) const
211  {
212  axpy( x, valueFactor, dofs );
213  axpy( x, jacobianFactor, dofs );
214  }
215 
218  template< class Point, class DofVector >
219  void axpy ( const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs ) const
220  {
221  typedef typename Geometry::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
222  const Geometry &geo = geometry();
223  const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed( coordinate( x ) );
224  LocalHessianRangeType tmpHessianFactor( RangeFieldType(0) );
225  // don't know how to work directly with the DiagonalMatrix
226  // returned from YaspGrid since gjit[k][l] is not correctly
227  // implemented for k!=l so convert first into a dense matrix...
228  Dune::FieldMatrix<double,DomainType::dimension,DomainType::dimension>
229  G = gjit;
230  DomainType Hg;
231  for( int r = 0; r < FunctionSpaceType::dimRange; ++r )
232  for( int j = 0; j < gjit.cols; ++j )
233  for( int k = 0; k < gjit.cols; ++k )
234  {
235  hessianFactor[r].mv(G[k],Hg);
236  tmpHessianFactor[r][j][k] += Hg * G[j];
237  }
239  shapeFunctionSet().hessianEach( x, f );
240  }
241 
242 
244  template< class QuadratureType, class DofVector, class RangeArray >
245  void evaluateAll ( const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges ) const
246  {
247  assert( ranges.size() >= quad.nop() );
248 
249  // if shape function set supports codegen and quadrature supports caching
250  if constexpr ( codegenShapeFunctionSet && std::is_base_of< CachingInterface, QuadratureType > :: value)
251  {
252  typedef Codegen :: EvaluateCallerInterfaceTraits< QuadratureType, RangeArray, DofVector > Traits;
253  typedef Codegen :: EvaluateCallerInterface< Traits > BaseEvaluationType;
254 
255  // get base function evaluate caller
256  const auto& baseEval = BaseEvaluationType::storage( *this, rangeCache( quad ), quad );
257 
258  // true if implementation exists, otherwise this is a nullptr
259  if( baseEval )
260  {
261  baseEval->evaluateRanges( quad, dofs, ranges );
262  return ;
263  }
264  }
265 
266  {
267  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
268  const unsigned int nop = quad.nop();
269  for( unsigned int qp = 0; qp < nop; ++qp )
270  {
271  evaluateAll( quad[ qp ], dofs, ranges[ qp ] );
272  }
273  }
274  }
275 
277  template< class Point, class DofVector >
278  void evaluateAll ( const Point &x, const DofVector &dofs, RangeType &value ) const
279  {
280  value = RangeType( 0 );
282  shapeFunctionSet().evaluateEach( x, f );
283  }
284 
286  template< class Point, class RangeArray >
287  void evaluateAll ( const Point &x, RangeArray &values ) const
288  {
289  assert( values.size() >= size() );
290  std::fill( values.begin(), values.end(), RangeType(0) );
291  AssignFunctor< RangeArray > f( values );
292  shapeFunctionSet().evaluateEach( x, f );
293  }
294 
296  template< class QuadratureType, class DofVector, class JacobianArray >
297  void jacobianAll ( const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians ) const
298  {
299  assert( jacobians.size() >= quad.nop() );
300 
301  // if shape function set supports codegen and quadrature supports caching
302  if constexpr ( codegenShapeFunctionSet && std::is_base_of< CachingInterface, QuadratureType > :: value)
303  {
304  typedef Codegen :: EvaluateCallerInterfaceTraits< QuadratureType, JacobianArray, DofVector, Geometry > Traits;
305  typedef Codegen :: EvaluateCallerInterface< Traits > BaseEvaluationType;
306 
307  // get base function evaluate caller (calls axpyRanges)
308  const auto& baseEval = BaseEvaluationType::storage( *this, jacobianCache( quad ), quad );
309 
310  // true if implementation exists
311  if( baseEval )
312  {
313  // call appropriate axpyJacobian method
314  const Geometry &geo = geometry();
315  baseEval->evaluateJacobians( quad, geo, dofs, jacobians );
316  return ;
317  }
318  }
319 
320  {
321  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
322  const unsigned int nop = quad.nop();
323  for( unsigned int qp = 0; qp < nop; ++qp )
324  {
325  jacobianAll( quad[ qp ], dofs, jacobians[ qp ] );
326  }
327  }
328  }
329 
331  template< class Point, class DofVector >
332  void jacobianAll ( const Point &x, const DofVector &dofs, JacobianRangeType &jacobian ) const
333  {
334  LocalJacobianRangeType localJacobian( RangeFieldType( 0 ) );
336  shapeFunctionSet().jacobianEach( x, f );
337 
338  typedef JacobianTransformation< Geometry > Transformation;
339  const Geometry &geo = geometry();
340  Transformation transformation( geo, coordinate( x ) );
341  transformation( localJacobian, jacobian );
342  }
343 
345  template< class Point, class JacobianRangeArray >
346  void jacobianAll ( const Point &x, JacobianRangeArray &jacobians ) const
347  {
348  assert( jacobians.size() >= size() );
349  typedef JacobianTransformation< Geometry > Transformation;
350  const Geometry &geo = geometry();
351 
352  Transformation transformation( geo, coordinate( x ) );
353  AssignFunctor< JacobianRangeArray, Transformation > f( jacobians, transformation );
354  shapeFunctionSet().jacobianEach( x, f );
355  }
356 
358  template< class QuadratureType, class DofVector, class HessianArray >
359  void hessianAll ( const QuadratureType &quad, const DofVector &dofs, HessianArray &hessians ) const
360  {
361  assert( hessians.size() >= quad.nop() );
362  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
363  const unsigned int nop = quad.nop();
364  for( unsigned int qp = 0; qp < nop; ++qp )
365  {
366  hessianAll( quad[ qp ], dofs, hessians[ qp ] );
367  }
368  }
369 
371  template< class Point, class DofVector >
372  void hessianAll ( const Point &x, const DofVector &dofs, HessianRangeType &hessian ) const
373  {
374  LocalHessianRangeType localHessian( typename LocalHessianRangeType::value_type( RangeFieldType( 0 ) ) );
376  shapeFunctionSet().hessianEach( x, f );
377 
378  typedef HessianTransformation< Geometry > Transformation;
379  const Geometry &geo = geometry();
380  Transformation transformation( geo, coordinate( x ) );
381  transformation( localHessian, hessian );
382  }
383 
385  template< class Point, class HessianRangeArray >
386  void hessianAll ( const Point &x, HessianRangeArray &hessians ) const
387  {
388  assert( hessians.size() >= size() );
389  typedef HessianTransformation< Geometry > Transformation;
390  const Geometry &geo = geometry();
391  Transformation transformation( geo, coordinate( x ) );
392  AssignFunctor< HessianRangeArray, Transformation > f( hessians, transformation );
393  shapeFunctionSet().hessianEach( x, f );
394  }
395 
397  const Entity &entity () const
398  {
399  assert( valid() );
400  return *entity_;
401  }
402 
404  bool valid () const { return bool(entity_); }
405 
407  Dune::GeometryType type () const { return entity().type(); }
408 
409  // Non-interface methods
410  // ---------------------
411 
414 
415  protected:
417  template< class QuadratureType, class RangeArray, class DofVector >
418  void axpyImpl ( const QuadratureType &quad, const RangeArray &rangeFactors, DofVector &dofs, const RangeType& ) const
419  {
420  assert( rangeFactors.size() >= quad.nop() );
421 
422  // if shape function set supports codegen and quadrature supports caching
423  if constexpr ( codegenShapeFunctionSet && std::is_base_of< CachingInterface, QuadratureType > :: value)
424  {
425  typedef Codegen :: EvaluateCallerInterfaceTraits< QuadratureType, RangeArray, DofVector > Traits;
426  typedef Codegen :: EvaluateCallerInterface< Traits > BaseEvaluationType;
427 
428  // get base function evaluate caller (calls axpyRanges)
429  const auto& baseEval = BaseEvaluationType::storage( *this, rangeCache( quad ), quad );
430 
431  // true if implementation exists
432  if( baseEval )
433  {
434  baseEval->axpyRanges( quad, rangeFactors, dofs );
435  return ;
436  }
437  }
438 
439  {
440  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
441  const unsigned int nop = quad.nop();
442  for( unsigned int qp = 0; qp < nop; ++qp )
443  {
444  axpy( quad[ qp ], rangeFactors[ qp ], dofs );
445  }
446  }
447  }
448 
450  template< class QuadratureType, class JacobianArray, class DofVector >
451  void axpyImpl ( const QuadratureType &quad, const JacobianArray &jacobianFactors, DofVector &dofs, const JacobianRangeType& ) const
452  {
453  assert( jacobianFactors.size() >= quad.nop() );
454  // if shape function set supports codegen and quadrature supports caching
455  if constexpr ( codegenShapeFunctionSet && std::is_base_of< CachingInterface, QuadratureType > :: value)
456  {
457  typedef Codegen :: EvaluateCallerInterfaceTraits< QuadratureType, JacobianArray, DofVector, Geometry > Traits;
458  typedef Codegen :: EvaluateCallerInterface< Traits > BaseEvaluationType;
459 
460  // get base function evaluate caller (calls axpyRanges)
461  const auto& baseEval = BaseEvaluationType::storage( *this, jacobianCache( quad ), quad );
462 
463  // true if implementation exists
464  if( baseEval )
465  {
466  // call appropriate axpyRanges method
467  const Geometry &geo = geometry();
468  baseEval->axpyJacobians( quad, geo, jacobianFactors, dofs );
469  return ;
470  }
471  }
472 
473  {
474  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
475  const unsigned int nop = quad.nop();
476  for( unsigned int qp = 0; qp < nop; ++qp )
477  {
478  axpy( quad[ qp ], jacobianFactors[ qp ], dofs );
479  }
480  }
481  }
482 
483  template <class QuadratureType>
484  const auto& rangeCache( const QuadratureType& quad ) const
485  {
486  return shapeFunctionSet().scalarShapeFunctionSet().impl().rangeCache( quad );
487  }
488 
489  template <class QuadratureType>
490  const auto& jacobianCache( const QuadratureType& quad ) const
491  {
492  return shapeFunctionSet().scalarShapeFunctionSet().impl().jacobianCache( quad );
493  }
494 
495  protected:
496  Geometry geometry () const { return geometry_.value(); }
497 
498  const EntityType *entity_ = nullptr;
500  std::optional< Geometry > geometry_;
501  };
502 
503  } // namespace Fem
504 
505 } // namespace Dune
506 
507 #endif // #ifndef DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH
Definition: bindguard.hh:11
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:14
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
Definition: explicitfieldvector.hh:75
Definition: misc/functor.hh:31
Definition: space/basisfunctionset/default.hh:52
EntityType::Geometry Geometry
Definition: space/basisfunctionset/default.hh:72
const auto & rangeCache(const QuadratureType &quad) const
Definition: space/basisfunctionset/default.hh:484
FunctionSpaceType::DomainType DomainType
domain type
Definition: space/basisfunctionset/default.hh:83
void evaluateAll(const Point &x, const DofVector &dofs, RangeType &value) const
Definition: space/basisfunctionset/default.hh:278
LocalFunctionSpaceType::HessianRangeType LocalHessianRangeType
Definition: space/basisfunctionset/default.hh:68
DefaultBasisFunctionSet()
constructor
Definition: space/basisfunctionset/default.hh:100
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: space/basisfunctionset/default.hh:87
const EntityType * entity_
Definition: space/basisfunctionset/default.hh:498
void axpy(const QuadratureType &quad, const VectorA &valuesA, const VectorB &valuesB, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:169
auto referenceElement() const -> decltype(Dune::ReferenceElements< ctype, Geometry::coorddimension >::general(std::declval< const Dune::GeometryType & >()))
return reference element
Definition: space/basisfunctionset/default.hh:147
const Entity & entity() const
return entity
Definition: space/basisfunctionset/default.hh:397
ToNewDimDomainFunctionSpace< LocalFunctionSpaceType, Geometry::coorddimension >::Type FunctionSpaceType
type of function space
Definition: space/basisfunctionset/default.hh:78
void axpy(const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor, DofVector &dofs) const
evaluate all basis function and derivatives and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:209
static constexpr bool codegenShapeFunctionSet
Definition: space/basisfunctionset/default.hh:63
std::optional< Geometry > geometry_
Definition: space/basisfunctionset/default.hh:500
FunctionSpaceType::RangeType RangeType
range type
Definition: space/basisfunctionset/default.hh:81
ShapeFunctionSetType::FunctionSpaceType LocalFunctionSpaceType
Definition: space/basisfunctionset/default.hh:66
LocalFunctionSpaceType::RangeFieldType RangeFieldType
Definition: space/basisfunctionset/default.hh:70
const auto & jacobianCache(const QuadratureType &quad) const
Definition: space/basisfunctionset/default.hh:490
LocalFunctionSpaceType::JacobianRangeType LocalJacobianRangeType
Definition: space/basisfunctionset/default.hh:67
void axpy(const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs) const
evaluate all derivatives of all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:192
void jacobianAll(const Point &x, JacobianRangeArray &jacobians) const
Definition: space/basisfunctionset/default.hh:346
bool valid() const
return true if entity pointer is set
Definition: space/basisfunctionset/default.hh:404
DefaultBasisFunctionSet(const EntityType &entity, const ShapeFunctionSet &shapeFunctionSet=ShapeFunctionSet())
constructor
Definition: space/basisfunctionset/default.hh:103
const ShapeFunctionSetType & shapeFunctionSet() const
return shape function set
Definition: space/basisfunctionset/default.hh:413
void axpy(const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs) const
Add H:D^2phi to each dof.
Definition: space/basisfunctionset/default.hh:219
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: space/basisfunctionset/default.hh:85
void hessianAll(const Point &x, HessianRangeArray &hessians) const
Definition: space/basisfunctionset/default.hh:386
void evaluateAll(const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges) const
Definition: space/basisfunctionset/default.hh:245
ShapeFunctionSetType shapeFunctionSet_
Definition: space/basisfunctionset/default.hh:499
ScalarFunctionSpaceType::RangeType ScalarRangeType
Definition: space/basisfunctionset/default.hh:91
DefaultBasisFunctionSet & operator=(const DefaultBasisFunctionSet &other)
Definition: space/basisfunctionset/default.hh:124
void hessianAll(const Point &x, const DofVector &dofs, HessianRangeType &hessian) const
Definition: space/basisfunctionset/default.hh:372
void hessianAll(const QuadratureType &quad, const DofVector &dofs, HessianArray &hessians) const
Definition: space/basisfunctionset/default.hh:359
DefaultBasisFunctionSet(const DefaultBasisFunctionSet &other)
Definition: space/basisfunctionset/default.hh:113
void evaluateAll(const Point &x, RangeArray &values) const
Definition: space/basisfunctionset/default.hh:287
ShapeFunctionSet ShapeFunctionSetType
shape function set type
Definition: space/basisfunctionset/default.hh:59
ScalarFunctionSpaceType::JacobianRangeType ScalarJacobianRangeType
Definition: space/basisfunctionset/default.hh:92
static const int pointSetId
Definition: space/basisfunctionset/default.hh:97
std::size_t size() const
return size of basis function set
Definition: space/basisfunctionset/default.hh:144
void jacobianAll(const Point &x, const DofVector &dofs, JacobianRangeType &jacobian) const
Definition: space/basisfunctionset/default.hh:332
Geometry geometry() const
Definition: space/basisfunctionset/default.hh:496
std::decay_t< decltype(Dune::ReferenceElements< ctype, Geometry::coorddimension >::general(std::declval< const Dune::GeometryType & >))) > ReferenceElementType
type of reference element
Definition: space/basisfunctionset/default.hh:95
void axpyImpl(const QuadratureType &quad, const RangeArray &rangeFactors, DofVector &dofs, const RangeType &) const
evaluate all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:418
Entity EntityType
entity type
Definition: space/basisfunctionset/default.hh:57
void axpyImpl(const QuadratureType &quad, const JacobianArray &jacobianFactors, DofVector &dofs, const JacobianRangeType &) const
evaluate all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:451
Dune::GeometryType type() const
return geometry type
Definition: space/basisfunctionset/default.hh:407
void axpy(const Point &x, const RangeType &valueFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:182
Geometry::ctype ctype
Definition: space/basisfunctionset/default.hh:74
FunctionSpaceType::ScalarFunctionSpaceType ScalarFunctionSpaceType
Definition: space/basisfunctionset/default.hh:89
void jacobianAll(const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians) const
Definition: space/basisfunctionset/default.hh:297
void axpy(const QuadratureType &quad, const Vector &values, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:157
int order() const
return order of basis function set
Definition: space/basisfunctionset/default.hh:141
Definition: space/basisfunctionset/functor.hh:108
Definition: space/basisfunctionset/functor.hh:132
Definition: transformation.hh:36
Definition: transformation.hh:92
A vector valued function space.
Definition: functionspace.hh:60
convert functions space to space with new dim domain
Definition: functionspace.hh:246
FunctionSpaceTraits::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:75
FunctionSpaceTraits::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:63
Interface class for shape function sets.
Definition: shapefunctionset/shapefunctionset.hh:33
void hessianEach(const Point &x, Functor functor) const
evalute hessian of each shape function
void evaluateEach(const Point &x, Functor functor) const
evalute each shape function
std::size_t size() const
return number of shape functions
int order() const
return order of shape functions
void jacobianEach(const Point &x, Functor functor) const
evalute jacobian of each shape function