1 #ifndef DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH
2 #define DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH
12 #include <dune/geometry/referenceelements.hh>
13 #include <dune/geometry/type.hh>
45 template<
class Entity,
class ShapeFunctionSet >
65 typedef typename GeometryType::ctype
ctype;
73 typedef typename FunctionSpaceType::DomainType
DomainType;
76 typedef typename FunctionSpaceType::RangeType
RangeType;
83 typedef std::decay_t< decltype( Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) ) >
ReferenceElementType;
110 return Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(
type() );
116 template<
class QuadratureType,
class Vector,
class DofVector >
117 void axpy (
const QuadratureType &quad,
const Vector &values, DofVector &
dofs )
const
120 const unsigned int nop = quad.nop();
121 for(
unsigned int qp = 0; qp < nop; ++qp )
123 axpy( quad[ qp ], values[ qp ],
dofs );
133 template<
class QuadratureType,
class VectorA,
class VectorB,
class DofVector >
134 void axpy (
const QuadratureType &quad,
const VectorA &valuesA,
const VectorB &valuesB, DofVector &
dofs )
const
137 const unsigned int nop = quad.nop();
138 for(
unsigned int qp = 0; qp < nop; ++qp )
140 axpy( quad[ qp ], valuesA[ qp ],
dofs );
141 axpy( quad[ qp ], valuesB[ qp ],
dofs );
148 template<
class Po
int,
class DofVector >
158 template<
class Po
int,
class DofVector >
161 typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
163 const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed(
coordinate( x ) );
165 for(
int r = 0; r < FunctionSpaceType::dimRange; ++r )
166 gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
173 template<
class Po
int,
class DofVector >
176 typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
178 const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed(
coordinate( x ) );
183 Dune::FieldMatrix<double,DomainType::dimension,DomainType::dimension>
186 for(
int r = 0; r < FunctionSpaceType::dimRange; ++r )
187 for(
int j = 0; j < gjit.cols; ++j )
188 for(
int k = 0; k < gjit.cols; ++k )
190 hessianFactor[r].mv(G[k],Hg);
191 tmpHessianFactor[r][j][k] += Hg * G[j];
200 template<
class Po
int,
class DofVector >
202 DofVector &
dofs )
const
209 template<
class QuadratureType,
class DofVector,
class RangeArray >
210 void evaluateAll (
const QuadratureType &quad,
const DofVector &
dofs, RangeArray &ranges )
const
213 const unsigned int nop = quad.nop();
214 for(
unsigned int qp = 0; qp < nop; ++qp )
221 template<
class Po
int,
class DofVector >
230 template<
class Po
int,
class RangeArray >
233 assert( values.size() >=
size() );
239 template<
class QuadratureType,
class DofVector,
class JacobianArray >
240 void jacobianAll (
const QuadratureType &quad,
const DofVector &
dofs, JacobianArray &jacobians )
const
243 const unsigned int nop = quad.nop();
244 for(
unsigned int qp = 0; qp < nop; ++qp )
251 template<
class Po
int,
class DofVector >
259 Transformation transformation( geo,
coordinate( x ) );
260 transformation( localJacobian, jacobian );
264 template<
class Po
int,
class JacobianRangeArray >
265 void jacobianAll (
const Point &x, JacobianRangeArray &jacobians )
const
267 assert( jacobians.size() >=
size() );
270 Transformation transformation( geo,
coordinate( x ) );
276 template<
class Po
int,
class DofVector >
284 Transformation transformation( geo,
coordinate( x ) );
285 transformation( localHessian, hessian );
289 template<
class Po
int,
class HessianRangeArray >
290 void hessianAll (
const Point &x, HessianRangeArray &hessians )
const
292 assert( hessians.size() >=
size() );
295 Transformation transformation( geo,
coordinate( x ) );
Definition: bindguard.hh:11
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:11
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
Definition: explicitfieldvector.hh:75
Definition: misc/functor.hh:31
Definition: basisfunctionset/default.hh:47
FunctionSpaceType::DomainType DomainType
domain type
Definition: basisfunctionset/default.hh:73
void evaluateAll(const Point &x, const DofVector &dofs, RangeType &value) const
Definition: basisfunctionset/default.hh:222
LocalFunctionSpaceType::HessianRangeType LocalHessianRangeType
Definition: basisfunctionset/default.hh:59
DefaultBasisFunctionSet()
constructor
Definition: basisfunctionset/default.hh:86
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: basisfunctionset/default.hh:80
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: basisfunctionset/default.hh:134
const Entity & entity() const
return entity
Definition: basisfunctionset/default.hh:301
GeometryType geometry() const
Definition: basisfunctionset/default.hh:318
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: basisfunctionset/default.hh:201
FunctionSpaceType::RangeType RangeType
range type
Definition: basisfunctionset/default.hh:76
ShapeFunctionSetType::FunctionSpaceType LocalFunctionSpaceType
Definition: basisfunctionset/default.hh:57
GeometryType::ctype ctype
Definition: basisfunctionset/default.hh:65
LocalFunctionSpaceType::RangeFieldType RangeFieldType
Definition: basisfunctionset/default.hh:61
LocalFunctionSpaceType::JacobianRangeType LocalJacobianRangeType
Definition: basisfunctionset/default.hh:58
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: basisfunctionset/default.hh:159
void jacobianAll(const Point &x, JacobianRangeArray &jacobians) const
Definition: basisfunctionset/default.hh:265
DefaultBasisFunctionSet(const EntityType &entity, const ShapeFunctionSet &shapeFunctionSet=ShapeFunctionSet())
constructor
Definition: basisfunctionset/default.hh:91
const ShapeFunctionSetType & shapeFunctionSet() const
return shape function set
Definition: basisfunctionset/default.hh:315
void axpy(const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs) const
Add H:D^2phi to each dof.
Definition: basisfunctionset/default.hh:174
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: basisfunctionset/default.hh:78
EntityType::Geometry GeometryType
Definition: basisfunctionset/default.hh:63
void hessianAll(const Point &x, HessianRangeArray &hessians) const
Definition: basisfunctionset/default.hh:290
void evaluateAll(const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges) const
Definition: basisfunctionset/default.hh:210
ToNewDimDomainFunctionSpace< LocalFunctionSpaceType, EntityType ::Geometry ::coorddimension >::Type FunctionSpaceType
type of function space
Definition: basisfunctionset/default.hh:70
void hessianAll(const Point &x, const DofVector &dofs, HessianRangeType &hessian) const
Definition: basisfunctionset/default.hh:277
void evaluateAll(const Point &x, RangeArray &values) const
Definition: basisfunctionset/default.hh:231
ShapeFunctionSet ShapeFunctionSetType
shape function set type
Definition: basisfunctionset/default.hh:54
std::size_t size() const
return size of basis function set
Definition: basisfunctionset/default.hh:104
void jacobianAll(const Point &x, const DofVector &dofs, JacobianRangeType &jacobian) const
Definition: basisfunctionset/default.hh:252
Entity EntityType
entity type
Definition: basisfunctionset/default.hh:52
std::decay_t< decltype(Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(std::declval< const Dune::GeometryType & >))) > ReferenceElementType
type of reference element
Definition: basisfunctionset/default.hh:83
Dune::GeometryType type() const
return geometry type
Definition: basisfunctionset/default.hh:308
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: basisfunctionset/default.hh:149
void jacobianAll(const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians) const
Definition: basisfunctionset/default.hh:240
auto referenceElement() const -> decltype(Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(std::declval< const Dune::GeometryType & >()))
return reference element
Definition: basisfunctionset/default.hh:107
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: basisfunctionset/default.hh:117
int order() const
return order of basis function set
Definition: basisfunctionset/default.hh:101
Definition: space/basisfunctionset/functor.hh:108
Definition: space/basisfunctionset/functor.hh:132
Definition: transformation.hh:36
Definition: transformation.hh:92
A vector valued function space.
Definition: functionspace.hh:60
convert functions space to space with new dim domain
Definition: functionspace.hh:246
FunctionSpaceTraits::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:75
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 hessianEach(const Point &x, Functor functor) const
evalute hessian of each shape function
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