dune-fem  2.4.1-rc
default.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH
2 #define DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH
3 
4 // C++ includes
5 #include <cassert>
6 #include <cstddef>
7 
8 // dune-common includes
9 #include <dune/common/nullptr.hh>
10 
11 // dune-geometry includes
12 #include <dune/geometry/referenceelements.hh>
13 #include <dune/geometry/type.hh>
14 
15 // dune-fem includes
20 #include <dune/fem/version.hh>
21 
22 
23 namespace Dune
24 {
25 
26  namespace Fem
27  {
28 
29  // DefaultBasisFunctionSet
30  // -----------------------
31 
45  template< class Entity, class ShapeFunctionSet >
47  {
49 
50  public:
52  typedef Entity EntityType;
55 
56  protected:
60 
62 
63  typedef typename EntityType::Geometry GeometryType;
64 
65  typedef typename GeometryType::ctype ctype;
66 
67  public:
68  // slight misuse of struct ToLocalFunctionSpace!!!
71 
73  typedef typename FunctionSpaceType::DomainType DomainType;
74 
76  typedef typename FunctionSpaceType::RangeType RangeType;
78  typedef typename FunctionSpaceType::JacobianRangeType JacobianRangeType;
80  typedef typename FunctionSpaceType::HessianRangeType HessianRangeType;
81 
83  typedef Dune::ReferenceElement< ctype, GeometryType::coorddimension > ReferenceElementType;
84 
87  : entity_( nullptr )
88  {}
89 
92  : entity_( &entity ),
93  shapeFunctionSet_( shapeFunctionSet )
94  {}
95 
96 
97  // Basis Function Set Interface Methods
98  // ------------------------------------
99 
101  int order () const { return shapeFunctionSet().order(); }
102 
104  std::size_t size () const { return shapeFunctionSet().size(); }
105 
107  const ReferenceElementType &referenceElement () const
108  {
109  return Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( type() );
110  }
111 
115  template< class QuadratureType, class Vector, class DofVector >
116  void axpy ( const QuadratureType &quad, const Vector &values, DofVector &dofs ) const
117  {
118  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
119  const unsigned int nop = quad.nop();
120  for( unsigned int qp = 0; qp < nop; ++qp )
121  {
122  axpy( quad[ qp ], values[ qp ], dofs );
123  }
124  }
125 
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
134  {
135  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
136  const unsigned int nop = quad.nop();
137  for( unsigned int qp = 0; qp < nop; ++qp )
138  {
139  axpy( quad[ qp ], valuesA[ qp ], dofs );
140  axpy( quad[ qp ], valuesB[ qp ], dofs );
141  }
142  }
143 
147  template< class Point, class DofVector >
148  void axpy ( const Point &x, const RangeType &valueFactor, DofVector &dofs ) const
149  {
150  FunctionalAxpyFunctor< RangeType, DofVector > f( valueFactor, dofs );
151  shapeFunctionSet().evaluateEach( x, f );
152  }
153 
157  template< class Point, class DofVector >
158  void axpy ( const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs ) const
159  {
160  typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
161  const GeometryType &geo = geometry();
162  const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed( coordinate( x ) );
163  LocalJacobianRangeType tmpJacobianFactor( RangeFieldType(0) );
164  for( int r = 0; r < FunctionSpaceType::dimRange; ++r )
165  gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
166 
168  shapeFunctionSet().jacobianEach( x, f );
169  }
170 
174  template< class Point, class DofVector >
175  void axpy ( const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor,
176  DofVector &dofs ) const
177  {
178  axpy( x, valueFactor, dofs );
179  axpy( x, jacobianFactor, dofs );
180  }
181 
183  template< class QuadratureType, class DofVector, class RangeArray >
184  void evaluateAll ( const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges ) const
185  {
186  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
187  const unsigned int nop = quad.nop();
188  for( unsigned int qp = 0; qp < nop; ++qp )
189  {
190  evaluateAll( quad[ qp ], dofs, ranges[ qp ] );
191  }
192  }
193 
195  template< class Point, class DofVector >
196  void evaluateAll ( const Point &x, const DofVector &dofs, RangeType &value ) const
197  {
198  value = RangeType( 0 );
199  AxpyFunctor< DofVector, RangeType > f( dofs, value );
200  shapeFunctionSet().evaluateEach( x, f );
201  }
202 
204  template< class Point, class RangeArray >
205  void evaluateAll ( const Point &x, RangeArray &values ) const
206  {
207  assert( values.size() >= size() );
208  AssignFunctor< RangeArray > f( values );
209  shapeFunctionSet().evaluateEach( x, f );
210  }
211 
213  template< class QuadratureType, class DofVector, class JacobianArray >
214  void jacobianAll ( const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians ) const
215  {
216  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
217  const unsigned int nop = quad.nop();
218  for( unsigned int qp = 0; qp < nop; ++qp )
219  {
220  jacobianAll( quad[ qp ], dofs, jacobians[ qp ] );
221  }
222  }
223 
225  template< class Point, class DofVector >
226  void jacobianAll ( const Point &x, const DofVector &dofs, JacobianRangeType &jacobian ) const
227  {
228  LocalJacobianRangeType localJacobian( RangeFieldType( 0 ) );
229  AxpyFunctor< DofVector, LocalJacobianRangeType > f( dofs, localJacobian );
230  shapeFunctionSet().jacobianEach( x, f );
231  const GeometryType &geo = geometry();
232  typedef JacobianTransformation< GeometryType > Transformation;
233  Transformation transformation( geo, coordinate( x ) );
234  transformation( localJacobian, jacobian );
235  }
236 
238  template< class Point, class JacobianRangeArray >
239  void jacobianAll ( const Point &x, JacobianRangeArray &jacobians ) const
240  {
241  assert( jacobians.size() >= size() );
242  const GeometryType &geo = geometry();
243  typedef JacobianTransformation< GeometryType > Transformation;
244  Transformation transformation( geo, coordinate( x ) );
245  AssignFunctor< JacobianRangeArray, Transformation > f( jacobians, transformation );
246  shapeFunctionSet().jacobianEach( x, f );
247  }
248 
250  template< class Point, class DofVector >
251  void hessianAll ( const Point &x, const DofVector &dofs, HessianRangeType &hessian ) const
252  {
253  LocalHessianRangeType localHessian( typename LocalHessianRangeType::value_type( RangeFieldType( 0 ) ) );
254  AxpyFunctor< DofVector, LocalHessianRangeType > f( dofs, localHessian );
255  shapeFunctionSet().hessianEach( x, f );
256  const GeometryType &geo = geometry();
257  typedef HessianTransformation< GeometryType > Transformation;
258  Transformation transformation( geo, coordinate( x ) );
259  transformation( localHessian, hessian );
260  }
261 
263  template< class Point, class HessianRangeArray >
264  void hessianAll ( const Point &x, HessianRangeArray &hessians ) const
265  {
266  assert( hessians.size() >= size() );
267  const GeometryType &geo = geometry();
268  typedef HessianTransformation< GeometryType > Transformation;
269  Transformation transformation( geo, coordinate( x ) );
270  AssignFunctor< HessianRangeArray, Transformation > f( hessians, transformation );
271  shapeFunctionSet().hessianEach( x, f );
272  }
273 
275  const Entity &entity () const
276  {
277  assert( entity_ );
278  return *entity_;
279  }
280 
282  Dune::GeometryType type () const { return entity().type(); }
283 
284 
285  // Non-interface methods
286  // ---------------------
287 
289  const ShapeFunctionSetType &shapeFunctionSet () const { return shapeFunctionSet_; }
290 
291  protected:
292  GeometryType geometry () const { return entity().geometry(); }
293 
294  private:
295  const EntityType *entity_;
296  ShapeFunctionSetType shapeFunctionSet_;
297  };
298 
299  } // namespace Fem
300 
301 } // namespace Dune
302 
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