1#ifndef DUNE_FEM_SPACE_SHAPEFUNCTIONSET_ORTHONORMAL_HH
2#define DUNE_FEM_SPACE_SHAPEFUNCTIONSET_ORTHONORMAL_HH
7#include <dune/geometry/type.hh>
28 template<
int dimension >
29 struct OrthonormalShapeFunctions;
32 struct OrthonormalShapeFunctions< 1 >
34 static constexpr std::size_t size (
int order )
36 return static_cast< std::size_t
>( order+1 );
41 struct OrthonormalShapeFunctions< 2 >
43 static constexpr std::size_t size (
int order )
45 return static_cast< std::size_t
>( (order+2)*(order+1)/2 );
50 struct OrthonormalShapeFunctions< 3 >
52 static constexpr std::size_t size (
int order )
54 return static_cast< std::size_t
>( ((order+1)*(order+2)*(2*order+3)/6+(order+1)*(order+2)/2)/2 );
65 template<
class FunctionSpace >
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" );
81 typedef RangeFieldType (*Evaluate) (
const int,
const DomainFieldType * );
82 typedef void (*Jacobian) (
const int i,
const DomainFieldType *, RangeFieldType * );
85 typedef RangeFieldType Array[ 3 ];
86 typedef void (*Hessian) (
const int i,
const DomainFieldType *, Array & );
92 static void setFunctionPointers(
const Dune::GeometryType& geomType,
93 Evaluate &evaluate, Jacobian &jacobian )
95 if( geomType.isLine() )
101 else if( geomType.isQuadrilateral() )
107 else if( geomType.isTriangle() )
113 else if( geomType.isHexahedron() )
119 else if ( geomType.isTetrahedron() )
125 else if( geomType.isPrism() )
131 else if ( geomType.isPyramid() )
138 DUNE_THROW(InvalidStateException,
"Invalid geometry type " << geomType );
162 evaluate_( nullptr ),
163 jacobian_( nullptr ),
168 type = Dune::GeometryTypes::cube(type.dim());
172 setFunctionPointers( type, evaluate_, jacobian_ );
178 if( type.isTriangle() )
180 else if( type.isQuadrilateral() )
206 int order ()
const {
return order_; }
214 return OrthonormalShapeFunctions< dimension >::size(
order );
218 template<
class Po
int,
class Functor >
223 const std::size_t
size = this->
size();
224 for( std::size_t i = 0; i <
size; ++i )
226 evaluate( i, y, value );
232 template<
class Po
int,
class Functor >
237 const std::size_t
size = this->
size();
238 for( std::size_t i = 0; i <
size; ++i )
240 this->jacobian( i, y, jacobian );
241 functor( i, jacobian );
246 template<
class Po
int,
class Functor >
251 const std::size_t
size = this->
size();
252 for( std::size_t i = 0; i <
size; ++i )
254 this->hessian( i, y, hessian );
255 functor( i, hessian );
264 value[ 0 ] = evaluate_( i, &x[ 0 ] );
269 jacobian_( i, &x[ 0 ], &jacobian[ 0 ][ 0 ] );
276 RangeFieldType values[] = { 0, 0, 0 };
277 hessian_( i , &x[ 0 ], values );
280 hessian[ 0 ][ j ][ k ] = values[ j + k ];
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: space/shapefunctionset/orthonormal.hh:67
FunctionSpaceType::RangeType RangeType
range type
Definition: space/shapefunctionset/orthonormal.hh:148
std::size_t constexpr size() const
return number of shape functions
Definition: space/shapefunctionset/orthonormal.hh:209
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
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: space/shapefunctionset/orthonormal.hh:150
FunctionSpace FunctionSpaceType
function space type
Definition: space/shapefunctionset/orthonormal.hh:75
FunctionSpaceType::DomainType DomainType
Definition: space/shapefunctionset/orthonormal.hh:146
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
static std::size_t constexpr size(int order)
please doc me
Definition: space/shapefunctionset/orthonormal.hh:212
OrthonormalShapeFunctionSet(const ThisType &)=default
OrthonormalShapeFunctionSet & operator=(const ThisType &)=default
OrthonormalShapeFunctionSet()=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
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