1#ifndef DUNE_FEM_L2NORM_HH
2#define DUNE_FEM_L2NORM_HH
18 template<
class Gr
idPart >
30 template<
class Function >
33 template<
class UFunction,
class VFunction >
50 const unsigned int order = 0,
51 const bool communicate =
true );
54 template<
class DiscreteFunctionType,
class PartitionSet >
55 typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
56 norm (
const DiscreteFunctionType &u,
const PartitionSet& partitionSet )
const;
59 template<
class DiscreteFunctionType >
60 typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
61 norm (
const DiscreteFunctionType &u )
const
63 return norm( u, Partitions::interior );
67 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class PartitionSet >
68 typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
70 const VDiscreteFunctionType &v,
71 const PartitionSet& partitionSet )
const;
74 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType >
75 typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
76 distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v )
const
78 return distance( u, v, Partitions::interior );
81 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
82 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &sum )
const;
84 template<
class LocalFunctionType,
class ReturnType >
85 void normLocal (
const EntityType &entity,
unsigned int order,
const LocalFunctionType &uLocal, ReturnType &sum )
const;
92 template<
class Gr
idPart >
96 communicate_(
BaseType::checkCommunicateFlag( communicate ) )
101 template<
class Gr
idPart >
102 template<
class DiscreteFunctionType,
class PartitionSet >
103 typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
106 typedef typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type RealType;
107 typedef FieldVector< RealType, 1 > ReturnType ;
110 ReturnType sum = BaseType :: forEach( u, ReturnType(0), partitionSet, order_ );
115 sum[ 0 ] = comm().sum( sum[ 0 ] );
123 template<
class Gr
idPart >
124 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class PartitionSet >
125 typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
128 const VDiscreteFunctionType &v,
129 const PartitionSet& partitionSet )
const
131 typedef typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type RealType;
132 typedef FieldVector< RealType, 1 > ReturnType ;
135 ReturnType sum = BaseType :: forEach( u, v, ReturnType(0), partitionSet, order_ );
140 sum[ 0 ] = comm().sum( sum[ 0 ] );
147 template<
class Gr
idPart >
148 template<
class LocalFunctionType,
class ReturnType >
159 template<
class Gr
idPart >
160 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
169 LocalDistanceType dist( uLocal, vLocal );
176 template<
class Gr
idPart >
177 template<
class Function >
183 typedef typename Dune::FieldTraits< RangeFieldType >::real_type
RealType;
187 : function_( function )
190 template<
class Po
int >
194 function_.evaluate( x, phi );
195 ret = phi.two_norm2();
199 const FunctionType &function_;
203 template<
class Gr
idPart >
204 template<
class UFunction,
class VFunction >
218 template<
class Po
int >
222 u_.evaluate( x, ret );
223 v_.evaluate( x, phi );
227 template<
class Po
int >
231 u_.jacobian( x, ret );
232 v_.jacobian( x, phi );
237 const UFunctionType &u_;
238 const VFunctionType &v_;
246 template<
class WeightFunction >
248 :
public L2Norm< typename WeightFunction::DiscreteFunctionSpaceType::GridPartType >
260 template<
class Function >
279 template<
class LocalFunctionType,
class ReturnType >
280 void normLocal (
const EntityType &entity,
unsigned int order,
const LocalFunctionType &uLocal, ReturnType &sum )
const;
282 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
283 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &sum )
const;
295 template<
class WeightFunction >
298 :
BaseType( weightFunction.space().gridPart(), order, communicate ),
299 wfLocal_( weightFunction )
301 static_assert( (WeightFunctionSpaceType::dimRange == 1),
302 "Wight function must be scalar." );
306 template<
class WeightFunction >
307 template<
class LocalFunctionType,
class ReturnType >
314 wfLocal_.bind( entity );
321 template<
class WeightFunction >
322 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
326 typedef typename BaseType::template FunctionDistance< ULocalFunctionType, VLocalFunctionType > LocalDistanceType;
331 wfLocal_.bind( entity );
332 LocalDistanceType dist( uLocal, vLocal );
340 template<
class WeightFunction >
341 template<
class Function >
347 typedef typename Dune::FieldTraits< RangeFieldType >::real_type
RealType;
352 : weightFunction_( weightFunction ),
353 function_( function )
356 template<
class Po
int >
360 weightFunction_.evaluate( x, weight );
363 function_.evaluate( x, phi );
364 ret = weight[ 0 ] * (phi * phi);
369 const FunctionType &function_;
double sqrt(const Dune::Fem::Double &v)
Definition: double.hh:977
Definition: bindguard.hh:11
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:604
Abstract class representing a function.
Definition: common/function.hh:50
FunctionSpaceType::RangeType RangeType
range type
Definition: common/function.hh:68
FunctionSpaceType::RangeFieldType RangeFieldType
field type of range
Definition: common/function.hh:64
Definition: domainintegral.hh:33
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: domainintegral.hh:44
const GridPartType::CollectiveCommunicationType & comm() const
Definition: domainintegral.hh:244
const GridPartType & gridPart() const
Definition: domainintegral.hh:242
Dune::FieldTraits< typenameDiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u) const
|| u ||_L2 on interior partition entities
Definition: l2norm.hh:61
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: l2norm.hh:149
GridPart GridPartType
Definition: l2norm.hh:25
const unsigned int order_
Definition: l2norm.hh:41
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: l2norm.hh:162
Integrator< QuadratureType > IntegratorType
Definition: l2norm.hh:39
BaseType::EntityType EntityType
Definition: l2norm.hh:37
const bool communicate_
Definition: l2norm.hh:42
Dune::FieldTraits< typenameUDiscreteFunctionType::RangeFieldType >::real_type distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const PartitionSet &partitionSet) const
|| u - v ||_L2 on given set of entities (partition set)
Definition: l2norm.hh:127
Dune::FieldTraits< typenameUDiscreteFunctionType::RangeFieldType >::real_type distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v) const
|| u - v ||_L2 on interior partition entities
Definition: l2norm.hh:76
Dune::FieldTraits< typenameDiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u, const PartitionSet &partitionSet) const
|| u ||_L2 on given set of entities (partition set)
Definition: l2norm.hh:104
CachingQuadrature< GridPartType, 0 > QuadratureType
Definition: l2norm.hh:38
L2Norm(const GridPartType &gridPart, const unsigned int order=0, const bool communicate=true)
constructor
Definition: l2norm.hh:93
Definition: l2norm.hh:179
FunctionSquare(const FunctionType &function)
Definition: l2norm.hh:186
Function FunctionType
Definition: l2norm.hh:180
FieldVector< RealType, 1 > RangeType
Definition: l2norm.hh:184
void evaluate(const Point &x, RangeType &ret) const
Definition: l2norm.hh:191
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: l2norm.hh:183
FunctionType::RangeFieldType RangeFieldType
Definition: l2norm.hh:182
Definition: l2norm.hh:206
UFunctionType::JacobianRangeType JacobianRangeType
Definition: l2norm.hh:212
void evaluate(const Point &x, RangeType &ret) const
Definition: l2norm.hh:219
VFunction VFunctionType
Definition: l2norm.hh:208
void jacobian(const Point &x, JacobianRangeType &ret) const
Definition: l2norm.hh:228
UFunctionType::RangeFieldType RangeFieldType
Definition: l2norm.hh:210
FunctionDistance(const UFunctionType &u, const VFunctionType &v)
Definition: l2norm.hh:214
UFunctionType::RangeType RangeType
Definition: l2norm.hh:211
UFunction UFunctionType
Definition: l2norm.hh:207
Definition: l2norm.hh:249
BaseType::IntegratorType IntegratorType
Definition: l2norm.hh:267
WeightFunction WeightFunctionType
Definition: l2norm.hh:254
WeightFunctionType::DiscreteFunctionSpaceType WeightFunctionSpaceType
Definition: l2norm.hh:256
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: l2norm.hh:309
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: l2norm.hh:324
ConstLocalFunction< WeightFunctionType > LocalWeightFunctionType
Definition: l2norm.hh:263
BaseType::EntityType EntityType
Definition: l2norm.hh:266
WeightFunctionSpaceType::GridPartType GridPartType
Definition: l2norm.hh:257
WeightFunctionType::RangeType WeightType
Definition: l2norm.hh:264
Definition: l2norm.hh:343
WeightedFunctionSquare(const LocalWeightFunctionType &weightFunction, const FunctionType &function)
Definition: l2norm.hh:350
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: l2norm.hh:347
Function FunctionType
Definition: l2norm.hh:344
FieldVector< RealType, 1 > RangeType
Definition: l2norm.hh:348
FunctionType::RangeFieldType RangeFieldType
Definition: l2norm.hh:346
void evaluate(const Point &x, RangeType &ret) const
Definition: l2norm.hh:357
integrator for arbitrary functions providing evaluate
Definition: integrator.hh:28
void integrateAdd(const EntityType &entity, const Function &function, typename Function ::RangeType &ret) const
add the integral over an entity to a variable
Definition: integrator.hh:67