1 #ifndef DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH 2 #define DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH 9 #include <dune/common/nullptr.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;
109 return Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(
type() );
115 template<
class QuadratureType,
class Vector,
class DofVector >
116 void axpy (
const QuadratureType &quad,
const Vector &values, DofVector &dofs )
const 119 const unsigned int nop = quad.nop();
120 for(
unsigned int qp = 0; qp < nop; ++qp )
122 axpy( quad[ qp ], values[ qp ], dofs );
132 template<
class QuadratureType,
class VectorA,
class VectorB,
class DofVector >
133 void axpy (
const QuadratureType &quad,
const VectorA &valuesA,
const VectorB &valuesB, DofVector &dofs )
const 136 const unsigned int nop = quad.nop();
137 for(
unsigned int qp = 0; qp < nop; ++qp )
139 axpy( quad[ qp ], valuesA[ qp ], dofs );
140 axpy( quad[ qp ], valuesB[ qp ], dofs );
147 template<
class Po
int,
class DofVector >
148 void axpy (
const Point &x,
const RangeType &valueFactor, DofVector &dofs )
const 157 template<
class Po
int,
class DofVector >
158 void axpy (
const Point &x,
const JacobianRangeType &jacobianFactor, DofVector &dofs )
const 160 typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
161 const GeometryType &geo =
geometry();
162 const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed(
coordinate( x ) );
164 for(
int r = 0; r < FunctionSpaceType::dimRange; ++r )
165 gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
174 template<
class Po
int,
class DofVector >
175 void axpy (
const Point &x,
const RangeType &valueFactor,
const JacobianRangeType &jacobianFactor,
176 DofVector &dofs )
const 178 axpy( x, valueFactor, dofs );
179 axpy( x, jacobianFactor, dofs );
183 template<
class QuadratureType,
class DofVector,
class RangeArray >
184 void evaluateAll (
const QuadratureType &quad,
const DofVector &dofs, RangeArray &ranges )
const 187 const unsigned int nop = quad.nop();
188 for(
unsigned int qp = 0; qp < nop; ++qp )
195 template<
class Po
int,
class DofVector >
196 void evaluateAll (
const Point &x,
const DofVector &dofs, RangeType &value )
const 204 template<
class Po
int,
class RangeArray >
207 assert( values.size() >=
size() );
213 template<
class QuadratureType,
class DofVector,
class JacobianArray >
214 void jacobianAll (
const QuadratureType &quad,
const DofVector &dofs, JacobianArray &jacobians )
const 217 const unsigned int nop = quad.nop();
218 for(
unsigned int qp = 0; qp < nop; ++qp )
225 template<
class Po
int,
class DofVector >
226 void jacobianAll (
const Point &x,
const DofVector &dofs, JacobianRangeType &jacobian )
const 231 const GeometryType &geo =
geometry();
233 Transformation transformation( geo,
coordinate( x ) );
234 transformation( localJacobian, jacobian );
238 template<
class Po
int,
class JacobianRangeArray >
239 void jacobianAll (
const Point &x, JacobianRangeArray &jacobians )
const 241 assert( jacobians.size() >=
size() );
242 const GeometryType &geo =
geometry();
244 Transformation transformation( geo,
coordinate( x ) );
250 template<
class Po
int,
class DofVector >
251 void hessianAll (
const Point &x,
const DofVector &dofs, HessianRangeType &hessian )
const 253 LocalHessianRangeType localHessian(
typename LocalHessianRangeType::value_type(
RangeFieldType( 0 ) ) );
256 const GeometryType &geo =
geometry();
258 Transformation transformation( geo,
coordinate( x ) );
259 transformation( localHessian, hessian );
263 template<
class Po
int,
class HessianRangeArray >
264 void hessianAll (
const Point &x, HessianRangeArray &hessians )
const 266 assert( hessians.size() >=
size() );
267 const GeometryType &geo =
geometry();
269 Transformation transformation( geo,
coordinate( x ) );
295 const EntityType *entity_;
296 ShapeFunctionSetType shapeFunctionSet_;
303 #endif // #ifndef DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH int order() const
return order of shape functions
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: default.hh:80
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:62
std::size_t size() const
return number of shape functions
Definition: transformation.hh:99
ShapeFunctionSet ShapeFunctionSetType
shape function set type
Definition: default.hh:54
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: default.hh:158
void jacobianAll(const Point &x, JacobianRangeArray &jacobians) const
Definition: default.hh:239
ToNewDimDomainFunctionSpace< LocalFunctionSpaceType, EntityType::Geometry::coorddimension >::Type FunctionSpaceType
type of function space
Definition: default.hh:70
A vector valued function space.
Definition: functionspace.hh:16
const ShapeFunctionSetType & shapeFunctionSet() const
return shape function set
Definition: default.hh:289
void evaluateAll(const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges) const
Definition: default.hh:184
FunctionSpaceType::RangeType RangeType
range type
Definition: default.hh:76
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:74
FieldVector< FieldMatrix< RangeFieldType, dimDomain, dimDomain >, dimRange > HessianRangeType
Intrinsic type used for the hessian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:78
GeometryType::ctype ctype
Definition: default.hh:65
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: default.hh:133
Definition: shapefunctionset/shapefunctionset.hh:33
void jacobianAll(const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians) const
Definition: default.hh:214
LocalFunctionSpaceType::HessianRangeType LocalHessianRangeType
Definition: default.hh:59
EntityType::Geometry GeometryType
Definition: default.hh:63
void evaluateEach(const Point &x, Functor functor) const
evalute each shape function
Definition: transformation.hh:34
Definition: space/basisfunctionset/functor.hh:87
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: default.hh:116
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: default.hh:175
Definition: misc/functor.hh:30
int order() const
return order of basis function set
Definition: default.hh:101
void evaluateAll(const Point &x, const DofVector &dofs, RangeType &value) const
Definition: default.hh:196
FunctionSpaceType::DomainType DomainType
domain type
Definition: default.hh:73
Definition: coordinate.hh:4
Definition: space/basisfunctionset/functor.hh:111
DefaultBasisFunctionSet(const EntityType &entity, const ShapeFunctionSet &shapeFunctionSet=ShapeFunctionSet())
constructor
Definition: default.hh:91
DefaultBasisFunctionSet()
constructor
Definition: default.hh:86
void hessianAll(const Point &x, const DofVector &dofs, HessianRangeType &hessian) const
Definition: default.hh:251
GeometryType geometry() const
Definition: default.hh:292
const Entity & entity() const
return entity
Definition: default.hh:275
void hessianEach(const Point &x, Functor functor) const
evalute hessian of each shape function
void hessianAll(const Point &x, HessianRangeArray &hessians) const
Definition: default.hh:264
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: default.hh:148
Entity EntityType
entity type
Definition: default.hh:52
void jacobianEach(const Point &x, Functor functor) const
evalute jacobian of each shape function
ShapeFunctionSetType::FunctionSpaceType LocalFunctionSpaceType
Definition: default.hh:57
void evaluateAll(const Point &x, RangeArray &values) const
Definition: default.hh:205
Dune::GeometryType type() const
return geometry type
Definition: default.hh:282
LocalFunctionSpaceType::RangeFieldType RangeFieldType
Definition: default.hh:61
const ReferenceElementType & referenceElement() const
return reference element
Definition: default.hh:107
LocalFunctionSpaceType::JacobianRangeType LocalJacobianRangeType
Definition: default.hh:58
void jacobianAll(const Point &x, const DofVector &dofs, JacobianRangeType &jacobian) const
Definition: default.hh:226
Dune::ReferenceElement< ctype, GeometryType::coorddimension > ReferenceElementType
type of reference element
Definition: default.hh:83
convert functions space to space with new dim domain
Definition: functionspace.hh:242
std::size_t size() const
return size of basis function set
Definition: default.hh:104
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:11
implementation of a basis function set for given entity
Definition: default.hh:46
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: default.hh:78