1#ifndef DUNE_FEM_SPACE_SHAPEFUNCTIONSET_LEGENDRE_HH
2#define DUNE_FEM_SPACE_SHAPEFUNCTIONSET_LEGENDRE_HH
30 template<
class FunctionSpace >
63 template<
class MultiIndex >
69 auto first = begin( multiIndex );
70 auto last = end( multiIndex );
71 assert( std::distance( first, last ) == dimDomain );
72 std::copy( first, last, multiIndex_.begin() );
80 return *std::max_element( multiIndex_.begin(), multiIndex_.end() );
84 const std::array< int, FunctionSpaceType::dimDomain > &
orders () const noexcept
93 for(
int i = 0; i < dimDomain; ++i )
101 for(
int k = 0; k < dimDomain; ++k )
105 for(
int i = 0; i < dimDomain; ++i )
106 jacobian[ 0 ][ i ] *= ( k == i ) ? dphi : phi;
114 for(
int k = 0; k < dimDomain; ++k )
118 for(
int i = 0; i < dimDomain; ++i )
121 for(
int j = i+1; j < dimDomain; ++j )
132 std::array< int, dimDomain > multiIndex_;
139 namespace __LegendreShapeFunctionSet
145 template<
class FunctionSpace >
148 typedef LegendreShapeFunction< FunctionSpace > ShapeFunctionType;
150 static const int dimDomain = FunctionSpace::dimDomain;
151 typedef std::array< int, dimDomain > MultiIndex;
154 explicit DefaultFactory (
int order )
158 int order () const noexcept {
return order_; }
160 std::size_t size () const noexcept
162 std::size_t size = 1;
163 for(
int i = 0; i < dimDomain; ++i )
168 template<
class InputIterator >
169 void operator() ( InputIterator first )
const noexcept
171 auto function = [&first](
const MultiIndex &multiIndex )
173 *first = ShapeFunctionType( multiIndex );
177 MultiIndex multiIndex;
178 fill( multiIndex, function, &multiIndex[ 0 ], dimDomain, order() );
182 template<
class Function >
183 static void fill ( MultiIndex &multiIndex, Function function,
184 int *begin, std::size_t n,
int order )
188 for( *begin = 0; *begin <= order; *begin += 1 )
189 fill( multiIndex, function, begin+1, n-1, order );
192 function( multiIndex );
216 template<
class FunctionSpace,
bool hierarchicalOrdering = false >
240 if( lhs.order() != rhs.order() )
241 return lhs.order() < rhs.order();
244 const auto &a = lhs.orders();
245 const auto &b = rhs.orders();
246 return std::lexicographical_compare( a.begin(), a.end(), b.begin(), b.end() );
287 template<
class Factory >
296 if( hierarchicalOrdering )
311 template<
class Po
int,
class Functor >
318 shapeFunction.evaluate(
coordinate( x ), value );
319 functor( i++, value );
324 template<
class Po
int,
class Functor >
331 shapeFunction.jacobian(
coordinate( x ), jacobian );
332 functor( i++, jacobian );
337 template<
class Po
int,
class Functor >
338 void hessianEach (
const Point &x, Functor functor )
const noexcept
344 shapeFunction.hessian(
coordinate( x ), hessian );
345 functor( i++, hessian );
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
implementation of a single scalar-valued Legendre shape function
Definition: shapefunctionset/legendre.hh:32
FunctionSpaceType::HessianRangeType HessianRangeType
hessian type
Definition: shapefunctionset/legendre.hh:51
void evaluate(const DomainType &x, RangeType &value) const noexcept
evaluate the function
Definition: shapefunctionset/legendre.hh:90
void jacobian(const DomainType &x, JacobianRangeType &jacobian) const noexcept
evaluate the Jacobian of the function
Definition: shapefunctionset/legendre.hh:98
void hessian(const DomainType &x, HessianRangeType &hessian) const noexcept
evaluate the hessian of the function
Definition: shapefunctionset/legendre.hh:111
FunctionSpaceType::RangeFieldType RangeFieldType
field type of range
Definition: shapefunctionset/legendre.hh:42
LegendreShapeFunction()=default
LegendreShapeFunction(const MultiIndex &multiIndex)
Definition: shapefunctionset/legendre.hh:64
FunctionSpace FunctionSpaceType
type of function space this function belongs to
Definition: shapefunctionset/legendre.hh:37
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian type
Definition: shapefunctionset/legendre.hh:49
int order() const noexcept
return polynomial order of this function
Definition: shapefunctionset/legendre.hh:78
FunctionSpaceType::DomainType DomainType
domain type
Definition: shapefunctionset/legendre.hh:45
FunctionSpaceType::RangeType RangeType
range type
Definition: shapefunctionset/legendre.hh:47
FunctionSpaceType::DomainFieldType DomainFieldType
field type of domain
Definition: shapefunctionset/legendre.hh:40
const std::array< int, FunctionSpaceType::dimDomain > & orders() const noexcept
return monomial orders of this function
Definition: shapefunctionset/legendre.hh:84
a Dune::Fem::ShapeFunctionSet of Legendre ansatz polynomials
Definition: shapefunctionset/legendre.hh:218
int order_
Definition: shapefunctionset/legendre.hh:351
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: shapefunctionset/legendre.hh:231
FunctionSpaceType::DomainType DomainType
domain type
Definition: shapefunctionset/legendre.hh:227
void hessianEach(const Point &x, Functor functor) const noexcept
evalute hessian of each shape function
Definition: shapefunctionset/legendre.hh:338
void evaluateEach(const Point &x, Functor functor) const noexcept
evalute each shape function
Definition: shapefunctionset/legendre.hh:312
FunctionSpaceType::RangeType RangeType
range type
Definition: shapefunctionset/legendre.hh:229
LegendreShapeFunctionSet()=default
default constructor resulting in uninitialized shape function set
int order() const noexcept
return order of shape functions
Definition: shapefunctionset/legendre.hh:305
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: shapefunctionset/legendre.hh:233
void jacobianEach(const Point &x, Functor functor) const noexcept
evalute jacobian of each shape function
Definition: shapefunctionset/legendre.hh:325
std::size_t size() const noexcept
return number of shape functions
Definition: shapefunctionset/legendre.hh:308
FunctionSpace FunctionSpaceType
function space type
Definition: shapefunctionset/legendre.hh:224
LegendreShapeFunctionSet(const Factory &factory)
initialize from user-defined factory object
Definition: shapefunctionset/legendre.hh:288
LegendreShapeFunction< FunctionSpace > ShapeFunctionType
Definition: shapefunctionset/legendre.hh:220
LegendreShapeFunctionSet(int order)
initialize with polynomial order
Definition: shapefunctionset/legendre.hh:264
std::vector< ShapeFunctionType > shapeFunctions_
Definition: shapefunctionset/legendre.hh:350
Definition: shapefunctionset/legendre.hh:237
bool operator()(const ShapeFunctionType &lhs, const ShapeFunctionType &rhs) noexcept
Definition: shapefunctionset/legendre.hh:238
static double hessian(const int num, const double x)
Definition: legendrepolynomials.hh:48
static double jacobian(const int num, const double x)
Definition: legendrepolynomials.hh:34
static double evaluate(const int num, const double x)
Definition: legendrepolynomials.hh:24