1 #ifndef DUNE_FEM_SPACE_SHAPEFUNCTIONSET_ORTHONORMAL_HH 2 #define DUNE_FEM_SPACE_SHAPEFUNCTIONSET_ORTHONORMAL_HH 7 #include <dune/common/exceptions.hh> 9 #include <dune/geometry/genericgeometry/topologytypes.hh> 10 #include <dune/geometry/type.hh> 33 template<
class FunctionSpace,
int polOrder >
41 template<
class FunctionSpace,
int polOrder >
46 template<
int order,
int dimension >
47 struct ShapeFunctionSetSize;
50 struct ShapeFunctionSetSize< order, 1 >
52 static const std::size_t
v = order + 1;
56 struct ShapeFunctionSetSize< order, 2 >
58 static const std::size_t
v = (order + 2) * (order + 1) / 2;
62 struct ShapeFunctionSetSize< order, 3 >
64 static const size_t v = ((order+1)*(order+2)*(2*order+3)/6
65 + (order+1)*(order+2)/2)/2;
69 static const std::size_t
v = ShapeFunctionSetSize< polOrder, FunctionSpace::dimDomain >::v;
77 template<
class FunctionSpace,
int polOrder >
81 "FunctionSpace must be scalar (i.e., dimRange = 1)." );
95 typedef Dune::GenericGeometry::Prism< Dune::GenericGeometry::Point >
Line;
99 typedef Dune::GenericGeometry::Pyramid< Line >
Triangle;
101 typedef Dune::GenericGeometry::Prism< Quadrilateral >
Hexahedron;
103 typedef Dune::GenericGeometry::Pyramid< Quadrilateral >
Pyramid;
105 typedef Dune::GenericGeometry::Prism< Triangle >
Prism;
110 template<
class Topology >
113 template<
unsigned int topologyId >
116 static const unsigned int v = topologyId | (
unsigned int)Dune::GenericGeometry::prismConstruction;
120 typedef typename Dune::GenericGeometry::Topology< UniqueId< Topology::id >::v, Topology::dimension >::type
Type;
124 template<
class Topology >
131 template<
class Topology >
134 template<
class Topology >
137 template<
class Topology >
146 template<
class FunctionSpace,
int polOrder >
147 template<
class Topology >
150 template<
class Functor >
154 for( std::size_t i = 0; i < size; ++i )
156 assert( x.size() == Topology::dimension );
157 functor( i, evaluate( basicGeometryType< Topology>(), i, x ) );
169 return OrthonormalBase1d::eval_line( i, &x[ 0 ] );
174 return OrthonormalBase2d::eval_quadrilateral_2d( i, &x[ 0 ] );
179 return OrthonormalBase2d::eval_triangle_2d( i, &x[ 0 ] );
184 return OrthonormalBase3d::eval_pyramid_3d( i, &x[ 0 ] );
189 return OrthonormalBase3d::eval_hexahedron_3d( i, &x[ 0 ] );
194 return OrthonormalBase3d::eval_prism_3d( i, &x[ 0 ] );
199 return OrthonormalBase3d::eval_tetrahedron_3d( i, &x[ 0 ] );
208 template<
class FunctionSpace,
int polOrder >
209 template<
class Topology >
212 template<
class Functor >
217 for( std::size_t i = 0; i < size; ++i )
219 assert( x.size() == Topology::dimension );
220 evaluate( basicGeometryType< Topology >(), i, x, jacobian );
221 functor( i, jacobian );
232 OrthonormalBase1d::grad_line( i , &x[ 0 ], &jacobian[ 0 ][ 0 ] );
237 OrthonormalBase2d::grad_quadrilateral_2d( i , &x[ 0 ], &jacobian[ 0 ][ 0 ] );
242 OrthonormalBase2d::grad_triangle_2d( i , &x[ 0 ], &jacobian[ 0 ][ 0 ] );
247 OrthonormalBase3d::grad_pyramid_3d( i , &x[ 0 ], &jacobian[ 0 ][ 0 ] );
252 OrthonormalBase3d::grad_hexahedron_3d( i , &x[ 0 ], &jacobian[ 0 ][ 0 ] );
257 OrthonormalBase3d::grad_prism_3d( i , &x[ 0 ], &jacobian[ 0 ][ 0 ] );
262 OrthonormalBase3d::grad_tetrahedron_3d( i , &x[ 0 ], &jacobian[ 0 ][ 0 ] );
271 template<
class FunctionSpace,
int polOrder >
272 template<
class Topology >
275 template<
class Functor >
280 for( std::size_t i = 0; i < size; ++i )
282 assert( x.size() == Topology::dimension );
283 evaluate( basicGeometryType< Topology >(), i, x, hessian );
284 functor( i, hessian );
293 OrthonormalBase2d::hess_quadrilateral_2d( i , &x[ 0 ], values );
297 hessian[ 0 ][ j ][ k ] = values[ j + k ];
304 OrthonormalBase2d::hess_triangle_2d( i , &x[ 0 ], values );
308 hessian[ 0 ][ j ][ k ] = values[ j + k ];
311 template<
class Other >
314 DUNE_THROW( NotImplemented,
"On orthonormal shape function set HessianAll() is only implemented for triangles" );
323 template<
class FunctionSpace,
int polOrder >
327 "FunctionSpace must be scalar (i.e., dimRange = 1)." );
340 static const int dimension = FunctionSpaceType::dimDomain;
352 : topologyId_( type.id() )
356 int order ()
const {
return polOrder; }
366 template<
class Po
int,
class Functor >
370 Dune::GenericGeometry::IfTopology< ShapeFunctionSetHelperType::template EvaluateEach, dimension >
371 ::apply( topologyId_, y, functor );
375 template<
class Po
int,
class Functor >
379 Dune::GenericGeometry::IfTopology< ShapeFunctionSetHelperType::template JacobianEach, dimension >
380 ::apply( topologyId_, y, functor );
384 template<
class Po
int,
class Functor >
388 Dune::GenericGeometry::IfTopology< ShapeFunctionSetHelperType::template HessianEach, dimension >
389 ::apply( topologyId_, y, functor );
393 unsigned int topologyId_;
400 #endif // #ifndef DUNE_FEM_SPACE_SHAPEFUNCTIONSET_ORTHONORMAL_HH static RangeType evaluate(const Hexahedron &, std::size_t i, const DomainType &x)
Definition: space/shapefunctionset/orthonormal.hh:187
void jacobianEach(const Point &x, Functor functor) const
evalute jacobian of each shape function
Definition: space/shapefunctionset/orthonormal.hh:376
OrthonormalBase_1D< DomainFieldType, RangeFieldType > OrthonormalBase1d
Definition: space/shapefunctionset/orthonormal.hh:226
static BasicGeometryType< Topology >::Type basicGeometryType()
Definition: space/shapefunctionset/orthonormal.hh:125
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:62
Definition: space/shapefunctionset/orthonormal.hh:138
static void apply(const DomainType &x, Functor functor)
Definition: space/shapefunctionset/orthonormal.hh:151
OrthonormalBase_2D< DomainFieldType, RangeFieldType > OrthonormalBase2d
Definition: space/shapefunctionset/orthonormal.hh:227
A vector valued function space.
Definition: functionspace.hh:16
Dune::GenericGeometry::Prism< Dune::GenericGeometry::Point > Line
Definition: space/shapefunctionset/orthonormal.hh:95
Dune::GenericGeometry::Prism< Line > Quadrilateral
Definition: space/shapefunctionset/orthonormal.hh:97
Definition: space/shapefunctionset/orthonormal.hh:42
static RangeType evaluate(const Tetrahedron &, std::size_t i, const DomainType &x)
Definition: space/shapefunctionset/orthonormal.hh:197
static RangeType evaluate(const Line &, std::size_t i, const DomainType &x)
Definition: space/shapefunctionset/orthonormal.hh:167
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:74
Definition: space/shapefunctionset/orthonormal.hh:78
FieldVector< FieldMatrix< RangeFieldType, dimDomain, dimDomain >, dimRange > HessianRangeType
Intrinsic type used for the hessian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:78
std::size_t size() const
return number of shape functions
Definition: space/shapefunctionset/orthonormal.hh:360
OrthonormalBase_1D< DomainFieldType, RangeFieldType > OrthonormalBase1d
Definition: space/shapefunctionset/orthonormal.hh:163
FunctionSpaceType::JacobianRangeType JacobianRangeType
Definition: space/shapefunctionset/orthonormal.hh:90
FunctionSpaceType::RangeFieldType RangeFieldType
Definition: space/shapefunctionset/orthonormal.hh:86
FunctionSpaceType::RangeType RangeType
range type
Definition: space/shapefunctionset/orthonormal.hh:345
Definition: orthonormalbase_2d.hh:16
dimension of domain vector space
Definition: functionspaceinterface.hh:45
FunctionSpaceType::DomainType DomainType
domain type
Definition: space/shapefunctionset/orthonormal.hh:343
static RangeType evaluate(const Quadrilateral &, std::size_t i, const DomainType &x)
Definition: space/shapefunctionset/orthonormal.hh:172
Dune::GenericGeometry::Pyramid< Line > Triangle
Definition: space/shapefunctionset/orthonormal.hh:99
dimension of range vector space
Definition: functionspaceinterface.hh:47
void evaluateEach(const Point &x, Functor functor) const
evalute each shape function
Definition: space/shapefunctionset/orthonormal.hh:367
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::DomainFieldType DomainFieldType
Intrinsic type used for values in the domain field (usually a double)
Definition: functionspaceinterface.hh:59
static void apply(const DomainType &x, Functor functor)
Definition: space/shapefunctionset/orthonormal.hh:276
static const std::size_t v
Definition: space/shapefunctionset/orthonormal.hh:69
Dune::GenericGeometry::Topology< UniqueId< Topology::id >::v, Topology::dimension >::type Type
Definition: space/shapefunctionset/orthonormal.hh:120
static void evaluate(const Pyramid &, std::size_t i, const DomainType &x, JacobianRangeType &jacobian)
Definition: space/shapefunctionset/orthonormal.hh:245
Definition: orthonormalbase_1d.hh:14
FunctionSpaceType::DomainFieldType DomainFieldType
Definition: space/shapefunctionset/orthonormal.hh:85
int order() const
return order of shape functions
Definition: space/shapefunctionset/orthonormal.hh:356
Definition: coordinate.hh:4
FunctionSpace FunctionSpaceType
Definition: space/shapefunctionset/orthonormal.hh:81
OrthonormalBase_2D< DomainFieldType, RangeFieldType > OrthonormalBase2d
Definition: space/shapefunctionset/orthonormal.hh:164
FunctionSpaceType::HessianRangeType HessianRangeType
Definition: space/shapefunctionset/orthonormal.hh:91
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:66
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: space/shapefunctionset/orthonormal.hh:347
Definition: space/shapefunctionset/orthonormal.hh:111
OrthonormalBase_3D< DomainFieldType, RangeFieldType > OrthonormalBase3d
Definition: space/shapefunctionset/orthonormal.hh:165
static RangeType evaluate(const Prism &, std::size_t i, const DomainType &x)
Definition: space/shapefunctionset/orthonormal.hh:192
FunctionSpaceType::RangeType RangeType
Definition: space/shapefunctionset/orthonormal.hh:89
static void evaluate(const Line &, std::size_t i, const DomainType &x, JacobianRangeType &jacobian)
Definition: space/shapefunctionset/orthonormal.hh:230
static void evaluate(const Triangle &, std::size_t i, const DomainType &x, HessianRangeType &hessian)
Definition: space/shapefunctionset/orthonormal.hh:300
OrthonormalShapeFunctionSet(const GeometryType &type)
Definition: space/shapefunctionset/orthonormal.hh:351
static RangeType evaluate(const Triangle &, std::size_t i, const DomainType &x)
Definition: space/shapefunctionset/orthonormal.hh:177
Dune::GenericGeometry::Prism< Quadrilateral > Hexahedron
Definition: space/shapefunctionset/orthonormal.hh:101
Definition: space/shapefunctionset/orthonormal.hh:132
Dune::GenericGeometry::Pyramid< Triangle > Tetrahedron
Definition: space/shapefunctionset/orthonormal.hh:107
static void apply(const DomainType &x, Functor functor)
Definition: space/shapefunctionset/orthonormal.hh:213
static RangeType evaluate(const Pyramid &, std::size_t i, const DomainType &x)
Definition: space/shapefunctionset/orthonormal.hh:182
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:70
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: space/shapefunctionset/orthonormal.hh:349
Dune::GenericGeometry::Pyramid< Quadrilateral > Pyramid
Definition: space/shapefunctionset/orthonormal.hh:103
FunctionSpace FunctionSpaceType
function space type
Definition: space/shapefunctionset/orthonormal.hh:337
FunctionSpaceType::DomainType DomainType
Definition: space/shapefunctionset/orthonormal.hh:88
static void evaluate(const Triangle &, std::size_t i, const DomainType &x, JacobianRangeType &jacobian)
Definition: space/shapefunctionset/orthonormal.hh:240
static void evaluate(const Other &, std::size_t i, const DomainType &x, HessianRangeType &hessian)
Definition: space/shapefunctionset/orthonormal.hh:312
Definition: space/shapefunctionset/orthonormal.hh:135
static void evaluate(const Tetrahedron &, std::size_t i, const DomainType &x, JacobianRangeType &jacobian)
Definition: space/shapefunctionset/orthonormal.hh:260
static void evaluate(const Quadrilateral &, std::size_t i, const DomainType &x, HessianRangeType &hessian)
Definition: space/shapefunctionset/orthonormal.hh:289
Definition: orthonormalbase_3d.hh:16
static void evaluate(const Quadrilateral &, std::size_t i, const DomainType &x, JacobianRangeType &jacobian)
Definition: space/shapefunctionset/orthonormal.hh:235
Dune::GenericGeometry::Prism< Triangle > Prism
Definition: space/shapefunctionset/orthonormal.hh:105
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:11
static void evaluate(const Hexahedron &, std::size_t i, const DomainType &x, JacobianRangeType &jacobian)
Definition: space/shapefunctionset/orthonormal.hh:250
OrthonormalBase_3D< DomainFieldType, RangeFieldType > OrthonormalBase3d
Definition: space/shapefunctionset/orthonormal.hh:228
void hessianEach(const Point &x, Functor functor) const
evalute hessian of each shape function
Definition: space/shapefunctionset/orthonormal.hh:385
static void evaluate(const Prism &, std::size_t i, const DomainType &x, JacobianRangeType &jacobian)
Definition: space/shapefunctionset/orthonormal.hh:255
Definition: space/shapefunctionset/orthonormal.hh:34