1#ifndef DUNE_FEM_SHAPEFUNCTIONSET_TENSORPRODUCT_HH
2#define DUNE_FEM_SHAPEFUNCTIONSET_TENSORPRODUCT_HH
10#include <dune/common/fmatrix.hh>
12#include <dune/common/fvector.hh>
13#include <dune/common/hybridutilities.hh>
14#include <dune/common/tupleutility.hh>
26 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
32 "dimDomain of FunctionSpace must coincide with length of ShapeFunctionSetTuple." );
34 "FunctionSpace must be scalar (i.e., dimRange = 1)." );
38 template<
int i >
struct Size;
39 template<
int i >
struct EvaluateAll;
40 template<
int i >
struct JacobianAll;
41 template<
int i >
struct HessianAll;
67 std::size_t
size ()
const;
69 template<
class Po
int,
class Functor >
70 void evaluateEach (
const Point &x, Functor functor )
const;
72 template<
class Po
int,
class Functor >
73 void jacobianEach (
const Point &x, Functor functor )
const;
75 template<
class Po
int,
class Functor >
76 void hessianEach (
const Point &x, Functor functor )
const;
79 template<
class Functor >
80 void doEvaluateEach (
int d,
RangeType value, std::size_t &index,
const RangeFieldType *buffer, Functor functor )
const;
81 template<
class Functor >
83 template<
class Functor >
86 ShapeFunctionSetTuple shapeFunctionSetTuple_;
87 std::array< std::size_t, dimension > sizes_;
96 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
103 buffer_[ i ] = value;
107 void operator() (
const std::size_t i,
const FieldVector< T, 1 > &value )
109 (*this)( i, value[ 0 ] );
113 void operator() (
const std::size_t i,
const FieldMatrix< T, 1, 1 > &value )
115 (*this)( i, value[ 0 ][ 0 ] );
127 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
131 static void apply (
const ShapeFunctionSetTuple &tuple, std::array< std::size_t, FunctionSpace::dimDomain > &
size )
133 size[ i ] = std::get< i >( tuple ).size();
142 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
144 struct TensorProductShapeFunctionSet< FunctionSpace, ShapeFunctionSetTuple >::EvaluateAll
148 Dune::FieldVector< DomainFieldType, 1 > xi( x[ i ] );
149 std::get< i >( tuple ).evaluateEach( xi, Assign( it ) );
150 it += std::get< i >( tuple ).size();
159 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
161 struct TensorProductShapeFunctionSet< FunctionSpace, ShapeFunctionSetTuple >::JacobianAll
165 Dune::FieldVector< DomainFieldType, 1 > xi( x[ i ] );
166 const std::size_t
size = std::get< i >( tuple ).size();
167 std::get< i >( tuple ).evaluateEach( xi, Assign( it ) );
168 std::get< i >( tuple ).jacobianEach( xi, Assign( it+
size ) );
178 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
180 struct TensorProductShapeFunctionSet< FunctionSpace, ShapeFunctionSetTuple >::HessianAll
184 Dune::FieldVector< DomainFieldType, 1 > xi( x[ i ] );
185 const std::size_t
size = std::get< i >( tuple ).size();
186 std::get< i >( tuple ).evaluateEach( xi, Assign( it ) );
187 std::get< i >( tuple ).jacobianEach( xi, Assign( it+
size ) );
188 std::get< i >( tuple ).hessianEach( xi, Assign( it+2*
size ) );
198 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
199 inline TensorProductShapeFunctionSet< FunctionSpace, ShapeFunctionSetTuple >
201 : shapeFunctionSetTuple_( shapeFunctionSetTuple )
204 std::size_t buffer_size = 0;
205 for(
int i = 0; i < dimension; ++i )
206 buffer_size += sizes_[ i ];
211 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
219 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
222 : shapeFunctionSetTuple_( other.shapeFunctionSetTuple_ )
224 std::size_t buffer_size = 0;
225 for(
int i = 0; i < dimension; ++i )
227 sizes_[ i ] = other.sizes_[ i ];
228 buffer_size += sizes_[ i ];
234 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
243 shapeFunctionSetTuple_ = other.shapeFunctionSetTuple_;
244 std::size_t buffer_size = 0;
245 for(
int i = 0; i < dimension; ++i )
247 sizes_[ i ] = other.sizes_[ i ];
248 buffer_size += sizes_[ i ];
255 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
259 Hybrid::forEach( std::make_index_sequence< std::tuple_size< ShapeFunctionSetTupleType >::value >{},
260 [ & ](
auto i ){ value =
std::max( value, std::get< i >( shapeFunctionSetTuple_ ).order() ); } );
265 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
269 std::size_t size( 1 );
270 for(
int i = 0; i < dimension; ++i )
276 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
277 template<
class Po
int,
class Functor >
284 std::size_t index = 0;
289 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
290 template<
class Po
int,
class Functor >
297 std::size_t index = 0;
299 for(
int i = 0; i < dimension; ++i )
301 doJacobianEach( 0, jacobian, index, buffer_, functor );
305 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
306 template<
class Po
int,
class Functor >
313 std::size_t index = 0;
315 for(
int i = 0; i < dimension; ++i )
316 for(
int j = 0; j < dimension; ++j )
318 doHessianEach( 0, hessian, index, buffer_, functor );
322 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
323 template<
class Functor >
325 ::doEvaluateEach (
int d, RangeType value, std::size_t &index,
const RangeFieldType *buffer, Functor functor )
const
329 for( std::size_t i = 0; i < sizes_[ d ]; ++i )
331 RangeType v( value );
332 v[ 0 ] *= buffer[ i ];
333 doEvaluateEach( d+1, v, index, buffer+sizes_[ d ], functor );
337 functor( index++, value );
341 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
342 template<
class Functor >
343 inline void TensorProductShapeFunctionSet< FunctionSpace, ShapeFunctionSetTuple >
344 ::doJacobianEach (
int d, JacobianRangeType jacobian, std::size_t &index,
const RangeFieldType *buffer, Functor functor )
const
348 for( std::size_t i = 0; i < sizes_[ d ]; ++i )
350 JacobianRangeType j( jacobian );
351 j[ 0 ][ d ] *= buffer[ i + sizes_[ d ] ];
352 for(
int k = 1; k < dimension; ++k )
353 j[ 0 ][ (d+k)%dimension ] *= buffer[ i ];
354 doJacobianEach( d+1, j, index, buffer+2*sizes_[ d ], functor );
358 functor( index++, jacobian );
362 template<
class FunctionSpace,
class ShapeFunctionSetTuple >
363 template<
class Functor >
364 inline void TensorProductShapeFunctionSet< FunctionSpace, ShapeFunctionSetTuple >
365 ::doHessianEach (
int d, HessianRangeType hessian, std::size_t &index,
const RangeFieldType *buffer, Functor functor )
const
369 for( std::size_t i = 0; i < sizes_[ d ]; ++i )
371 HessianRangeType h( hessian );
372 h[ 0 ][ d ][ d ] *= buffer[ i + 2*sizes_[ d ] ];
373 for(
int j = 1; j < dimension; ++j )
375 h[ 0 ][ (d+j)%dimension ][ d ] *= buffer[ i * sizes_[ d ] ];
376 h[ 0 ][ d ][ (d+j)%dimension ] *= buffer[ i * sizes_[ d ] ];
377 for(
int k = 1; k < dimension; ++k )
378 h[ 0 ][ (d+j)%dimension ][ (d+k)%dimension ] *= buffer[ i ];
380 doHessianEach( d+1, h, index, buffer+3*sizes_[ d ], functor );
384 functor( index++, hessian );
double max(const Dune::Fem::Double &v, const double p)
Definition: double.hh:965
Definition: bindguard.hh:11
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:14
static void forEach(IndexRange< T, sz > range, F &&f)
Definition: hybrid.hh:129
static DUNE_PRIVATE void apply(Args &&... args)
Definition: forloop.hh:23
A vector valued function space.
Definition: functionspace.hh:60
ExplicitFieldVector< FieldMatrix< RangeFieldType, dimDomain, dimDomain >, dimRange > HessianRangeType
Intrinsic type used for the hessian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:79
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: tensorproduct.hh:28
void evaluateEach(const Point &x, Functor functor) const
Definition: tensorproduct.hh:279
FunctionSpaceType::RangeType RangeType
Definition: tensorproduct.hh:53
void jacobianEach(const Point &x, Functor functor) const
Definition: tensorproduct.hh:292
FunctionSpaceType::HessianRangeType HessianRangeType
Definition: tensorproduct.hh:55
FunctionSpaceType::DomainType DomainType
Definition: tensorproduct.hh:52
TensorProductShapeFunctionSet()=default
void hessianEach(const Point &x, Functor functor) const
Definition: tensorproduct.hh:308
int order() const
Definition: tensorproduct.hh:256
std::size_t size() const
Definition: tensorproduct.hh:267
FunctionSpaceType::JacobianRangeType JacobianRangeType
Definition: tensorproduct.hh:54
const ThisType & operator=(const ThisType &other)
Definition: tensorproduct.hh:237
FunctionSpace FunctionSpaceType
Definition: tensorproduct.hh:46
FunctionSpaceType::RangeFieldType RangeFieldType
Definition: tensorproduct.hh:50
FunctionSpaceType::DomainFieldType DomainFieldType
Definition: tensorproduct.hh:49
~TensorProductShapeFunctionSet()
Definition: tensorproduct.hh:213
ShapeFunctionSetTuple ShapeFunctionSetTupleType
Definition: tensorproduct.hh:47
Definition: tensorproduct.hh:98
Assign(RangeFieldType *buffer)
Definition: tensorproduct.hh:99