1 #ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_ORTHONORMAL_HH 2 #define DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_ORTHONORMAL_HH 9 #include <dune/geometry/referenceelements.hh> 22 template<
class BasisFunctionSet >
24 :
public LocalRieszProjection< BasisFunctionSet, OrthonormalLocalRieszProjection< BasisFunctionSet > >
42 : basisFunctionSet_(
std::forward< BasisFunctionSetType >( basisFunctionSet ) ),
43 factor_( ratio( basisFunctionSet.entity().geometry() ) )
60 : basisFunctionSet_(
std::
move( other.basisFunctionSet_ ) ),
61 factor_( other.factor_ )
64 ThisType &
operator= (
const ThisType & ) =
default;
68 basisFunctionSet_ =
std::move( other.basisFunctionSet_ );
69 factor_( other.factor_ );
82 return basisFunctionSet_;
86 template<
class F,
class LocalDofVector >
87 void apply (
const F &f, LocalDofVector &dofs )
const 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 ];
98 template<
class Geometry >
99 static RangeFieldType ratio (
const Geometry &geometry )
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();
107 BasisFunctionSetType basisFunctionSet_;
108 RangeFieldType factor_;
115 #endif // #ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_RIESZ_ORTHONORMAL_HH OrthonormalLocalRieszProjection(BasisFunctionSetType &&basisFunctionSet)
Definition: operator/projection/local/riesz/orthonormal.hh:46
OrthonormalLocalRieszProjection(ThisType &&other)
Definition: operator/projection/local/riesz/orthonormal.hh:59
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:62
BasisFunctionSetType basisFunctionSet() const
return basis function set
Definition: operator/projection/local/riesz/orthonormal.hh:80
interface documentation of a local Riesz projection
Definition: localrieszprojection.hh:22
BaseType::BasisFunctionSetType BasisFunctionSetType
basis function set
Definition: operator/projection/local/riesz/orthonormal.hh:31
Definition: coordinate.hh:4
Definition: operator/projection/local/riesz/orthonormal.hh:23
void move(ArrayInterface< T > &array, const unsigned int oldOffset, const unsigned int newOffset, const unsigned int length)
Definition: array_inline.hh:38
OrthonormalLocalRieszProjection(const BasisFunctionSetType &basisFunctionSet)
Definition: operator/projection/local/riesz/orthonormal.hh:41
BasisFunctionSet BasisFunctionSetType
basis function set
Definition: localrieszprojection.hh:26
void apply(const F &f, LocalDofVector &dofs) const
compute Riesz representation
Definition: operator/projection/local/riesz/orthonormal.hh:87
ThisType & operator=(const ThisType &)=default