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 > >
39 : basisFunctionSet_( basisFunctionSet ),
40 inverse_( basisFunctionSet.size(), basisFunctionSet.size() )
42 assemble( basisFunctionSet, inverse_ );
63 : basisFunctionSet_(
std::
move( other.basisFunctionSet_ ) ),
64 inverse_(
std::
move( other.inverse_ ) )
67 ThisType &
operator= (
const ThisType &other ) =
default;
71 basisFunctionSet_ =
std::move( other.basisFunctionSet_ );
86 template<
class F,
class LocalDofVector >
87 void apply (
const F&f, LocalDofVector &dofs )
const 89 inverse_.mv( f, dofs );
95 template<
class DenseMatrix >
98 const std::size_t size = basisFunctionSet.size();
99 std::vector< typename BasisFunctionSetType::RangeType > phi( size );
101 matrix =
static_cast< typename DenseMatrix::value_type
>( 0 );
102 assert( matrix.N() == basisFunctionSet.size() );
103 assert( matrix.M() == basisFunctionSet.size() );
105 const auto &geometry = basisFunctionSet.entity().geometry();
106 const Quadrature quadrature( geometry.type(), 2*basisFunctionSet.order() );
107 const std::size_t nop = quadrature.
nop();
108 for( std::size_t qp = 0; qp < nop; ++qp )
110 basisFunctionSet.evaluateAll( quadrature[ qp ], phi );
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;
126 BasisFunctionSetType basisFunctionSet_;
127 Dune::DynamicMatrix< typename BasisFunctionSetType::FunctionSpaceType::RangeFieldType > inverse_;
134 #endif // #ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_DENSE_HH
DenseLocalRieszProjection(BasisFunctionSetType &&basisFunctionSet)
Definition: dense.hh:46
int nop() const
obtain the number of integration points
Definition: quadrature.hh:226
DenseMatrix based on std::vector< std::vector< T > >
Definition: blockmatrix.hh:23
interface documentation of a local Riesz projection
Definition: localrieszprojection.hh:22
Definition: coordinate.hh:4
void apply(const F &f, LocalDofVector &dofs) const
compute Riesz representation
Definition: dense.hh:87
BaseType::BasisFunctionSetType BasisFunctionSetType
Definition: dense.hh:32
BasisFunctionSetType basisFunctionSet() const
return basis function set
Definition: dense.hh:83
void move(ArrayInterface< T > &array, const unsigned int oldOffset, const unsigned int newOffset, const unsigned int length)
Definition: array_inline.hh:38
DenseLocalRieszProjection(const BasisFunctionSetType &basisFunctionSet)
Definition: dense.hh:38
ThisType & operator=(const ThisType &other)=default
DenseLocalRieszProjection(ThisType &&other)
Definition: dense.hh:62
Definition: basisfunctionset/basisfunctionset.hh:31
actual interface class for quadratures
Definition: quadrature.hh:320