dune-fem  2.6-git
basisfunctionset/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 #include <type_traits>
9 #include <utility>
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 std::decay_t< decltype( Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) ) > 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  auto referenceElement () const
108  -> decltype( Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) )
109  {
110  return Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( type() );
111  }
112 
116  template< class QuadratureType, class Vector, class DofVector >
117  void axpy ( const QuadratureType &quad, const Vector &values, DofVector &dofs ) const
118  {
119  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
120  const unsigned int nop = quad.nop();
121  for( unsigned int qp = 0; qp < nop; ++qp )
122  {
123  axpy( quad[ qp ], values[ qp ], dofs );
124  }
125  }
126 
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
135  {
136  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
137  const unsigned int nop = quad.nop();
138  for( unsigned int qp = 0; qp < nop; ++qp )
139  {
140  axpy( quad[ qp ], valuesA[ qp ], dofs );
141  axpy( quad[ qp ], valuesB[ qp ], dofs );
142  }
143  }
144 
148  template< class Point, class DofVector >
149  void axpy ( const Point &x, const RangeType &valueFactor, DofVector &dofs ) const
150  {
152  shapeFunctionSet().evaluateEach( x, f );
153  }
154 
158  template< class Point, class DofVector >
159  void axpy ( const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs ) const
160  {
161  typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
162  const GeometryType &geo = geometry();
163  const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed( coordinate( x ) );
164  LocalJacobianRangeType tmpJacobianFactor( RangeFieldType(0) );
165  for( int r = 0; r < FunctionSpaceType::dimRange; ++r )
166  gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
167 
169  shapeFunctionSet().jacobianEach( x, f );
170  }
173  template< class Point, class DofVector >
174  void axpy ( const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs ) const
175  {
176  typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
177  const GeometryType &geo = geometry();
178  const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed( coordinate( x ) );
179  LocalHessianRangeType tmpHessianFactor( RangeFieldType(0) );
180  // don't know how to work directly with the DiagonalMatrix
181  // returned from YaspGrid since gjit[k][l] is not correctly
182  // implemented for k!=l so convert first into a dense matrix...
183  Dune::FieldMatrix<double,DomainType::dimension,DomainType::dimension>
184  G = gjit;
185  DomainType Hg;
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 )
189  {
190  hessianFactor[r].mv(G[k],Hg);
191  tmpHessianFactor[r][j][k] += Hg * G[j];
192  }
194  shapeFunctionSet().hessianEach( x, f );
195  }
196 
200  template< class Point, class DofVector >
201  void axpy ( const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor,
202  DofVector &dofs ) const
203  {
204  axpy( x, valueFactor, dofs );
205  axpy( x, jacobianFactor, dofs );
206  }
207 
209  template< class QuadratureType, class DofVector, class RangeArray >
210  void evaluateAll ( const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges ) const
211  {
212  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
213  const unsigned int nop = quad.nop();
214  for( unsigned int qp = 0; qp < nop; ++qp )
215  {
216  evaluateAll( quad[ qp ], dofs, ranges[ qp ] );
217  }
218  }
219 
221  template< class Point, class DofVector >
222  void evaluateAll ( const Point &x, const DofVector &dofs, RangeType &value ) const
223  {
224  value = RangeType( 0 );
226  shapeFunctionSet().evaluateEach( x, f );
227  }
228 
230  template< class Point, class RangeArray >
231  void evaluateAll ( const Point &x, RangeArray &values ) const
232  {
233  assert( values.size() >= size() );
234  AssignFunctor< RangeArray > f( values );
235  shapeFunctionSet().evaluateEach( x, f );
236  }
237 
239  template< class QuadratureType, class DofVector, class JacobianArray >
240  void jacobianAll ( const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians ) const
241  {
242  // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
243  const unsigned int nop = quad.nop();
244  for( unsigned int qp = 0; qp < nop; ++qp )
245  {
246  jacobianAll( quad[ qp ], dofs, jacobians[ qp ] );
247  }
248  }
249 
251  template< class Point, class DofVector >
252  void jacobianAll ( const Point &x, const DofVector &dofs, JacobianRangeType &jacobian ) const
253  {
254  LocalJacobianRangeType localJacobian( RangeFieldType( 0 ) );
256  shapeFunctionSet().jacobianEach( x, f );
257  const GeometryType &geo = geometry();
258  typedef JacobianTransformation< GeometryType > Transformation;
259  Transformation transformation( geo, coordinate( x ) );
260  transformation( localJacobian, jacobian );
261  }
262 
264  template< class Point, class JacobianRangeArray >
265  void jacobianAll ( const Point &x, JacobianRangeArray &jacobians ) const
266  {
267  assert( jacobians.size() >= size() );
268  const GeometryType &geo = geometry();
269  typedef JacobianTransformation< GeometryType > Transformation;
270  Transformation transformation( geo, coordinate( x ) );
271  AssignFunctor< JacobianRangeArray, Transformation > f( jacobians, transformation );
272  shapeFunctionSet().jacobianEach( x, f );
273  }
274 
276  template< class Point, class DofVector >
277  void hessianAll ( const Point &x, const DofVector &dofs, HessianRangeType &hessian ) const
278  {
279  LocalHessianRangeType localHessian( typename LocalHessianRangeType::value_type( RangeFieldType( 0 ) ) );
281  shapeFunctionSet().hessianEach( x, f );
282  const GeometryType &geo = geometry();
283  typedef HessianTransformation< GeometryType > Transformation;
284  Transformation transformation( geo, coordinate( x ) );
285  transformation( localHessian, hessian );
286  }
287 
289  template< class Point, class HessianRangeArray >
290  void hessianAll ( const Point &x, HessianRangeArray &hessians ) const
291  {
292  assert( hessians.size() >= size() );
293  const GeometryType &geo = geometry();
294  typedef HessianTransformation< GeometryType > Transformation;
295  Transformation transformation( geo, coordinate( x ) );
296  AssignFunctor< HessianRangeArray, Transformation > f( hessians, transformation );
297  shapeFunctionSet().hessianEach( x, f );
298  }
299 
301  const Entity &entity () const
302  {
303  assert( entity_ );
304  return *entity_;
305  }
306 
308  Dune::GeometryType type () const { return entity().type(); }
309 
310 
311  // Non-interface methods
312  // ---------------------
313 
315  const ShapeFunctionSetType &shapeFunctionSet () const { return shapeFunctionSet_; }
316 
317  protected:
318  GeometryType geometry () const { return entity().geometry(); }
319 
320  private:
321  const EntityType *entity_;
322  ShapeFunctionSetType shapeFunctionSet_;
323  };
324 
325  } // namespace Fem
326 
327 } // namespace Dune
328 
329 #endif // #ifndef DUNE_FEM_BASISFUNCTIONSET_DEFAULT_HH
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