DUNE-FEM (unstable)

Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction > Class Template Reference

A constant local function carrying values for one entity. More...

#include <dune/fem/function/localfunction/const.hh>

Public Types

typedef BaseType::SizeType SizeType
 type of SizeType
 
typedef EntityType::Geometry Geometry
 type of the geometry, the local function lives on is given by the space
 
typedef FunctionSpaceType::DomainFieldType DomainFieldType
 field type of the domain
 
typedef FunctionSpaceType::RangeFieldType RangeFieldType
 field type of the range
 
typedef Geometry::LocalCoordinate LocalCoordinateType
 type of local coordinates
 
typedef Traits::derived_type derived_type
 type of derived vector class
 
typedef Traits::value_type value_type
 export the type representing the field
 
typedef FieldTraits< value_type >::field_type field_type
 export the type representing the field
 
typedef Traits::value_type block_type
 export the type representing the components
 
typedef Traits::size_type size_type
 The type used for the index access and size operation.
 
typedef DenseIterator< DenseVector, value_typeIterator
 Iterator class for sequential access.
 
typedef Iterator iterator
 typedef for stl compliant access
 
typedef DenseIterator< const DenseVector, const value_typeConstIterator
 ConstIterator class for sequential access.
 
typedef ConstIterator const_iterator
 typedef for stl compliant access
 

Public Member Functions

 ConstLocalDiscreteFunction (const DiscreteFunctionType &df)
 constructor creating a local function without binding it to an entity More...
 
 ConstLocalDiscreteFunction (const typename DiscreteFunctionType::LocalFunctionType &localFunction)
 cast a MutableLocalFunction into this one !!! expensive !!!
 
 ConstLocalDiscreteFunction (const DiscreteFunctionType &df, const EntityType &entity)
 constructor creating a local function and binding it to an entity More...
 
 ConstLocalDiscreteFunction (const ThisType &other)
 copy constructor
 
 ConstLocalDiscreteFunction (ThisType &&other)
 move constructor
 
template<class Point >
RangeType evaluate (const Point &p) const
 evaluate the local function More...
 
template<class Point >
JacobianRangeType jacobian (const Point &p) const
 evaluate Jacobian of the local function More...
 
template<class Point >
HessianRangeType hessian (const Point &p) const
 evaluate Hessian of the local function More...
 
void init (const EntityType &entity)
 interface for local functions :: init More...
 
const LocalDofVectorType & localDofVector () const
 return const reference to local Dof Vector

 
LocalDofVectorType & localDofVector ()
 return mutable reference to local Dof Vector

 
template<class PointType >
void evaluate (const PointType &x, RangeType &ret) const
 evaluate the local function More...
 
template<class PointType >
void jacobian (const PointType &x, JacobianRangeType &ret) const
 evaluate Jacobian of the local function More...
 
template<class PointType >
void hessian (const PointType &x, HessianRangeType &ret) const
 evaluate Hessian of the local function More...
 
value_typeoperator[] (size_type i)
 random access
 
derived_typeaxpy (const field_type &a, const DenseVector< Other > &x)
 vector space axpy operation ( *this += a x )
 
int order () const
 obtain the order of this local function More...
 
const BasisFunctionSetTypebasisFunctionSet () const
 obtain the basis function set for this local function More...
 
const EntityTypeentity () const
 obtain the entity, this local function lives on More...
 
const Geometrygeometry () const
 obtain the geometry, this local function lives on More...
 
int numDofs () const
 obtain the number of local DoFs More...
 
SizeType size () const
 obtain the number of local DoFs More...
 
template<class QuadratureType , class ... Vectors>
void axpyQuadrature (const QuadratureType &quad, const Vectors &... values)
 evaluate all basisfunctions for all quadrature points, multiply with the given factor and add the result to the local coefficients

 
template<class QuadratureType , class RangeVectorType , class JacobianRangeVectorType >
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 result to the local coefficients
 
template<class QuadratureType , class ... Vectors>
void evaluateQuadrature (const QuadratureType &quad, Vectors &... vec) const
 evaluate all basisfunctions for all quadrature points and store the results in the result vector
 
template<class QuadratureType , class ... Vectors>
void jacobianQuadrature (const QuadratureType &quad, Vectors &... vec) const
 evaluate all Jacobians for all basis functions for all quadrature points and store the results in the result vector
 
template<class QuadratureType , class ... Vectors>
void hessianQuadrature (const QuadratureType &quad, Vectors &... vec) const
 evaluate all hessians of all basis functions for all quadrature points and store the results in the result vector
 
bool valid () const
 Returns true if local function if bind or init was previously called.
 
value_typefront ()
 return reference to first element
 
const value_typefront () const
 return reference to first element
 
value_typeback ()
 return reference to last element
 
const value_typeback () const
 return reference to last element
 
bool empty () const
 checks whether the container is empty
 
Iterator begin ()
 begin iterator
 
ConstIterator begin () const
 begin ConstIterator
 
Iterator end ()
 end iterator
 
ConstIterator end () const
 end ConstIterator
 
Iterator beforeEnd ()
 
ConstIterator beforeEnd () const
 
Iterator beforeBegin ()
 
ConstIterator beforeBegin () const
 
Iterator find (size_type i)
 return iterator to given element or end()
 
ConstIterator find (size_type i) const
 return iterator to given element or end()
 
derived_typeoperator+= (const DenseVector< Other > &x)
 vector space addition
 
std::enable_if< std::is_convertible< ValueType, value_type >::value, derived_type >::type & operator+= (const ValueType &kk)
 vector space add scalar to all comps More...
 
derived_typeoperator-= (const DenseVector< Other > &x)
 vector space subtraction
 
std::enable_if< std::is_convertible< ValueType, value_type >::value, derived_type >::type & operator-= (const ValueType &kk)
 vector space subtract scalar from all comps More...
 
derived_type operator+ (const DenseVector< Other > &b) const
 Binary vector addition.
 
derived_type operator- (const DenseVector< Other > &b) const
 Binary vector subtraction.
 
derived_type operator- () const
 Vector negation.
 
std::enable_if< std::is_convertible< FieldType, field_type >::value, derived_type >::type & operator*= (const FieldType &kk)
 vector space multiplication with scalar More...
 
std::enable_if< std::is_convertible< FieldType, field_type >::value, derived_type >::type & operator/= (const FieldType &kk)
 vector space division by scalar More...
 
bool operator== (const DenseVector< Other > &x) const
 Binary vector comparison.
 
bool operator!= (const DenseVector< Other > &x) const
 Binary vector incomparison.
 
PromotionTraits< field_type, typenameDenseVector< Other >::field_type >::PromotedType operator* (const DenseVector< Other > &x) const
 indefinite vector dot product \(\left (x^T \cdot y \right)\) which corresponds to Petsc's VecTDot More...
 
PromotionTraits< field_type, typenameDenseVector< Other >::field_type >::PromotedType dot (const DenseVector< Other > &x) const
 vector dot product \(\left (x^H \cdot y \right)\) which corresponds to Petsc's VecDot More...
 
FieldTraits< value_type >::real_type one_norm () const
 one norm (sum over absolute values of entries)
 
FieldTraits< value_type >::real_type one_norm_real () const
 simplified one norm (uses Manhattan norm for complex values)
 
FieldTraits< value_type >::real_type two_norm () const
 two norm sqrt(sum over squared values of entries)
 
FieldTraits< value_type >::real_type two_norm2 () const
 square of two norm (sum over squared values of entries), need for block recursion
 
FieldTraits< vt >::real_type infinity_norm () const
 infinity norm (maximum of absolute values of entries)
 
FieldTraits< vt >::real_type infinity_norm () const
 infinity norm (maximum of absolute values of entries)
 
FieldTraits< vt >::real_type infinity_norm_real () const
 simplified infinity norm (uses Manhattan norm for complex values)
 
FieldTraits< vt >::real_type infinity_norm_real () const
 simplified infinity norm (uses Manhattan norm for complex values)
 
size_type N () const
 number of blocks in the vector (are of size 1 here)
 
size_type dim () const
 dimension of the vector space
 

Static Public Attributes

static const int dimDomain = FunctionSpaceType::dimDomain
 dimension of the domain
 
static const int dimRange = FunctionSpaceType::dimRange
 dimension of the range
 
static constexpr int blocklevel
 The number of block levels we contain. This is the leaf, that is, 1.
 

Protected Member Functions

void clear ()
 set all DoFs to zero
 
template<class T >
void assign (const LocalFunction< BasisFunctionSet, T > &other)
 assign all DoFs of this local function More...
 
template<class PointType >
void axpy (const PointType &x, const RangeType &factor)
 axpy operation for local function More...
 
template<class PointType >
void axpy (const PointType &x, const JacobianRangeType &factor)
 axpy operation for local function More...
 
template<class PointType >
void axpy (const PointType &x, const RangeType &factor1, const JacobianRangeType &factor2)
 axpy operation for local function More...
 
void init (const BasisFunctionSetType &basisFunctionSet)
 initialize the local function for an basisFunctionSet More...
 

Related Functions

(Note that these are not member functions.)

std::ostream & operator<< (std::ostream &s, const DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > > &v)
 Write a DenseVector to an output stream. More...
 

Detailed Description

template<class DiscreteFunction>
class Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >

A constant local function carrying values for one entity.

A ConstLocalDiscreteFunction is a LocalFunction which is basically doing the same as the LocalFunction of a discrete function. The difference is that the local dofs are not kept as references but are copied to a local storage. Therefore, this is a const local function and any modification of dofs is not allowed.

Note
Local DoF numbers correspond directly to array indices. Hence it may be more cache efficient to generate a ConstLocalFunction when only a const access to the local function is needed.
Parameters
DiscreteFunctiontype of the discrete function, the local function shall belong to

Constructor & Destructor Documentation

◆ ConstLocalDiscreteFunction() [1/2]

template<class DiscreteFunction >
Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::ConstLocalDiscreteFunction ( const DiscreteFunctionType &  df)
inlineexplicit

constructor creating a local function without binding it to an entity

Creates the local function without initializing the fields depending on the current entity.

Note
Before using the local function it must be initilized by
localFunction.init( entity );
const EntityType & entity() const
obtain the entity, this local function lives on
Definition: localfunction.hh:305
Parameters
[in]dfdiscrete function the local function shall belong to

◆ ConstLocalDiscreteFunction() [2/2]

template<class DiscreteFunction >
Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::ConstLocalDiscreteFunction ( const DiscreteFunctionType &  df,
const EntityType entity 
)
inline

constructor creating a local function and binding it to an entity

Creates the local function and initilizes the fields depending on the current entity. It is not necessary, though allowed, to call init before using the discrete function.

Note
The degrees of freedom are not initialized by this function.
Parameters
[in]dfdiscrete function the local function shall belong to
[in]entityentity for initialize the local function to

References Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::entity(), and Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::localDofVector().

Member Function Documentation

◆ assign()

template<class BasisFunctionSet , class LocalDofVector >
template<class T >
void Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::assign ( const LocalFunction< BasisFunctionSet, T > &  other)
inlineprotectedinherited

assign all DoFs of this local function

Parameters
[in]lflocal function to assign DoFs from

◆ axpy() [1/3]

template<class BasisFunctionSet , class LocalDofVector >
template<class PointType >
void Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::axpy ( const PointType &  x,
const JacobianRangeType factor 
)
inlineprotectedinherited

axpy operation for local function

Denoting the DoFs of the local function by \(u_i\) and the basis functions by \(\varphi_i\), this function performs the following operation:

\[ u_i = u_i + factor \cdot \nabla\varphi_i( x ) \]

Parameters
[in]xpoint to evaluate jacobian of basis functions in
[in]factoraxpy factor

◆ axpy() [2/3]

template<class BasisFunctionSet , class LocalDofVector >
template<class PointType >
void Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::axpy ( const PointType &  x,
const RangeType factor 
)
inlineprotectedinherited

axpy operation for local function

Denoting the DoFs of the local function by \(u_i\) and the basis functions by \(\varphi_i\), this function performs the following operation:

\[ u_i = u_i + factor \cdot \varphi_i( x ) \]

Parameters
[in]xpoint to evaluate basis functions in
[in]factoraxpy factor

◆ axpy() [3/3]

template<class BasisFunctionSet , class LocalDofVector >
template<class PointType >
void Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::axpy ( const PointType &  x,
const RangeType factor1,
const JacobianRangeType factor2 
)
inlineprotectedinherited

axpy operation for local function

Denoting the DoFs of the local function by \(u_i\) and the basis functions by \(\varphi_i\), this function performs the following operation:

\[ u_i = u_i + factor1 \cdot \varphi_i( x ) + factor2 \cdot \nabla\varphi_i( x ) \]

Parameters
[in]xpoint to evaluate basis functions in
[in]factor1axpy factor for \(\varphi( x )\)
[in]factor2axpy factor for \(\nabla\varphi( x )\)

◆ basisFunctionSet()

◆ beforeBegin() [1/2]

Iterator Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::beforeBegin ( )
inlineinherited
Returns
an iterator that is positioned before the first entry of the vector.

◆ beforeBegin() [2/2]

ConstIterator Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::beforeBegin ( ) const
inlineinherited
Returns
an iterator that is positioned before the first entry of the vector.

◆ beforeEnd() [1/2]

Iterator Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::beforeEnd ( )
inlineinherited
Returns
an iterator that is positioned before the end iterator of the vector, i.e. at the last entry.

◆ beforeEnd() [2/2]

ConstIterator Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::beforeEnd ( ) const
inlineinherited
Returns
an iterator that is positioned before the end iterator of the vector. i.e. at the last element

◆ dot()

PromotionTraits< field_type, typenameDenseVector< Other >::field_type >::PromotedType Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::dot ( const DenseVector< Other > &  x) const
inlineinherited

vector dot product \(\left (x^H \cdot y \right)\) which corresponds to Petsc's VecDot

http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Vec/VecDot.html

Parameters
xother vector
Returns

◆ entity()

◆ evaluate() [1/2]

template<class DiscreteFunction >
template<class Point >
RangeType Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::evaluate ( const Point &  p) const
inline

evaluate the local function

Parameters
[in]xevaluation point in local coordinates
Returns
value of the function in the given point

References Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::evaluate().

Referenced by Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::evaluate().

◆ evaluate() [2/2]

template<class DiscreteFunction >
template<class PointType >
void Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::evaluate ( const PointType &  x,
RangeType &  ret 
) const
inline

evaluate the local function

Parameters
[in]xevaluation point in local coordinates
[out]retvalue of the function in the given point

◆ geometry()

template<class BasisFunctionSet , class LocalDofVector >
const Geometry & Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::geometry ( ) const
inlineinherited

obtain the geometry, this local function lives on

Returns
reference to the geometry

References Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::basisFunctionSet().

◆ hessian() [1/2]

template<class DiscreteFunction >
template<class Point >
HessianRangeType Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::hessian ( const Point &  p) const
inline

evaluate Hessian of the local function

Note
Though the Hessian is evaluated on the reference element, the return value is the Hessian with respect to the actual entity.
Parameters
[in]xevaluation point in local coordinates
Returns
Hessian of the function in the evaluation point

References Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::hessian().

Referenced by Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::hessian().

◆ hessian() [2/2]

template<class DiscreteFunction >
template<class PointType >
void Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::hessian ( const PointType &  x,
HessianRangeType ret 
) const
inline

evaluate Hessian of the local function

Note
Though the Hessian is evaluated on the reference element, the return value is the Hessian with respect to the actual entity.
Parameters
[in]xevaluation point in local coordinates
[out]retHessian of the function in the evaluation point

◆ init() [1/2]

template<class BasisFunctionSet , class LocalDofVector >
void Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::init ( const BasisFunctionSetType basisFunctionSet)
inlineprotectedinherited

initialize the local function for an basisFunctionSet

Binds the local function to an basisFunctionSet and entity.

Note
A local function must be initialized to an entity before it can be used.
This function can be called multiple times to use the local function for more than one entity.
Parameters
[in]basisFunctionSetto bind the local function to

References Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::basisFunctionSet(), and Dune::Fem::BasisFunctionSet< Entity, Range >::size().

◆ init() [2/2]

template<class DiscreteFunction >
void Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::init ( const EntityType entity)
inline

interface for local functions :: init

Local functions are used to represend a discrete function on one entity. The LocalFunctionInterface defines the functionality that can be expected from such a local function. :: init

References Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::basisFunctionSet(), Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::entity(), and Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::localDofVector().

◆ jacobian() [1/2]

template<class DiscreteFunction >
template<class Point >
JacobianRangeType Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::jacobian ( const Point &  p) const
inline

evaluate Jacobian of the local function

Note
Though the Jacobian is evaluated on the reference element, the return value is the Jacobian with respect to the actual entity.
Parameters
[in]xevaluation point in local coordinates
Returns
Jacobian of the function in the evaluation point

References Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::jacobian().

Referenced by Dune::Fem::ConstLocalDiscreteFunction< DiscreteFunction >::jacobian().

◆ jacobian() [2/2]

template<class DiscreteFunction >
template<class PointType >
void Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::jacobian ( const PointType &  x,
JacobianRangeType &  ret 
) const
inline

evaluate Jacobian of the local function

Note
Though the Jacobian is evaluated on the reference element, the return value is the Jacobian with respect to the actual entity.
Parameters
[in]xevaluation point in local coordinates
[out]retJacobian of the function in the evaluation point

◆ numDofs()

template<class BasisFunctionSet , class LocalDofVector >
int Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::numDofs ( ) const
inlineinherited

obtain the number of local DoFs

Obtain the number of local DoFs of this local function. The value is identical to the number of basis functons on the entity.

Returns
number of local DoFs

References Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::localDofVector().

◆ operator*()

PromotionTraits< field_type, typenameDenseVector< Other >::field_type >::PromotedType Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::operator* ( const DenseVector< Other > &  x) const
inlineinherited

indefinite vector dot product \(\left (x^T \cdot y \right)\) which corresponds to Petsc's VecTDot

http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Vec/VecTDot.html

Parameters
xother vector
Returns

◆ operator*=()

std::enable_if< std::is_convertible< FieldType, field_type >::value, derived_type >::type & Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::operator*= ( const FieldType &  kk)
inlineinherited

vector space multiplication with scalar

we use enable_if to avoid an ambiguity, if the function parameter can be converted to field_type implicitly. (see FS#1457)

The function is only enabled, if the parameter is directly convertible to field_type.

◆ operator+=()

std::enable_if< std::is_convertible< ValueType, value_type >::value, derived_type >::type & Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::operator+= ( const ValueType &  kk)
inlineinherited

vector space add scalar to all comps

we use enable_if to avoid an ambiguity, if the function parameter can be converted to value_type implicitly. (see FS#1457)

The function is only enabled, if the parameter is directly convertible to value_type.

◆ operator-=()

std::enable_if< std::is_convertible< ValueType, value_type >::value, derived_type >::type & Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::operator-= ( const ValueType &  kk)
inlineinherited

vector space subtract scalar from all comps

we use enable_if to avoid an ambiguity, if the function parameter can be converted to value_type implicitly. (see FS#1457)

The function is only enabled, if the parameter is directly convertible to value_type.

◆ operator/=()

std::enable_if< std::is_convertible< FieldType, field_type >::value, derived_type >::type & Dune::DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > >::operator/= ( const FieldType &  kk)
inlineinherited

vector space division by scalar

we use enable_if to avoid an ambiguity, if the function parameter can be converted to field_type implicitly. (see FS#1457)

The function is only enabled, if the parameter is directly convertible to field_type.

◆ order()

template<class BasisFunctionSet , class LocalDofVector >
int Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::order ( ) const
inlineinherited

obtain the order of this local function

The order of a local function refers to the polynomial order required to integrate it exactly.

Note
It is not completely clear what this value should be, e.g., for bilinear basis functions.
Returns
order of the local function

References Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::basisFunctionSet(), and Dune::Fem::BasisFunctionSet< Entity, Range >::order().

Referenced by Dune::Fem::LocalOrthonormalL2Projection< GridPart, BasisFunctionSet, Quadrature >::apply().

◆ size()

template<class BasisFunctionSet , class LocalDofVector >
SizeType Dune::Fem::LocalFunction< BasisFunctionSet, LocalDofVector >::size ( ) const
inlineinherited

Friends And Related Function Documentation

◆ operator<<()

std::ostream & operator<< ( std::ostream &  s,
const DenseVector< LocalFunction< BasisFunctionSet, LocalDofVector > > &  v 
)
related

Write a DenseVector to an output stream.

Parameters
[in]sstd :: ostream to write to
[in]vDenseVector to write
Returns
the output stream (s)

The documentation for this class was generated from the following file:
  • dune/fem/function/localfunction/const.hh
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 20, 23:30, 2024)