dune-fem  2.8-git
operator/projection/local/riesz/orthonormal.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_ORTHONORMAL_HH
2 #define DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_ORTHONORMAL_HH
3 
4 #include <cassert>
5 #include <cstddef>
6 
7 #include <utility>
8 
9 #include <dune/geometry/referenceelements.hh>
10 
11 #include "localrieszprojection.hh"
12 
13 namespace Dune
14 {
15 
16  namespace Fem
17  {
18 
19  // OrthonormalLocalRieszProjection
20  // -------------------------------
21 
22  template< class BasisFunctionSet >
24  : public LocalRieszProjection< BasisFunctionSet, OrthonormalLocalRieszProjection< BasisFunctionSet > >
25  {
28 
29  public:
32 
33  private:
34  typedef typename BasisFunctionSet::FunctionSpaceType::RangeFieldType RangeFieldType;
35 
36  public:
42  : basisFunctionSet_( std::forward< BasisFunctionSetType >( basisFunctionSet ) ),
43  factor_( ratio( basisFunctionSet.entity().geometry() ) )
44  {}
45 
47  : basisFunctionSet_( std::forward< BasisFunctionSetType >( basisFunctionSet ) ),
48  factor_( ratio( basisFunctionSet.entity().geometry() ) )
49  {}
50 
58 
60  : basisFunctionSet_( std::move( other.basisFunctionSet_ ) ),
61  factor_( other.factor_ )
62  {}
63 
64  ThisType &operator= ( const ThisType & ) = default;
65 
67  {
68  basisFunctionSet_ = std::move( other.basisFunctionSet_ );
69  factor_( other.factor_ );
70  return *this;
71  }
72 
81  {
82  return basisFunctionSet_;
83  }
84 
86  template< class F, class LocalDofVector >
87  void apply ( const F &f, LocalDofVector &dofs ) const
88  {
89  assert( f.size() == dofs.size() );
90  const std::size_t size = dofs.size();
91  for( std::size_t i = 0u; i < size; ++i )
92  dofs[ i ] = factor_*f[ i ];
93  }
94 
97  private:
98  template< class Geometry >
99  static RangeFieldType ratio ( const Geometry &geometry )
100  {
101  assert( geometry.affine() );
102  const auto &referenceElement
103  = Dune::ReferenceElements< typename Geometry::ctype, Geometry::mydimension >::general( geometry.type() );
104  return referenceElement.volume()/geometry.volume();
105  }
106 
107  BasisFunctionSetType basisFunctionSet_;
108  RangeFieldType factor_;
109  };
110 
111  } // namespace Fem
112 
113 } // namepsace Dune
114 
115 #endif // #ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_ORTHONORMAL_HH
Definition: bindguard.hh:11
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
interface documentation of a local Riesz projection
Definition: localrieszprojection.hh:23
Definition: operator/projection/local/riesz/orthonormal.hh:25
OrthonormalLocalRieszProjection(const ThisType &)=default
OrthonormalLocalRieszProjection(ThisType &&other)
Definition: operator/projection/local/riesz/orthonormal.hh:59
OrthonormalLocalRieszProjection(BasisFunctionSetType &&basisFunctionSet)
Definition: operator/projection/local/riesz/orthonormal.hh:46
BaseType::BasisFunctionSetType BasisFunctionSetType
basis function set
Definition: operator/projection/local/riesz/orthonormal.hh:31
OrthonormalLocalRieszProjection(const BasisFunctionSetType &basisFunctionSet)
Definition: operator/projection/local/riesz/orthonormal.hh:41
ThisType & operator=(const ThisType &)=default
void apply(const F &f, LocalDofVector &dofs) const
compute Riesz representation
Definition: operator/projection/local/riesz/orthonormal.hh:87
BasisFunctionSetType basisFunctionSet() const
return basis function set
Definition: operator/projection/local/riesz/orthonormal.hh:80
Interface class for basis function sets.
Definition: basisfunctionset/basisfunctionset.hh:31
FunctionSpaceTraits::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:63