1#ifndef DUNE_FEM_SPACE_BASISFUNCTIONSET_TRANSFORMED_HH
2#define DUNE_FEM_SPACE_BASISFUNCTIONSET_TRANSFORMED_HH
9#include <dune/geometry/referenceelements.hh>
10#include <dune/geometry/type.hh>
42 template<
class Entity,
class ShapeFunctionSet,
class Transformation >
55 typedef typename GeometryType::ctype
ctype;
74 typedef std::decay_t<
decltype( Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) ) >
ReferenceElementType;
101 return Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(
type() );
107 template<
class QuadratureType,
class Vector,
class DofVector >
108 void axpy (
const QuadratureType &quad,
const Vector &values, DofVector &
dofs )
const
111 const unsigned int nop = quad.nop();
112 for(
unsigned int qp = 0; qp < nop; ++qp )
113 axpy( quad[ qp ], values[ qp ],
dofs );
122 template<
class QuadratureType,
class VectorA,
class VectorB,
class DofVector >
123 void axpy (
const QuadratureType &quad,
const VectorA &valuesA,
const VectorB &valuesB, DofVector &
dofs )
const
126 const unsigned int nop = quad.nop();
127 for(
unsigned int qp = 0; qp < nop; ++qp )
129 axpy( quad[ qp ], valuesA[ qp ],
dofs );
130 axpy( quad[ qp ], valuesB[ qp ],
dofs );
137 template<
class Po
int,
class DofVector >
148 template<
class Po
int,
class DofVector >
151 typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
153 const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed(
coordinate( x ) );
156 gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
166 template<
class Po
int,
class DofVector >
169 DUNE_THROW( NotImplemented,
"hessian axpy for TransformedBasisFunctionSet not implemented." );
175 template<
class Po
int,
class DofVector >
177 DofVector &
dofs )
const
184 template<
class QuadratureType,
class DofVector,
class RangeArray >
185 void evaluateAll (
const QuadratureType &quad,
const DofVector &
dofs, RangeArray &ranges )
const
188 const unsigned int nop = quad.nop();
189 for(
unsigned int qp = 0; qp < nop; ++qp )
194 template<
class Po
int,
class DofVector >
204 template<
class Po
int,
class RangeArray >
207 assert( values.size() >=
size() );
211 for(
auto &e : values )
212 e = transform.apply( e );
216 template<
class QuadratureType,
class DofVector,
class JacobianArray >
217 void jacobianAll (
const QuadratureType &quad,
const DofVector &
dofs, JacobianArray &jacobians )
const
220 const unsigned int nop = quad.nop();
221 for(
unsigned int qp = 0; qp < nop; ++qp )
226 template<
class Po
int,
class DofVector >
238 template<
class Po
int,
class JacobianRangeArray >
239 void jacobianAll (
const Point &x, JacobianRangeArray &jacobians )
const
241 assert( jacobians.size() >=
size() );
248 for(
auto &jacobian : jacobians )
249 jacobian = transform.apply( jacobian );
253 template<
class Po
int,
class DofVector >
256 DUNE_THROW( NotImplemented,
"hessianAll for TransformedBasisFunctionSet not implemented." );
260 template<
class Po
int,
class HessianRangeArray >
261 void hessianAll (
const Point &x, HessianRangeArray &hessians )
const
263 DUNE_THROW( NotImplemented,
"hessianAll for TransformedBasisFunctionSet not implemented." );
274 bool valid ()
const {
return bool(entity_); }
288 return Transformation(
geometry(), x );
Definition: bindguard.hh:11
void jacobianTransformation(const GeometryJacobianInverseTransposed &gjit, const FieldMatrix< K, ROWS, GeometryJacobianInverseTransposed::cols > &a, FieldMatrix< K, ROWS, GeometryJacobianInverseTransposed::rows > &b)
Definition: transformation.hh:21
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:14
Definition: explicitfieldvector.hh:75
Definition: misc/functor.hh:31
Definition: space/basisfunctionset/functor.hh:108
Definition: space/basisfunctionset/functor.hh:132
Definition: transformation.hh:36
implementation of a basis function set for given entity
Definition: transformed.hh:44
void hessianAll(const Point &x, HessianRangeArray &hessians) const
Definition: transformed.hh:261
const ShapeFunctionSetType & shapeFunctionSet() const
return shape function set
Definition: transformed.hh:284
void jacobianAll(const Point &x, JacobianRangeArray &jacobians) const
Definition: transformed.hh:239
ShapeFunctionSet ShapeFunctionSetType
shape function set type
Definition: transformed.hh:51
void hessianAll(const Point &x, const DofVector &dofs, HessianRangeType &hessian) const
Definition: transformed.hh:254
GeometryType geometry() const
Definition: transformed.hh:292
void evaluateAll(const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges) const
Definition: transformed.hh:185
EntityType::Geometry GeometryType
Definition: transformed.hh:54
int order() const
return order of basis function set
Definition: transformed.hh:92
void axpy(const QuadratureType &quad, const Vector &values, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:108
void evaluateAll(const Point &x, RangeArray &values) const
Definition: transformed.hh:205
GeometryType::JacobianTransposed JacobianTransposed
Definition: transformed.hh:56
std::size_t size() const
return size of basis function set
Definition: transformed.hh:95
void axpy(const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:149
TransformedBasisFunctionSet(const EntityType &entity, const ShapeFunctionSet &shapeFunctionSet=ShapeFunctionSet())
constructor
Definition: transformed.hh:82
void axpy(const QuadratureType &quad, const VectorA &valuesA, const VectorB &valuesB, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:123
GeometryType::ctype ctype
Definition: transformed.hh:55
bool valid() const
return true if entity pointer is set
Definition: transformed.hh:274
Entity EntityType
entity type
Definition: transformed.hh:49
Dune::GeometryType type() const
return geometry type
Definition: transformed.hh:277
auto referenceElement() const -> decltype(Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(std::declval< const Dune::GeometryType & >()))
return reference element
Definition: transformed.hh:98
FunctionSpaceType::RangeType RangeType
range type
Definition: transformed.hh:65
void axpy(const Point &x, const RangeType &valueFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:138
std::decay_t< decltype(Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(std::declval< const Dune::GeometryType & >())) > ReferenceElementType
type of reference element
Definition: transformed.hh:74
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: transformed.hh:67
ShapeFunctionSetType::FunctionSpaceType FunctionSpaceType
type of function space
Definition: transformed.hh:60
FunctionSpaceType::DomainType DomainType
domain type
Definition: transformed.hh:63
Transformation transformation(const DomainType &x) const
Definition: transformed.hh:286
const Entity & entity() const
return entity
Definition: transformed.hh:267
void evaluateAll(const Point &x, const DofVector &dofs, RangeType &value) const
Definition: transformed.hh:195
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: transformed.hh:69
void axpy(const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:176
TransformedBasisFunctionSet()
constructor
Definition: transformed.hh:77
void jacobianAll(const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians) const
Definition: transformed.hh:217
void axpy(const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:167
void jacobianAll(const Point &x, const DofVector &dofs, JacobianRangeType &jacobian) const
Definition: transformed.hh:227
FunctionSpaceType::RangeFieldType RangeFieldType
Definition: transformed.hh:71
A vector valued function space.
Definition: functionspace.hh:60
FunctionSpaceTraits::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:71
@ 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
Interface class for shape function sets.
Definition: shapefunctionset/shapefunctionset.hh:33
void evaluateEach(const Point &x, Functor functor) const
evalute each shape function
std::size_t size() const
return number of shape functions
int order() const
return order of shape functions
void jacobianEach(const Point &x, Functor functor) const
evalute jacobian of each shape function