dune-fem  2.8-git
space/shapefunctionset/orthonormal.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_SPACE_SHAPEFUNCTIONSET_ORTHONORMAL_HH
2 #define DUNE_FEM_SPACE_SHAPEFUNCTIONSET_ORTHONORMAL_HH
3 
4 #include <cassert>
5 #include <cstddef>
6 
7 #include <dune/geometry/type.hh>
8 
11 
16 
17 namespace Dune
18 {
19 
20  namespace Fem
21  {
22 
23 #ifndef DOXYGEN
24 
25  // OrthonormalShapeFunctions
26  // -------------------------
27 
28  template< int dimension >
29  struct OrthonormalShapeFunctions;
30 
31  template<>
32  struct OrthonormalShapeFunctions< 1 >
33  {
34  static constexpr std::size_t size ( int order )
35  {
36  return static_cast< std::size_t >( order+1 );
37  }
38  };
39 
40  template<>
41  struct OrthonormalShapeFunctions< 2 >
42  {
43  static constexpr std::size_t size ( int order )
44  {
45  return static_cast< std::size_t >( (order+2)*(order+1)/2 );
46  }
47  };
48 
49  template<>
50  struct OrthonormalShapeFunctions< 3 >
51  {
52  static constexpr std::size_t size ( int order )
53  {
54  return static_cast< std::size_t >( ((order+1)*(order+2)*(2*order+3)/6+(order+1)*(order+2)/2)/2 );
55  }
56  };
57 
58 #endif // #ifndef DOXYGEN
59 
60 
61 
62  // OrthonormalShapeFunctionSet
63  // ---------------------------
64 
65  template< class FunctionSpace >
67  {
69 
70  static_assert( FunctionSpace::dimDomain <= 3, "OrthonormalShapeFunctionSet only implemented up to domain dimension 3" );
71  static_assert( FunctionSpace::dimRange == 1, "OrthonormalShapeFunctionSet only implemented for scalar function spaces" );
72 
73  public:
76 
77  private:
78  typedef typename FunctionSpaceType::DomainFieldType DomainFieldType;
79  typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
80 
81  typedef RangeFieldType (*Evaluate) ( const int, const DomainFieldType * );
82  typedef void (*Jacobian) ( const int i, const DomainFieldType *, RangeFieldType * );
83 
84  // hessian is only available for 2d-simplices (aka triangles)
85  typedef RangeFieldType Array[ 3 ];
86  typedef void (*Hessian) ( const int i, const DomainFieldType *, Array & );
87 
91 
92  static void setFunctionPointers(const Dune::GeometryType& geomType,
93  Evaluate &evaluate, Jacobian &jacobian )
94  {
95  if( geomType.isLine() )
96  {
99  return ;
100  }
101  else if( geomType.isQuadrilateral() )
102  {
105  return ;
106  }
107  else if( geomType.isTriangle() )
108  {
111  return ;
112  }
113  else if( geomType.isHexahedron() )
114  {
117  return ;
118  }
119  else if ( geomType.isTetrahedron() )
120  {
123  return ;
124  }
125  else if( geomType.isPrism() )
126  {
129  return ;
130  }
131  else if ( geomType.isPyramid() )
132  {
135  return ;
136  }
137 
138  DUNE_THROW(InvalidStateException,"Invalid geometry type " << geomType );
139  }
140 
141  public:
144 
153 
159 
160  OrthonormalShapeFunctionSet ( GeometryType type, int order )
161  : order_( order ),
162  evaluate_( nullptr ),
163  jacobian_( nullptr ),
164  hessian_( nullptr )
165  {
166  // for type none simply create cube basis function set.
167  if( type.isNone() )
168  type = Dune::GeometryTypes::cube(type.dim());
169 
170  // set functions pointers for evaluate and jacobian
171  // depending on geometry type
172  setFunctionPointers( type, evaluate_, jacobian_ );
173  assert( evaluate_ );
174  assert( jacobian_ );
175 
176  if( dimension == 2 )
177  {
178  if( type.isTriangle() )
180  else if( type.isQuadrilateral() )
182  }
183  }
184 
185  OrthonormalShapeFunctionSet ( const ThisType & ) = default;
186 
188 
196 
198 
206  int order () const { return order_; }
207 
209  std::size_t constexpr size () const { return size( order() ); }
210 
212  static std::size_t constexpr size ( int order )
213  {
214  return OrthonormalShapeFunctions< dimension >::size( order );
215  }
216 
218  template< class Point, class Functor >
219  void evaluateEach ( const Point &x, Functor functor ) const
220  {
221  const DomainType y = Dune::Fem::coordinate( x );
222  RangeType value;
223  const std::size_t size = this->size();
224  for( std::size_t i = 0; i < size; ++i )
225  {
226  evaluate( i, y, value );
227  functor( i, value );
228  }
229  }
230 
232  template< class Point, class Functor >
233  void jacobianEach ( const Point &x, Functor functor ) const
234  {
235  const DomainType y = Dune::Fem::coordinate( x );
236  JacobianRangeType jacobian;
237  const std::size_t size = this->size();
238  for( std::size_t i = 0; i < size; ++i )
239  {
240  this->jacobian( i, y, jacobian );
241  functor( i, jacobian );
242  }
243  }
244 
246  template< class Point, class Functor >
247  void hessianEach ( const Point &x, Functor functor ) const
248  {
249  const DomainType y = Dune::Fem::coordinate( x );
250  HessianRangeType hessian;
251  const std::size_t size = this->size();
252  for( std::size_t i = 0; i < size; ++i )
253  {
254  this->hessian( i, y, hessian );
255  functor( i, hessian );
256  }
257  }
258 
261  private:
262  void evaluate ( std::size_t i, const DomainType &x, RangeType &value ) const
263  {
264  value[ 0 ] = evaluate_( i, &x[ 0 ] );
265  }
266 
267  void jacobian ( std::size_t i, const DomainType &x, JacobianRangeType &jacobian ) const
268  {
269  jacobian_( i, &x[ 0 ], &jacobian[ 0 ][ 0 ] );
270  }
271 
272  void hessian ( std::size_t i, const DomainType &x, HessianRangeType &hessian ) const
273  {
274  assert( hessian_ );
275 
276  RangeFieldType values[] = { 0, 0, 0 };
277  hessian_( i , &x[ 0 ], values );
278  for( unsigned int j = 0; j < FunctionSpaceType::dimDomain; ++j )
279  for( unsigned int k = 0; k < FunctionSpaceType::dimDomain; ++k )
280  hessian[ 0 ][ j ][ k ] = values[ j + k ];
281  }
282 
283  int order_;
284  Evaluate evaluate_;
285  Jacobian jacobian_;
286  Hessian hessian_;
287  };
288 
289 
290  } // namespace Fem
291 
292 } // namespace Dune
293 
294 #endif // #ifndef DUNE_FEM_SPACE_SHAPEFUNCTIONSET_ORTHONORMAL_HH
Definition: bindguard.hh:11
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:14
Definition: explicitfieldvector.hh:75
A vector valued function space.
Definition: functionspace.hh:60
FunctionSpaceTraits::DomainFieldType DomainFieldType
Intrinsic type used for values in the domain field (usually a double)
Definition: functionspaceinterface.hh:60
FunctionSpaceTraits::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:71
@ dimDomain
dimension of domain vector space
Definition: functionspaceinterface.hh:46
@ dimRange
dimension of range vector space
Definition: functionspaceinterface.hh:48
FunctionSpaceTraits::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:75
FunctionSpaceTraits::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:67
FunctionSpaceTraits::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:63
Definition: orthonormalbase_1d.hh:15
static void grad_line(const int i, DomainType xi, JacobianRangeType grad)
Definition: orthonormalbase_1d.hh:121
static RangeField eval_line(const int i, DomainType xi)
Definition: orthonormalbase_1d.hh:25
Definition: orthonormalbase_2d.hh:17
static void grad_quadrilateral_2d(const int i, DomainType xi, JacobianRangeType grad)
Definition: orthonormalbase_2d.hh:3084
static RangeField eval_quadrilateral_2d(const int i, DomainType xi)
Definition: orthonormalbase_2d.hh:2511
static void grad_triangle_2d(const int i, DomainType xi, JacobianRangeType grad)
Definition: orthonormalbase_2d.hh:927
static void hess_quadrilateral_2d(const int i, DomainType xi, HessianRangeType &h)
Definition: orthonormalbase_2d.hh:3882
static void hess_triangle_2d(const int i, DomainType xi, HessianRangeType &h)
Definition: orthonormalbase_2d.hh:2256
static RangeField eval_triangle_2d(const int i, DomainType xi)
Definition: orthonormalbase_2d.hh:25
Definition: orthonormalbase_3d.hh:17
static RangeField eval_tetrahedron_3d(const int i, DomainType xi)
Definition: orthonormalbase_3d.hh:24
static void grad_tetrahedron_3d(const int i, DomainType xi, JacobianRangeType grad)
Definition: orthonormalbase_3d.hh:9298
static RangeField eval_prism_3d(const int i, DomainType xi)
Definition: orthonormalbase_3d.hh:53195
static RangeField eval_pyramid_3d(const int i, DomainType xi)
Definition: orthonormalbase_3d.hh:27863
static void grad_pyramid_3d(const int i, DomainType xi, JacobianRangeType grad)
Definition: orthonormalbase_3d.hh:36345
static void grad_hexahedron_3d(const int i, DomainType xi, JacobianRangeType grad)
Definition: orthonormalbase_3d.hh:70486
static void grad_prism_3d(const int i, DomainType xi, JacobianRangeType grad)
Definition: orthonormalbase_3d.hh:58113
static RangeField eval_hexahedron_3d(const int i, DomainType xi)
Definition: orthonormalbase_3d.hh:67284
Definition: space/shapefunctionset/orthonormal.hh:67
FunctionSpaceType::RangeType RangeType
range type
Definition: space/shapefunctionset/orthonormal.hh:148
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: space/shapefunctionset/orthonormal.hh:152
OrthonormalShapeFunctionSet(ThisType &&)=default
OrthonormalShapeFunctionSet(GeometryType type, int order)
Definition: space/shapefunctionset/orthonormal.hh:160
int order() const
return order of shape functions
Definition: space/shapefunctionset/orthonormal.hh:206
static constexpr std::size_t size(int order)
please doc me
Definition: space/shapefunctionset/orthonormal.hh:212
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: space/shapefunctionset/orthonormal.hh:150
FunctionSpace FunctionSpaceType
function space type
Definition: space/shapefunctionset/orthonormal.hh:70
FunctionSpaceType::DomainType DomainType
Definition: space/shapefunctionset/orthonormal.hh:146
OrthonormalShapeFunctionSet & operator=(const ThisType &)=default
static const int dimension
Definition: space/shapefunctionset/orthonormal.hh:143
void evaluateEach(const Point &x, Functor functor) const
evalute each shape function
Definition: space/shapefunctionset/orthonormal.hh:219
constexpr std::size_t size() const
return number of shape functions
Definition: space/shapefunctionset/orthonormal.hh:209
OrthonormalShapeFunctionSet(const ThisType &)=default
void jacobianEach(const Point &x, Functor functor) const
evalute jacobian of each shape function
Definition: space/shapefunctionset/orthonormal.hh:233
void hessianEach(const Point &x, Functor functor) const
evalute hessian of each shape function
Definition: space/shapefunctionset/orthonormal.hh:247