1#ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_DENSE_HH
2#define DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_DENSE_HH
10#include <dune/common/dynmatrix.hh>
23 template<
class BasisFunctionSet,
class Quadrature >
25 :
public LocalRieszProjection< BasisFunctionSet, DenseLocalRieszProjection< BasisFunctionSet, Quadrature > >
63 : basisFunctionSet_(
std::move( other.basisFunctionSet_ ) ),
64 inverse_(
std::move( other.inverse_ ) )
71 basisFunctionSet_ = std::move( other.basisFunctionSet_ );
72 inverse_ = std::move( other.inverse_ );
86 template<
class F,
class LocalDofVector >
87 void apply (
const F&f, LocalDofVector &
dofs )
const
89 inverse_.mv( f,
dofs );
95 template<
class DenseMatrix >
99 std::vector< typename BasisFunctionSetType::RangeType > phi( size );
101 matrix =
static_cast< typename DenseMatrix::value_type
>( 0 );
107 const std::size_t nop = quadrature.nop();
108 for( std::size_t qp = 0; qp < nop; ++qp )
112 const auto weight = quadrature.weight( qp )*geometry.integrationElement( quadrature.point( qp ) );
113 for( std::size_t i = 0; i < size; ++i )
115 matrix[ i ][ i ] += weight*( phi[ i ] * phi[ i ] );
116 for( std::size_t j = i+1; j < size; ++j )
118 const auto value = weight*( phi[ i ]* phi[ j ] );
119 matrix[ i ][ j ] += value;
120 matrix[ j ][ i ] += value;
127 Dune::DynamicMatrix< typename BasisFunctionSetType::FunctionSpaceType::RangeFieldType > inverse_;
Definition: bindguard.hh:11
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
DenseMatrix based on std::vector< std::vector< T > >
Definition: blockmatrix.hh:24
BasisFunctionSetType basisFunctionSet() const
return basis function set
Definition: dense.hh:83
ThisType & operator=(const ThisType &other)=default
DenseLocalRieszProjection(const BasisFunctionSetType &basisFunctionSet)
Definition: dense.hh:38
void apply(const F &f, LocalDofVector &dofs) const
compute Riesz representation
Definition: dense.hh:87
DenseLocalRieszProjection(const ThisType &other)=default
DenseLocalRieszProjection(ThisType &&other)
Definition: dense.hh:62
DenseLocalRieszProjection(BasisFunctionSetType &&basisFunctionSet)
Definition: dense.hh:46
BaseType::BasisFunctionSetType BasisFunctionSetType
Definition: dense.hh:32
interface documentation of a local Riesz projection
Definition: localrieszprojection.hh:23
actual interface class for quadratures
Definition: quadrature.hh:405
Interface class for basis function sets.
Definition: basisfunctionset/basisfunctionset.hh:31
const EntityType & entity() const
return entity
std::size_t size() const
return size of basis function set
int order() const
return order of basis function set
void evaluateAll(const Quadrature &quad, const DofVector &dofs, RangeArray &ranges) const
evaluate all basis functions and store the result in the ranges array