1 #ifndef DUNE_FEM_LPNORM_HH 2 #define DUNE_FEM_LPNORM_HH 4 #include <dune/common/typetraits.hh> 6 #include <dune/grid/common/rangegenerators.hh> 13 #include <dune/common/bartonnackmanifcheck.hh> 25 template<
class Gr
idPart,
class NormImplementation >
40 typedef typename GridPartType::template Codim< 0 >::EntityType
EntityType;
42 template <
bool uDiscrete,
bool vDiscrete>
45 template <
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class ReturnType>
46 static ReturnType
forEach (
const ThisType &norm,
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
const ReturnType &initialValue,
unsigned int order )
48 static_assert( uDiscrete && vDiscrete,
"Distance can only be calculated between GridFunctions" );
50 order = (order == 0 ? 2*
std::max( u.space().order(), v.space().order() ) : order);
53 const auto end = norm.gridPart_.template end<0>();
54 for (
auto it = norm.gridPart_.template begin<0>(); it!=end; ++it)
56 const auto entity = *it;
57 const typename UDiscreteFunctionType::LocalFunctionType uLocal = u.localFunction( entity );
58 const typename VDiscreteFunctionType::LocalFunctionType vLocal = v.localFunction( entity );
66 template <
bool vDiscrete>
70 class VDiscreteFunctionType,
73 ReturnType
forEach (
const ThisType& norm,
74 const F& f,
const VDiscreteFunctionType& v,
75 const ReturnType& initialValue,
76 const unsigned int order )
79 GridFunction u(
"LPNorm::adapter" , f , v.space().gridPart(), v.space().order() );
82 forEach( norm, u, v, initialValue, order );
87 template <
bool uDiscrete>
90 template <
class UDiscreteFunctionType,
94 ReturnType
forEach (
const ThisType& norm,
95 const UDiscreteFunctionType& u,
97 const ReturnType& initialValue,
98 const unsigned int order )
101 forEach( norm, f, u, initialValue, order );
105 template<
class DiscreteFunctionType,
class ReturnType >
106 ReturnType
forEach (
const DiscreteFunctionType &u,
const ReturnType &initialValue,
unsigned int order = 0 )
const 108 static_assert( (IsBaseOf<Fem::HasLocalFunction, DiscreteFunctionType>::value),
109 "Norm only implemented for quantities implementing a local function!" );
111 order = (order == 0 ? 2*u.space().order() : order);
114 const auto end = gridPart_.template end<0>();
115 for (
auto it = gridPart_.template begin<0>(); it!=end; ++it)
118 const EntityType &entity = *it;
119 typename DiscreteFunctionType::LocalFunctionType uLocal = u.localFunction( entity );
125 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class ReturnType >
126 ReturnType
forEach (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
const ReturnType &initialValue,
unsigned int order = 0 )
const 128 enum { uDiscrete = std::is_convertible<UDiscreteFunctionType, HasLocalFunction>::value };
129 enum { vDiscrete = std::is_convertible<VDiscreteFunctionType, HasLocalFunction>::value };
134 forEach( *
this, u, v, initialValue, order );
139 : gridPart_( gridPart )
143 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
144 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &
sum )
const 146 CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
asImp().
distanceLocal( entity, order, uLocal, vLocal, sum ) );
149 template<
class LocalFunctionType,
class ReturnType >
150 void normLocal (
const EntityType &entity,
unsigned int order,
const LocalFunctionType &uLocal, ReturnType &
sum )
const 152 CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
asImp().
normLocal( entity, order, uLocal, sum ) );
155 const GridPartType &
gridPart ()
const {
return gridPart_; }
157 typename GridPartType::CollectiveCommunicationType
comm ()
const 163 const GridPartType &gridPart_;
179 virtual int operator() (
const double p)=0;
188 int operator() (
const double p)
191 const double q = p / (p-1);
199 template<
class Gr
idPart,
class OrderCalculator = DefaultOrderCalculator >
208 using BaseType::gridPart;
209 using BaseType::comm;
212 template<
class Function >
215 template<
class UFunction,
class VFunction >
216 struct FunctionDistance;
218 typedef typename BaseType::EntityType
EntityType;
225 template<
class DiscreteFunctionType >
226 typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
227 norm (
const DiscreteFunctionType &u )
const;
229 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType >
230 typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
231 distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v )
const;
233 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
234 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &
sum )
const;
236 template<
class LocalFunctionType,
class ReturnType >
237 void normLocal (
const EntityType &entity,
unsigned int order,
const LocalFunctionType &uLocal, ReturnType &
sum )
const;
239 int order (
const int spaceOrder )
const ;
252 template<
class WeightFunction,
class OrderCalculator = DefaultOrderCalculator >
254 :
public LPNorm< typename WeightFunction::DiscreteFunctionSpaceType::GridPartType,
267 template<
class Function >
276 using BaseType::gridPart;
277 using BaseType::comm;
280 using BaseType::norm;
281 using BaseType::distance;
283 explicit WeightedLPNorm (
const WeightFunctionType &weightFunction,
const double p );
285 template<
class LocalFunctionType,
class ReturnType >
286 void normLocal (
const EntityType &entity,
unsigned int order,
const LocalFunctionType &uLocal, ReturnType &
sum )
const;
288 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
289 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &
sum )
const;
292 const WeightFunctionType &weightFunction_;
300 template<
class Gr
idPart,
class OrderCalculator >
304 orderCalc_( new OrderCalculator() )
308 template<
class Gr
idPart,
class OrderCalculator>
315 template<
class Gr
idPart,
class OrderCalculator>
316 template<
class DiscreteFunctionType >
317 inline typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
320 typedef typename DiscreteFunctionType::RangeFieldType RangeFieldType;
321 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
322 typedef FieldVector< RealType, 1 > ReturnType ;
328 return std::pow (
comm().
sum( sum[ 0 ] ), (1.0 /
p_) );
331 template<
class Gr
idPart,
class OrderCalculator >
332 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType >
333 inline typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
335 ::distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v )
const 337 typedef typename UDiscreteFunctionType::RangeFieldType RangeFieldType;
338 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
339 typedef FieldVector< RealType, 1 > ReturnType ;
345 return std::pow(
comm().
sum( sum[ 0 ] ), (1.0/
p_) );
348 template<
class Gr
idPart,
class OrderCalculator >
349 template<
class LocalFunctionType,
class ReturnType >
360 template<
class Gr
idPart,
class OrderCalculator >
361 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
369 LocalDistanceType dist( uLocal, vLocal );
376 template<
class Gr
idPart,
class OrderCalculator >
377 template<
class Function >
383 typedef typename Dune::FieldTraits< RangeFieldType >::real_type
RealType;
387 : function_( function ),
391 template<
class Po
int >
392 void evaluate (
const Point &x, RangeType &ret )
const 395 function_.evaluate( x, phi );
396 ret = std :: pow ( phi.two_norm(),
p_);
400 const FunctionType &function_;
405 template<
class Gr
idPart,
class OrderCalculator >
406 template<
class UFunction,
class VFunction >
413 typedef typename Dune::FieldTraits< RangeFieldType >::real_type
RealType;
421 template<
class Po
int >
422 void evaluate (
const Point &x, RangeType &ret )
const 425 u_.evaluate( x, ret );
426 v_.evaluate( x, phi );
430 template<
class Po
int >
431 void jacobian (
const Point &x, JacobianRangeType &ret )
const 433 JacobianRangeType phi;
434 u_.jacobian( x, ret );
435 v_.jacobian( x, phi );
440 const UFunctionType &u_;
441 const VFunctionType &v_;
448 template<
class WeightFunction,
class OrderCalculator >
452 weightFunction_( weightFunction ),
455 static_assert( (WeightFunctionSpaceType::dimRange == 1),
456 "Weight function must be scalar." );
460 template<
class WeightFunction,
class OrderCalculator >
461 template<
class LocalFunctionType,
class ReturnType >
476 template<
class WeightFunction,
class OrderCalculator >
477 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
481 typedef typename BaseType::template FunctionDistance< ULocalFunctionType, VLocalFunctionType > LocalDistanceType;
488 LocalDistanceType dist( uLocal, vLocal );
495 template<
class WeightFunction,
class OrderCalculator>
496 template<
class Function >
502 typedef typename Dune::FieldTraits< RangeFieldType >::real_type
RealType;
506 const FunctionType &
function,
508 : weightFunction_( weightFunction ),
509 function_( function ),
513 template<
class Po
int >
514 void evaluate (
const Point &x, RangeType &ret )
const 517 weightFunction_.evaluate( x, weight );
520 function_.evaluate( x, phi );
521 ret = weight[ 0 ] * std::pow ( phi.two_norm(),
p_);
526 const FunctionType &function_;
534 #endif // #ifndef DUNE_FEM_LPNORM_HH BaseType::EntityType EntityType
Definition: lpnorm.hh:216
WeightFunctionType::LocalFunctionType LocalWeightFunctionType
Definition: lpnorm.hh:268
Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u) const
Definition: lpnorm.hh:318
Definition: lpnorm.hh:213
int order(const int spaceOrder) const
Definition: lpnorm.hh:309
Definition: lpnorm.hh:253
UFunctionType::RangeType RangeType
Definition: lpnorm.hh:414
LPNorm(const GridPartType &gridPart, const double p)
Definition: lpnorm.hh:301
static ReturnType forEach(const ThisType &norm, const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const ReturnType &initialValue, unsigned int order)
Definition: lpnorm.hh:46
Definition: lpnorm.hh:268
static constexpr T max(T a)
Definition: utility.hh:65
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: lpnorm.hh:479
LPNormBase(const GridPartType &gridPart)
Definition: lpnorm.hh:138
static ReturnType forEach(const ThisType &norm, const F &f, const VDiscreteFunctionType &v, const ReturnType &initialValue, const unsigned int order)
Definition: lpnorm.hh:73
GridPart GridPartType
Definition: lpnorm.hh:206
static constexpr T sum(T a)
Definition: utility.hh:33
GridPart GridPartType
Definition: lpnorm.hh:35
FunctionType::RangeFieldType RangeFieldType
Definition: lpnorm.hh:501
Function FunctionType
Definition: lpnorm.hh:380
const NormImplementation & asImp() const
Definition: bartonnackmaninterface.hh:37
WeightFunctionSpaceType::GridPartType GridPartType
Definition: lpnorm.hh:264
OrderCalculator * orderCalc_
Definition: lpnorm.hh:243
void evaluate(const Point &x, RangeType &ret) const
Definition: lpnorm.hh:422
static double max(const Double &v, const double p)
Definition: double.hh:387
FunctionType::RangeFieldType RangeFieldType
Definition: lpnorm.hh:382
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: lpnorm.hh:383
FunctionSpaceType::RangeType RangeType
range type
Definition: function.hh:68
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: lpnorm.hh:413
BaseType::IntegratorType IntegratorType
Definition: lpnorm.hh:274
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:71
BaseType::EntityType EntityType
Definition: lpnorm.hh:273
ReturnType forEach(const DiscreteFunctionType &u, const ReturnType &initialValue, unsigned int order=0) const
Definition: lpnorm.hh:106
FunctionDistance(const UFunctionType &u, const VFunctionType &v)
Definition: lpnorm.hh:417
WeightedLPNorm(const WeightFunctionType &weightFunction, const double p)
Definition: lpnorm.hh:450
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: lpnorm.hh:40
FieldVector< RealType, 1 > RangeType
Definition: lpnorm.hh:503
void jacobian(const Point &x, JacobianRangeType &ret) const
Definition: lpnorm.hh:431
Definition: coordinate.hh:4
WeightedFunctionMultiplicator(const LocalWeightFunctionType &weightFunction, const FunctionType &function, double p)
Definition: lpnorm.hh:505
UFunctionType::RangeFieldType RangeFieldType
Definition: lpnorm.hh:412
FunctionMultiplicator(const FunctionType &function, double p)
Definition: lpnorm.hh:386
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: lpnorm.hh:144
integrator for arbitrary functions providing evaluate
Definition: integrator.hh:28
Definition: lpnorm.hh:200
Definition: lpnorm.hh:216
void evaluate(const Point &x, RangeType &ret) const
Definition: lpnorm.hh:392
FieldVector< RealType, 1 > RangeType
Definition: lpnorm.hh:384
default Quadrature Order Calculator
Definition: lpnorm.hh:186
static ReturnType forEach(const ThisType &norm, const UDiscreteFunctionType &u, const F &f, const ReturnType &initialValue, const unsigned int order)
Definition: lpnorm.hh:94
VFunction VFunctionType
Definition: lpnorm.hh:410
ReturnType forEach(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const ReturnType &initialValue, unsigned int order=0) const
Definition: lpnorm.hh:126
UFunction UFunctionType
Definition: lpnorm.hh:409
GridFunctionAdapter provides local functions for a Function.
Definition: gridfunctionadapter.hh:36
Function FunctionType
Definition: lpnorm.hh:499
UFunctionType::JacobianRangeType JacobianRangeType
Definition: lpnorm.hh:415
CachingQuadrature< GridPartType, 0 > QuadratureType
Definition: lpnorm.hh:219
GridPartType::CollectiveCommunicationType comm() const
Definition: lpnorm.hh:157
WeightFunction WeightFunctionType
Definition: lpnorm.hh:261
double p_
Definition: lpnorm.hh:242
FunctionSpaceType::RangeFieldType RangeFieldType
field type of range
Definition: function.hh:64
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: lpnorm.hh:502
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: lpnorm.hh:363
Abstract class representing a function.
Definition: function.hh:43
WeightFunctionType::RangeType WeightType
Definition: lpnorm.hh:271
Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v) const
Definition: lpnorm.hh:335
Quadrature Order Interface.
Definition: lpnorm.hh:177
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: lpnorm.hh:463
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: lpnorm.hh:351
void evaluate(const Point &x, RangeType &ret) const
Definition: lpnorm.hh:514
Integrator< QuadratureType > IntegratorType
Definition: lpnorm.hh:220
const GridPartType & gridPart() const
Definition: lpnorm.hh:155
WeightFunctionType::DiscreteFunctionSpaceType WeightFunctionSpaceType
Definition: lpnorm.hh:263
Definition: bartonnackmaninterface.hh:15
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: lpnorm.hh:150