dune-fem 2.8.0
Loading...
Searching...
No Matches
dense.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_DENSE_HH
2#define DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_DENSE_HH
3
4#include <cassert>
5#include <cstddef>
6
7#include <utility>
8#include <vector>
9
10#include <dune/common/dynmatrix.hh>
11
13
14namespace Dune
15{
16
17 namespace Fem
18 {
19
20 // DenseLocalRieszProjection
21 // -------------------------
22
23 template< class BasisFunctionSet, class Quadrature >
25 : public LocalRieszProjection< BasisFunctionSet, DenseLocalRieszProjection< BasisFunctionSet, Quadrature > >
26 {
29
30 public:
33
39 : basisFunctionSet_( basisFunctionSet ),
40 inverse_( basisFunctionSet.size(), basisFunctionSet.size() )
41 {
42 assemble( basisFunctionSet, inverse_ );
43 inverse_.invert();
44 }
45
47 : basisFunctionSet_( std::forward< BasisFunctionSetType >( basisFunctionSet ) ),
48 inverse_( basisFunctionSet.size(), basisFunctionSet.size() )
49 {
50 assemble( basisFunctionSet, inverse_ );
51 inverse_.invert();
52 }
53
60 DenseLocalRieszProjection ( const ThisType &other ) = default;
61
63 : basisFunctionSet_( std::move( other.basisFunctionSet_ ) ),
64 inverse_( std::move( other.inverse_ ) )
65 {}
66
67 ThisType &operator= ( const ThisType &other ) = default;
68
70 {
71 basisFunctionSet_ = std::move( other.basisFunctionSet_ );
72 inverse_ = std::move( other.inverse_ );
73 return *this;
74 }
75
83 BasisFunctionSetType basisFunctionSet () const { return basisFunctionSet_; }
84
86 template< class F, class LocalDofVector >
87 void apply ( const F&f, LocalDofVector &dofs ) const
88 {
89 inverse_.mv( f, dofs );
90 }
91
94 private:
95 template< class DenseMatrix >
96 void assemble ( const BasisFunctionSetType &basisFunctionSet, DenseMatrix &matrix )
97 {
98 const std::size_t size = basisFunctionSet.size();
99 std::vector< typename BasisFunctionSetType::RangeType > phi( size );
100
101 matrix = static_cast< typename DenseMatrix::value_type >( 0 );
102 assert( matrix.N() == basisFunctionSet.size() );
103 assert( matrix.M() == basisFunctionSet.size() );
104
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 )
109 {
110 basisFunctionSet.evaluateAll( quadrature[ qp ], phi );
111
112 const auto weight = quadrature.weight( qp )*geometry.integrationElement( quadrature.point( qp ) );
113 for( std::size_t i = 0; i < size; ++i )
114 {
115 matrix[ i ][ i ] += weight*( phi[ i ] * phi[ i ] );
116 for( std::size_t j = i+1; j < size; ++j )
117 {
118 const auto value = weight*( phi[ i ]* phi[ j ] );
119 matrix[ i ][ j ] += value;
120 matrix[ j ][ i ] += value;
121 }
122 }
123 }
124 }
125
126 BasisFunctionSetType basisFunctionSet_;
127 Dune::DynamicMatrix< typename BasisFunctionSetType::FunctionSpaceType::RangeFieldType > inverse_;
128 };
129
130 } // namespace Fem
131
132} // namepsace Dune
133
134#endif // #ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_DENSE_HH
135
STL namespace.
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