1 #ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_L2PROJECTION_HH 2 #define DUNE_FEM_OPERATOR_PROJECTION_LOCAL_L2PROJECTION_HH 8 #include <dune/common/dynvector.hh> 9 #include <dune/common/bartonnackmanifcheck.hh> 26 template<
class BasisFunctionSet,
class Implementation >
38 #endif // #ifndef DOXYGEN 67 return impl().basisFunctionSet();
78 template<
class LocalFunction,
class LocalDofVector >
81 apply( localFunction, localDofVector );
92 template<
class LocalFunction,
class LocalDofVector >
95 CHECK_INTERFACE_IMPLEMENTATION(
impl().
apply( localFunction, localDofVector ) );
96 impl().apply( localFunction, localDofVector );
102 const Implementation &
impl ()
const 104 return static_cast< const Implementation &
>( *this );
113 template<
class LocalRieszProjection,
class Quadrature >
115 :
public LocalL2Projection< typename LocalRieszProjection::BasisFunctionSetType, DefaultLocalL2Projection< LocalRieszProjection, Quadrature > >
128 typedef typename RangeType::value_type RangeFieldType;
135 template<
class... Args >
137 : rieszProjection_(
std::forward< Args >( args )... )
154 #endif // #ifndef DOXYGEN 157 : rieszProjection_(
std::
move( other.rieszProjection_ ) )
164 rieszProjection_ =
std::move( other.rieszProjection_ );
177 return rieszProjection_.basisFunctionSet();
181 template<
class LocalFunction,
class LocalDofVector >
185 const auto geometry = basisFunctionSet.entity().geometry();
187 f_.resize( basisFunctionSet.size() );
188 f_ =
static_cast< RangeFieldType
>( 0 );
190 Quadrature quadrature( geometry.type(), localFunction.
order() + basisFunctionSet.order() );
191 const std::size_t nop = quadrature.
nop();
192 for( std::size_t qp = 0; qp < nop; ++qp )
195 localFunction.
evaluate( quadrature[ qp ], value );
196 value *= quadrature.weight( qp )*geometry.integrationElement( quadrature.point( qp ) );
197 basisFunctionSet.axpy( quadrature[ qp ], value, f_ );
200 rieszProjection_.apply( f_, dofs );
204 LocalRieszProjectionType rieszProjection_;
205 mutable Dune::DynamicVector< RangeFieldType > f_;
212 #endif // #ifndef DUNE_FEM_OPERATOR_PROJECTION_LOCAL_L2PROJECTION_HH LocalRieszProjection LocalRieszProjectionType
type of local Riesz type
Definition: local/l2projection.hh:122
DefaultLocalL2Projection(Args &&...args)
Definition: local/l2projection.hh:136
FunctionSpaceType::RangeType RangeType
range type
Definition: basisfunctionset/basisfunctionset.hh:45
Definition: local/l2projection.hh:114
void operator()(const LocalFunction &localFunction, LocalDofVector &localDofVector) const
please doc me
Definition: local/l2projection.hh:79
int nop() const
obtain the number of integration points
Definition: quadrature.hh:226
LocalL2Projection & operator=(const LocalL2Projection &)=default
assignment operator
BaseType::BasisFunctionSetType BasisFunctionSetType
basis function set type
Definition: local/l2projection.hh:124
interface documentation of a local Riesz projection
Definition: localrieszprojection.hh:22
interface for local functions
Definition: localfunction.hh:41
const Implementation & impl() const
Definition: local/l2projection.hh:102
LocalL2Projection(const LocalL2Projection &)=default
copy constructor
Definition: coordinate.hh:4
BasisFunctionSet BasisFunctionSetType
basis function set type
Definition: local/l2projection.hh:31
BasisFunctionSet basisFunctionSet() const
return basis function set
Definition: local/l2projection.hh:64
void move(ArrayInterface< T > &array, const unsigned int oldOffset, const unsigned int newOffset, const unsigned int length)
Definition: array_inline.hh:38
please doc me
Definition: local/l2projection.hh:27
LocalL2Projection(LocalL2Projection &&)
move constructor
Definition: local/l2projection.hh:48
void evaluate(const PointType &x, RangeType &ret) const
evaluate the local function
Definition: localfunction.hh:300
BasisFunctionSetType basisFunctionSet() const
return basis function set
Definition: local/l2projection.hh:175
DefaultLocalL2Projection(ThisType &&other)
Definition: local/l2projection.hh:156
int order() const
obtain the order of this local function
Definition: localfunction.hh:273
void apply(const LocalFunction &localFunction, LocalDofVector &dofs) const
please doc me
Definition: local/l2projection.hh:182
Definition: basisfunctionset/basisfunctionset.hh:31
void apply(const LocalFunction &localFunction, LocalDofVector &localDofVector) const
please doc me
Definition: local/l2projection.hh:93
actual interface class for quadratures
Definition: quadrature.hh:320