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>
50 template<
class Entity,
class ShapeFunctionSet >
72 typedef typename EntityType::Geometry
Geometry ;
74 typedef typename Geometry::ctype
ctype;
81 typedef typename FunctionSpaceType::RangeType
RangeType;
83 typedef typename FunctionSpaceType::DomainType
DomainType;
95 typedef std::decay_t<
decltype( Dune::ReferenceElements< ctype, Geometry::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) ) >
ReferenceElementType;
97 static const int pointSetId = detail::SelectPointSetId< ShapeFunctionSetType >::value;
148 -> decltype(
Dune::ReferenceElements<
ctype,
Geometry::coorddimension >::general(
std::declval< const
Dune::GeometryType & >() ) )
150 return Dune::ReferenceElements< ctype, Geometry::coorddimension >::general(
type() );
156 template<
class QuadratureType,
class Vector,
class DofVector >
157 void axpy (
const QuadratureType &quad,
const Vector &values, DofVector &
dofs )
const
168 template<
class QuadratureType,
class VectorA,
class VectorB,
class DofVector >
169 void axpy (
const QuadratureType &quad,
const VectorA &valuesA,
const VectorB &valuesB, DofVector &
dofs )
const
171 assert( valuesA.size() > 0 );
172 assert( valuesB.size() > 0 );
181 template<
class Po
int,
class DofVector >
191 template<
class Po
int,
class DofVector >
194 typedef typename Geometry::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
196 const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed(
coordinate( x ) );
198 for(
int r = 0; r < FunctionSpaceType::dimRange; ++r )
199 gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
208 template<
class Po
int,
class DofVector >
210 DofVector &
dofs )
const
218 template<
class Po
int,
class DofVector >
221 typedef typename Geometry::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
223 const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed(
coordinate( x ) );
228 Dune::FieldMatrix<double,DomainType::dimension,DomainType::dimension>
231 for(
int r = 0; r < FunctionSpaceType::dimRange; ++r )
232 for(
int j = 0; j < gjit.cols; ++j )
233 for(
int k = 0; k < gjit.cols; ++k )
235 hessianFactor[r].mv(G[k],Hg);
236 tmpHessianFactor[r][j][k] += Hg * G[j];
244 template<
class QuadratureType,
class DofVector,
class RangeArray >
245 void evaluateAll (
const QuadratureType &quad,
const DofVector &
dofs, RangeArray &ranges )
const
247 assert( ranges.size() >= quad.nop() );
252 typedef Codegen :: EvaluateCallerInterfaceTraits< QuadratureType, RangeArray, DofVector > Traits;
253 typedef Codegen :: EvaluateCallerInterface< Traits > BaseEvaluationType;
256 const auto& baseEval = BaseEvaluationType::storage( *
this,
rangeCache( quad ), quad );
261 baseEval->evaluateRanges( quad,
dofs, ranges );
268 const unsigned int nop = quad.nop();
269 for(
unsigned int qp = 0; qp < nop; ++qp )
277 template<
class Po
int,
class DofVector >
286 template<
class Po
int,
class RangeArray >
289 assert( values.size() >=
size() );
290 std::fill( values.begin(), values.end(),
RangeType(0) );
296 template<
class QuadratureType,
class DofVector,
class JacobianArray >
297 void jacobianAll (
const QuadratureType &quad,
const DofVector &
dofs, JacobianArray &jacobians )
const
299 assert( jacobians.size() >= quad.nop() );
304 typedef Codegen :: EvaluateCallerInterfaceTraits< QuadratureType, JacobianArray, DofVector, Geometry > Traits;
305 typedef Codegen :: EvaluateCallerInterface< Traits > BaseEvaluationType;
308 const auto& baseEval = BaseEvaluationType::storage( *
this,
jacobianCache( quad ), quad );
315 baseEval->evaluateJacobians( quad, geo,
dofs, jacobians );
322 const unsigned int nop = quad.nop();
323 for(
unsigned int qp = 0; qp < nop; ++qp )
331 template<
class Po
int,
class DofVector >
340 Transformation transformation( geo,
coordinate( x ) );
341 transformation( localJacobian, jacobian );
345 template<
class Po
int,
class JacobianRangeArray >
346 void jacobianAll (
const Point &x, JacobianRangeArray &jacobians )
const
348 assert( jacobians.size() >=
size() );
352 Transformation transformation( geo,
coordinate( x ) );
358 template<
class QuadratureType,
class DofVector,
class HessianArray >
359 void hessianAll (
const QuadratureType &quad,
const DofVector &
dofs, HessianArray &hessians )
const
361 assert( hessians.size() >= quad.nop() );
363 const unsigned int nop = quad.nop();
364 for(
unsigned int qp = 0; qp < nop; ++qp )
371 template<
class Po
int,
class DofVector >
380 Transformation transformation( geo,
coordinate( x ) );
381 transformation( localHessian, hessian );
385 template<
class Po
int,
class HessianRangeArray >
386 void hessianAll (
const Point &x, HessianRangeArray &hessians )
const
388 assert( hessians.size() >=
size() );
391 Transformation transformation( geo,
coordinate( x ) );
417 template<
class QuadratureType,
class RangeArray,
class DofVector >
418 void axpyImpl (
const QuadratureType &quad,
const RangeArray &rangeFactors, DofVector &
dofs,
const RangeType& )
const
420 assert( rangeFactors.size() >= quad.nop() );
425 typedef Codegen :: EvaluateCallerInterfaceTraits< QuadratureType, RangeArray, DofVector > Traits;
426 typedef Codegen :: EvaluateCallerInterface< Traits > BaseEvaluationType;
429 const auto& baseEval = BaseEvaluationType::storage( *
this,
rangeCache( quad ), quad );
434 baseEval->axpyRanges( quad, rangeFactors,
dofs );
441 const unsigned int nop = quad.nop();
442 for(
unsigned int qp = 0; qp < nop; ++qp )
444 axpy( quad[ qp ], rangeFactors[ qp ],
dofs );
450 template<
class QuadratureType,
class JacobianArray,
class DofVector >
453 assert( jacobianFactors.size() >= quad.nop() );
457 typedef Codegen :: EvaluateCallerInterfaceTraits< QuadratureType, JacobianArray, DofVector, Geometry > Traits;
458 typedef Codegen :: EvaluateCallerInterface< Traits > BaseEvaluationType;
461 const auto& baseEval = BaseEvaluationType::storage( *
this,
jacobianCache( quad ), quad );
468 baseEval->axpyJacobians( quad, geo, jacobianFactors,
dofs );
475 const unsigned int nop = quad.nop();
476 for(
unsigned int qp = 0; qp < nop; ++qp )
478 axpy( quad[ qp ], jacobianFactors[ qp ],
dofs );
483 template <
class QuadratureType>
486 return shapeFunctionSet().scalarShapeFunctionSet().impl().rangeCache( quad );
489 template <
class QuadratureType>
492 return shapeFunctionSet().scalarShapeFunctionSet().impl().jacobianCache( quad );
Definition: bindguard.hh:11
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/default.hh:52
EntityType::Geometry Geometry
Definition: space/basisfunctionset/default.hh:72
FunctionSpaceType::DomainType DomainType
domain type
Definition: space/basisfunctionset/default.hh:83
void evaluateAll(const Point &x, const DofVector &dofs, RangeType &value) const
Definition: space/basisfunctionset/default.hh:278
LocalFunctionSpaceType::HessianRangeType LocalHessianRangeType
Definition: space/basisfunctionset/default.hh:68
DefaultBasisFunctionSet()
constructor
Definition: space/basisfunctionset/default.hh:100
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: space/basisfunctionset/default.hh:87
const EntityType * entity_
Definition: space/basisfunctionset/default.hh:498
ToNewDimDomainFunctionSpace< LocalFunctionSpaceType, Geometry::coorddimension >::Type FunctionSpaceType
type of function space
Definition: space/basisfunctionset/default.hh:78
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: space/basisfunctionset/default.hh:169
const ShapeFunctionSetType & shapeFunctionSet() const
return shape function set
Definition: space/basisfunctionset/default.hh:413
auto referenceElement() const -> decltype(Dune::ReferenceElements< ctype, Geometry::coorddimension >::general(std::declval< const Dune::GeometryType & >()))
return reference element
Definition: space/basisfunctionset/default.hh:147
void axpy(const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor, DofVector &dofs) const
evaluate all basis function and derivatives and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:209
static constexpr bool codegenShapeFunctionSet
Definition: space/basisfunctionset/default.hh:63
std::optional< Geometry > geometry_
Definition: space/basisfunctionset/default.hh:500
FunctionSpaceType::RangeType RangeType
range type
Definition: space/basisfunctionset/default.hh:81
ShapeFunctionSetType::FunctionSpaceType LocalFunctionSpaceType
Definition: space/basisfunctionset/default.hh:66
LocalFunctionSpaceType::RangeFieldType RangeFieldType
Definition: space/basisfunctionset/default.hh:70
LocalFunctionSpaceType::JacobianRangeType LocalJacobianRangeType
Definition: space/basisfunctionset/default.hh:67
void axpy(const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs) const
evaluate all derivatives of all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:192
void jacobianAll(const Point &x, JacobianRangeArray &jacobians) const
Definition: space/basisfunctionset/default.hh:346
bool valid() const
return true if entity pointer is set
Definition: space/basisfunctionset/default.hh:404
DefaultBasisFunctionSet(const EntityType &entity, const ShapeFunctionSet &shapeFunctionSet=ShapeFunctionSet())
constructor
Definition: space/basisfunctionset/default.hh:103
const auto & jacobianCache(const QuadratureType &quad) const
Definition: space/basisfunctionset/default.hh:490
void axpy(const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs) const
Add H:D^2phi to each dof.
Definition: space/basisfunctionset/default.hh:219
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: space/basisfunctionset/default.hh:85
void hessianAll(const Point &x, HessianRangeArray &hessians) const
Definition: space/basisfunctionset/default.hh:386
void evaluateAll(const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges) const
Definition: space/basisfunctionset/default.hh:245
ShapeFunctionSetType shapeFunctionSet_
Definition: space/basisfunctionset/default.hh:499
ScalarFunctionSpaceType::RangeType ScalarRangeType
Definition: space/basisfunctionset/default.hh:91
std::decay_t< decltype(Dune::ReferenceElements< ctype, Geometry::coorddimension >::general(std::declval< const Dune::GeometryType & >())) > ReferenceElementType
type of reference element
Definition: space/basisfunctionset/default.hh:95
DefaultBasisFunctionSet & operator=(const DefaultBasisFunctionSet &other)
Definition: space/basisfunctionset/default.hh:124
void hessianAll(const Point &x, const DofVector &dofs, HessianRangeType &hessian) const
Definition: space/basisfunctionset/default.hh:372
void hessianAll(const QuadratureType &quad, const DofVector &dofs, HessianArray &hessians) const
Definition: space/basisfunctionset/default.hh:359
DefaultBasisFunctionSet(const DefaultBasisFunctionSet &other)
Definition: space/basisfunctionset/default.hh:113
void evaluateAll(const Point &x, RangeArray &values) const
Definition: space/basisfunctionset/default.hh:287
ShapeFunctionSet ShapeFunctionSetType
shape function set type
Definition: space/basisfunctionset/default.hh:59
ScalarFunctionSpaceType::JacobianRangeType ScalarJacobianRangeType
Definition: space/basisfunctionset/default.hh:92
const Entity & entity() const
return entity
Definition: space/basisfunctionset/default.hh:397
static const int pointSetId
Definition: space/basisfunctionset/default.hh:97
std::size_t size() const
return size of basis function set
Definition: space/basisfunctionset/default.hh:144
void jacobianAll(const Point &x, const DofVector &dofs, JacobianRangeType &jacobian) const
Definition: space/basisfunctionset/default.hh:332
Geometry geometry() const
Definition: space/basisfunctionset/default.hh:496
const auto & rangeCache(const QuadratureType &quad) const
Definition: space/basisfunctionset/default.hh:484
void axpyImpl(const QuadratureType &quad, const RangeArray &rangeFactors, DofVector &dofs, const RangeType &) const
evaluate all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:418
Entity EntityType
entity type
Definition: space/basisfunctionset/default.hh:57
void axpyImpl(const QuadratureType &quad, const JacobianArray &jacobianFactors, DofVector &dofs, const JacobianRangeType &) const
evaluate all basis function and multiply with given values and add to dofs
Definition: space/basisfunctionset/default.hh:451
Dune::GeometryType type() const
return geometry type
Definition: space/basisfunctionset/default.hh:407
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: space/basisfunctionset/default.hh:182
Geometry::ctype ctype
Definition: space/basisfunctionset/default.hh:74
FunctionSpaceType::ScalarFunctionSpaceType ScalarFunctionSpaceType
Definition: space/basisfunctionset/default.hh:89
void jacobianAll(const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians) const
Definition: space/basisfunctionset/default.hh:297
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: space/basisfunctionset/default.hh:157
int order() const
return order of basis function set
Definition: space/basisfunctionset/default.hh:141
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