dune-fem  2.4.1-rc
localfunction.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_FUNCTION_LOCALFUNCTION_LOCALFUNCTION_HH
2 #define DUNE_FEM_FUNCTION_LOCALFUNCTION_LOCALFUNCTION_HH
3 
4 #include <utility>
5 
6 // dune-common includes
7 #include <dune/common/fvector.hh>
8 
9 
10 namespace Dune
11 {
12 
13  namespace Fem
14  {
15 
30  // LocalFunction
31  // -------------
32 
40  template< class BasisFunctionSet, class LocalDofVector >
42  {
44  public:
45 
48 
50  typedef LocalDofVector LocalDofVectorType;
51 
53  typedef typename LocalDofVectorType::value_type DofType;
54 
56  typedef typename LocalDofVectorType::size_type SizeType;
57 
60 
63 
76 
78  typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
79 
84 
87 
89  explicit LocalFunction ( const BasisFunctionSetType &basisFunctionSet )
90  : basisFunctionSet_( basisFunctionSet )
91  {
92  localDofVector_.resize( basisFunctionSet.size() );
93  }
94 
96  explicit LocalFunction ( const LocalDofVectorType &localDofVector ) : localDofVector_( localDofVector ) {}
97 
99  LocalFunction ( const BasisFunctionSetType &basisFunctionSet, const LocalDofVector &localDofVector )
100  : basisFunctionSet_( basisFunctionSet ),
101  localDofVector_( localDofVector )
102  {
103  localDofVector_.resize( basisFunctionSet.size() );
104  }
105 
107  explicit LocalFunction ( LocalDofVectorType &&localDofVector ) : localDofVector_( localDofVector ) {}
108 
110  LocalFunction ( const BasisFunctionSetType &basisFunctionSet, LocalDofVector &&localDofVector )
111  : basisFunctionSet_( basisFunctionSet ),
113  {
114  localDofVector_.resize( basisFunctionSet.size() );
115  }
116 
118  LocalFunction ( ThisType && other )
121  {}
122 
124  LocalFunction ( const ThisType & other )
127  {}
128 
134  const DofType &operator[] ( SizeType num ) const { return localDofVector_[ num ]; }
135 
141  DofType &operator[] ( SizeType num ) { return localDofVector_[ num ]; }
142 
152  template< class T >
154  {
155  localDofVector() += other.localDofVector();
156  return *this;
157  }
158 
163  template< class T >
165  {
166  localDofVector() = other.localDofVector();
167  }
168 
170  void clear ()
171  {
172  std::fill( localDofVector().begin(), localDofVector().end(), DofType( 0 ) );
173  }
174 
184  template< class T >
186  {
187  localDofVector() -= other.localDofVector();
188  return *this;
189  }
190 
201  template< class T >
202  ThisType &axpy ( const RangeFieldType s, const LocalFunction< BasisFunctionSet, T > &other )
203  {
204  localDofVector().axpy( s, other.localDofVector() );
205  return *this;
206  }
207 
220  template< class PointType >
221  void axpy ( const PointType &x, const RangeType &factor )
222  {
223  basisFunctionSet().axpy( x, factor, localDofVector() );
224  }
225 
238  template< class PointType >
239  void axpy ( const PointType &x, const JacobianRangeType &factor)
240  {
241  basisFunctionSet().axpy( x, factor, localDofVector() );
242  }
243 
257  template< class PointType >
258  void axpy ( const PointType &x, const RangeType &factor1, const JacobianRangeType &factor2 )
259  {
260  basisFunctionSet().axpy( x, factor1, factor2, localDofVector() );
261  }
262 
273  int order () const { return basisFunctionSet().order(); }
274 
279  const BasisFunctionSetType &basisFunctionSet () const { return basisFunctionSet_; }
280 
285  const EntityType &entity () const { return basisFunctionSet().entity(); }
286 
288  void init ( const BasisFunctionSetType &basisFunctionSet )
289  {
291  localDofVector_.resize( basisFunctionSet.size() );
292  }
293 
299  template< class PointType >
300  void evaluate ( const PointType &x, RangeType &ret ) const
301  {
303  }
304 
313  template< class PointType >
314  void jacobian ( const PointType &x, JacobianRangeType &ret ) const
315  {
317  }
318 
327  template< class PointType >
328  void hessian ( const PointType &x, HessianRangeType &ret ) const
329  {
331  }
332 
340  int numDofs () const { return localDofVector().size(); }
341 
342  SizeType size () const { return localDofVector().size(); }
343 
352  int numScalarDofs () const
353  {
354  assert( numDofs() % dimRange == 0 );
355  return numDofs() / dimRange;
356  }
357 
360  template< class QuadratureType, class VectorType >
361  void axpyQuadrature ( const QuadratureType &quad, const VectorType &values )
362  {
363  basisFunctionSet().axpy( quad, values, localDofVector() );
364  }
365 
368  template< class QuadratureType, class RangeVectorType, class JacobianRangeVectorType >
369  void axpyQuadrature ( const QuadratureType &quad,
370  const RangeVectorType& rangeVector,
371  const JacobianRangeVectorType& jacobianVector )
372  {
373  basisFunctionSet().axpy( quad, rangeVector, jacobianVector, localDofVector() );
374  }
375 
377  template< class QuadratureType, class VectorType >
378  void evaluateQuadrature( const QuadratureType &quad, VectorType &result ) const
379  {
380  assert( result.size() > 0 );
381  evaluateQuadrature( quad, result, result[ 0 ] );
382  }
383 
384 
386  const LocalDofVectorType &localDofVector () const { return localDofVector_; }
387 
389  LocalDofVectorType &localDofVector () { return localDofVector_; }
390 
391  protected:
392  // evaluate local function and store results in vector of RangeTypes
393  // this method only helps to identify the correct method on
394  // the basis function set
395  template< class QuadratureType, class VectorType >
396  void evaluateQuadrature( const QuadratureType &quad, VectorType &result, const RangeType & ) const
397  {
398  basisFunctionSet().evaluateAll( quad, localDofVector(), result );
399  }
400 
401  // evaluate jacobian of local function and store result in vector of
402  // JacobianRangeTypes, this method only helps to identify the correct method on
403  // the basis function set
404  template< class QuadratureType, class VectorType >
405  void evaluateQuadrature( const QuadratureType &quad, VectorType &result, const JacobianRangeType & ) const
406  {
407  basisFunctionSet().jacobianAll( quad, localDofVector(), result );
408  }
409 
410  BasisFunctionSetType basisFunctionSet_;
411  LocalDofVectorType localDofVector_;
412  };
413 
416  } // end namespace Fem
417 
418 } // end namespace Dune
419 
420 #endif // #ifndef DUNE_FEM_FUNCTION_LOCALFUNCTION_LOCALFUNCTION_HH
LocalDofVectorType localDofVector_
Definition: localfunction.hh:411
int numDofs() const
obtain the number of local DoFs
Definition: localfunction.hh:340
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:62
LocalDofVectorType::value_type DofType
type of DoF use with the discrete function
Definition: localfunction.hh:53
BasisFunctionSetType basisFunctionSet_
Definition: localfunction.hh:410
LocalFunction(const BasisFunctionSetType &basisFunctionSet)
ctor taking a basisFunctionSet, calling default ctor for LocalDofVectorType, and resize ...
Definition: localfunction.hh:89
FunctionSpaceType::DomainFieldType DomainFieldType
field type of the domain
Definition: localfunction.hh:65
LocalFunction()
default constructor, calls default ctor of BasisFunctionSetType and LocalDofVectorType ...
Definition: localfunction.hh:86
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
FunctionSpaceType::RangeFieldType RangeFieldType
field type of the range
Definition: localfunction.hh:67
void jacobianAll(const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians) const
evaluate the jacobian of all basis functions and store the result in the jacobians array ...
dimension of domain vector space
Definition: functionspaceinterface.hh:45
LocalDofVectorType & localDofVector()
return mutable reference to local Dof Vector
Definition: localfunction.hh:389
void axpy(const Quadrature &quad, const Vector &values, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
dimension of range vector space
Definition: functionspaceinterface.hh:47
void evaluateQuadrature(const QuadratureType &quad, VectorType &result, const JacobianRangeType &) const
Definition: localfunction.hh:405
LocalFunction(const LocalDofVectorType &localDofVector)
ctor taking a localDofVector, calling default ctor for BasisFunctionSetType
Definition: localfunction.hh:96
static const int dimRange
dimension of the range
Definition: localfunction.hh:83
const DofType & operator[](SizeType num) const
access to local dofs (read-only)
Definition: localfunction.hh:134
const BasisFunctionSetType & basisFunctionSet() const
obtain the basis function set for this local function
Definition: localfunction.hh:279
FunctionSpaceType::JacobianRangeType JacobianRangeType
type of the Jacobian, i.e., type of evaluated Jacobian matrix
Definition: localfunction.hh:73
interface for local functions
Definition: localfunction.hh:41
LocalFunction(const ThisType &other)
copy constructor
Definition: localfunction.hh:124
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::DomainFieldType DomainFieldType
Intrinsic type used for values in the domain field (usually a double)
Definition: functionspaceinterface.hh:59
BasisFunctionSetType::EntityType EntityType
type of the entity, the local function lives on is given by the space
Definition: localfunction.hh:59
int order() const
return order of basis function set
LocalDofVectorType::size_type SizeType
type of index
Definition: localfunction.hh:56
Definition: coordinate.hh:4
void evaluateAll(const Quadrature &quad, const DofVector &dofs, RangeArray &ranges) const
evaluate all basis functions and store the result in the ranges array
ThisType & operator-=(const LocalFunction< BasisFunctionSet, T > &other)
subtract another local function to this one
Definition: localfunction.hh:185
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:66
EntityType::Geometry::LocalCoordinate LocalCoordinateType
type of local coordinates
Definition: localfunction.hh:78
void axpy(const PointType &x, const RangeType &factor1, const JacobianRangeType &factor2)
axpy operation for local function
Definition: localfunction.hh:258
void clear()
set all DoFs to zero
Definition: localfunction.hh:170
const EntityType & entity() const
obtain the entity, this local function lives on
Definition: localfunction.hh:285
void evaluateQuadrature(const QuadratureType &quad, VectorType &result) const
evaluate all basisfunctions for all quadrature points and store the results in the result vector ...
Definition: localfunction.hh:378
void axpyQuadrature(const QuadratureType &quad, const RangeVectorType &rangeVector, const JacobianRangeVectorType &jacobianVector)
evaluate all basisfunctions for all quadrature points, multiply with the given factor and add the res...
Definition: localfunction.hh:369
static const int dimDomain
dimension of the domain
Definition: localfunction.hh:81
int numScalarDofs() const
obtain the number of local DoFs in the scalar case
Definition: localfunction.hh:352
STL namespace.
LocalFunction(LocalDofVectorType &&localDofVector)
half move ctor
Definition: localfunction.hh:107
void axpy(const PointType &x, const JacobianRangeType &factor)
axpy operation for local function
Definition: localfunction.hh:239
LocalFunction(const BasisFunctionSetType &basisFunctionSet, const LocalDofVector &localDofVector)
copy given agruments
Definition: localfunction.hh:99
LocalDofVector LocalDofVectorType
type of local Dof Vector
Definition: localfunction.hh:50
BasisFunctionSetType::FunctionSpaceType FunctionSpaceType
type of functionspace
Definition: localfunction.hh:62
ThisType & operator+=(const LocalFunction< BasisFunctionSet, T > &other)
add another local function to this one
Definition: localfunction.hh:153
void evaluateQuadrature(const QuadratureType &quad, VectorType &result, const RangeType &) const
Definition: localfunction.hh:396
FunctionSpace< typename Entity::Geometry::ctype, typename Range::value_type, Entity::Geometry::coorddimension, Range::dimension > FunctionSpaceType
function space type
Definition: basisfunctionset/basisfunctionset.hh:40
void init(const BasisFunctionSetType &basisFunctionSet)
Definition: localfunction.hh:288
void move(ArrayInterface< T > &array, const unsigned int oldOffset, const unsigned int newOffset, const unsigned int length)
Definition: array_inline.hh:38
LocalFunction(ThisType &&other)
move constructor
Definition: localfunction.hh:118
void hessianAll(const Point &x, const DofVector &dofs, HessianRangeType &hessian) const
void jacobian(const PointType &x, JacobianRangeType &ret) const
evaluate Jacobian of the local function
Definition: localfunction.hh:314
void assign(const LocalFunction< BasisFunctionSet, T > &other)
assign all DoFs of this local function
Definition: localfunction.hh:164
void axpy(const PointType &x, const RangeType &factor)
axpy operation for local function
Definition: localfunction.hh:221
void evaluate(const PointType &x, RangeType &ret) const
evaluate the local function
Definition: localfunction.hh:300
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:70
Entity EntityType
entity type
Definition: basisfunctionset/basisfunctionset.hh:35
FunctionSpaceType::DomainType DomainType
type of domain vectors, i.e., type of coordinates
Definition: localfunction.hh:69
const EntityType & entity() const
return entity
void hessian(const PointType &x, HessianRangeType &ret) const
evaluate Hessian of the local function
Definition: localfunction.hh:328
void axpyQuadrature(const QuadratureType &quad, const VectorType &values)
evaluate all basisfunctions for all quadrature points, multiply with the given factor and add the res...
Definition: localfunction.hh:361
SizeType size() const
Definition: localfunction.hh:342
int order() const
obtain the order of this local function
Definition: localfunction.hh:273
BasisFunctionSet BasisFunctionSetType
type of basis function set
Definition: localfunction.hh:47
FunctionSpaceType::HessianRangeType HessianRangeType
type of the Hessian
Definition: localfunction.hh:75
Definition: basisfunctionset/basisfunctionset.hh:31
LocalFunction(const BasisFunctionSetType &basisFunctionSet, LocalDofVector &&localDofVector)
half move ctor
Definition: localfunction.hh:110
ThisType & axpy(const RangeFieldType s, const LocalFunction< BasisFunctionSet, T > &other)
add a multiple of another local function to this one
Definition: localfunction.hh:202
const LocalDofVectorType & localDofVector() const
return const reference to local Dof Vector
Definition: localfunction.hh:386
FunctionSpaceType::RangeType RangeType
type of range vectors, i.e., type of function values
Definition: localfunction.hh:71