dune-fem  2.4.1-rc
lpnorm.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_LPNORM_HH
2 #define DUNE_FEM_LPNORM_HH
3 
4 #include <dune/common/typetraits.hh>
5 
6 #include <dune/grid/common/rangegenerators.hh>
7 
10 
12 
13 #include <dune/common/bartonnackmanifcheck.hh>
15 
16 namespace Dune
17 {
18 
19  namespace Fem
20  {
21 
22  // LPNormBase
23  // ----------
24 
25  template< class GridPart, class NormImplementation >
26  class LPNormBase
27  : public BartonNackmanInterface< LPNormBase< GridPart, NormImplementation >,
28  NormImplementation >
29  {
31  NormImplementation > BaseType ;
33 
34  public:
35  typedef GridPart GridPartType;
36 
37  protected:
38  using BaseType :: asImp ;
39 
40  typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
41 
42  template <bool uDiscrete, bool vDiscrete>
44  {
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 )
47  {
48  static_assert( uDiscrete && vDiscrete, "Distance can only be calculated between GridFunctions" );
49 
50  order = (order == 0 ? 2*std::max( u.space().order(), v.space().order() ) : order);
51 
52  ReturnType sum( 0 );
53  const auto end = norm.gridPart_.template end<0>();
54  for (auto it = norm.gridPart_.template begin<0>(); it!=end; ++it)
55  {
56  const auto entity = *it;
57  const typename UDiscreteFunctionType::LocalFunctionType uLocal = u.localFunction( entity );
58  const typename VDiscreteFunctionType::LocalFunctionType vLocal = v.localFunction( entity );
59  norm.distanceLocal( entity, order, uLocal, vLocal, sum );
60  }
61  return sum;
62  }
63  };
64 
65  // this specialization creates a grid function adapter
66  template <bool vDiscrete>
67  struct ForEachCaller<false, vDiscrete>
68  {
69  template <class F,
70  class VDiscreteFunctionType,
71  class ReturnType>
72  static
73  ReturnType forEach ( const ThisType& norm,
74  const F& f, const VDiscreteFunctionType& v,
75  const ReturnType& initialValue,
76  const unsigned int order )
77  {
78  typedef GridFunctionAdapter< F, GridPartType> GridFunction ;
79  GridFunction u( "LPNorm::adapter" , f , v.space().gridPart(), v.space().order() );
80 
82  forEach( norm, u, v, initialValue, order );
83  }
84  };
85 
86  // this specialization simply switches arguments
87  template <bool uDiscrete>
88  struct ForEachCaller<uDiscrete, false>
89  {
90  template <class UDiscreteFunctionType,
91  class F,
92  class ReturnType>
93  static
94  ReturnType forEach ( const ThisType& norm,
95  const UDiscreteFunctionType& u,
96  const F& f,
97  const ReturnType& initialValue,
98  const unsigned int order )
99  {
101  forEach( norm, f, u, initialValue, order );
102  }
103  };
104 
105  template< class DiscreteFunctionType, class ReturnType >
106  ReturnType forEach ( const DiscreteFunctionType &u, const ReturnType &initialValue, unsigned int order = 0 ) const
107  {
108  static_assert( (IsBaseOf<Fem::HasLocalFunction, DiscreteFunctionType>::value),
109  "Norm only implemented for quantities implementing a local function!" );
110 
111  order = (order == 0 ? 2*u.space().order() : order);
112 
113  ReturnType sum( 0 );
114  const auto end = gridPart_.template end<0>();
115  for (auto it = gridPart_.template begin<0>(); it!=end; ++it)
116  // for( const EntityType &entity : elements( gridPart_, Partitions::interior ) )
117  {
118  const EntityType &entity = *it;
119  typename DiscreteFunctionType::LocalFunctionType uLocal = u.localFunction( entity );
120  normLocal( entity, order, uLocal, sum );
121  }
122  return sum;
123  }
124 
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
127  {
128  enum { uDiscrete = std::is_convertible<UDiscreteFunctionType, HasLocalFunction>::value };
129  enum { vDiscrete = std::is_convertible<VDiscreteFunctionType, HasLocalFunction>::value };
130 
131  // call forEach depending on which argument is a grid function,
132  // i.e. has a local function
134  forEach( *this, u, v, initialValue, order );
135  }
136 
137  public:
138  explicit LPNormBase ( const GridPartType &gridPart )
139  : gridPart_( gridPart )
140  {}
141 
142  protected:
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
145  {
146  CHECK_AND_CALL_INTERFACE_IMPLEMENTATION( asImp().distanceLocal( entity, order, uLocal, vLocal, sum ) );
147  }
148 
149  template< class LocalFunctionType, class ReturnType >
150  void normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const
151  {
152  CHECK_AND_CALL_INTERFACE_IMPLEMENTATION( asImp().normLocal( entity, order, uLocal, sum ) );
153  }
154 
155  const GridPartType &gridPart () const { return gridPart_; }
156 
157  typename GridPartType::CollectiveCommunicationType comm () const
158  {
159  return gridPart().comm();
160  }
161 
162  private:
163  const GridPartType &gridPart_;
164  };
165 
166 
167 
168  // TODO weighte LP norm might be adapted later
169  // LPNorm
170  //
171  // !!!!! It is not cleared which quadrature order have to be applied for p > 2!!!!
172  // !!!!! For p = 1 this norm does not work !!!!
173  // ------
174 
175 
178  {
179  virtual int operator() (const double p)=0;
180  };
181 
183  // can be re-implemented in order to use
184  // a different type of calculation
185  // which can be sepcified in the second template argument of LPNorm
187  {
188  int operator() (const double p)
189  {
190  int ret=0;
191  const double q = p / (p-1);
192  double max = std::max(p,q);
193  assert(max < std::numeric_limits<int>::max()/2. );
194  ret = max +1;
195  return ret;
196  }
197  };
198 
199  template< class GridPart, class OrderCalculator = DefaultOrderCalculator >
200  class LPNorm : public LPNormBase < GridPart, LPNorm< GridPart, OrderCalculator > >
201  {
204 
205  public:
206  typedef GridPart GridPartType;
207 
208  using BaseType::gridPart;
209  using BaseType::comm;
210 
211  protected:
212  template< class Function >
214 
215  template< class UFunction, class VFunction >
216  struct FunctionDistance;
217 
218  typedef typename BaseType::EntityType EntityType;
221 
222  public:
223  explicit LPNorm ( const GridPartType &gridPart, const double p );
224 
225  template< class DiscreteFunctionType >
226  typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
227  norm ( const DiscreteFunctionType &u ) const;
228 
229  template< class UDiscreteFunctionType, class VDiscreteFunctionType >
230  typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
231  distance ( const UDiscreteFunctionType &u, const VDiscreteFunctionType &v ) const;
232 
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;
235 
236  template< class LocalFunctionType, class ReturnType >
237  void normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const;
238 
239  int order ( const int spaceOrder ) const ;
240 
241  protected:
242  double p_ ;
243  OrderCalculator *orderCalc_;
244  };
245 
246 
247 
248 
249  // WeightedLPNorm
250  // --------------
251 
252  template< class WeightFunction, class OrderCalculator = DefaultOrderCalculator >
254  : public LPNorm< typename WeightFunction::DiscreteFunctionSpaceType::GridPartType,
255  OrderCalculator >
256  {
259 
260  public:
261  typedef WeightFunction WeightFunctionType;
262 
263  typedef typename WeightFunctionType::DiscreteFunctionSpaceType WeightFunctionSpaceType;
264  typedef typename WeightFunctionSpaceType::GridPartType GridPartType;
265 
266  protected:
267  template< class Function >
269 
270  typedef typename WeightFunctionType::LocalFunctionType LocalWeightFunctionType;
271  typedef typename WeightFunctionType::RangeType WeightType;
272 
275 
276  using BaseType::gridPart;
277  using BaseType::comm;
278 
279  public:
280  using BaseType::norm;
281  using BaseType::distance;
282 
283  explicit WeightedLPNorm ( const WeightFunctionType &weightFunction, const double p );
284 
285  template< class LocalFunctionType, class ReturnType >
286  void normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const;
287 
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;
290 
291  private:
292  const WeightFunctionType &weightFunction_;
293  const double p_;
294  };
295 
296 
297  // Implementation of LPNorm
298  // ------------------------
299 
300  template< class GridPart, class OrderCalculator >
302  : BaseType( gridPart ),
303  p_(p),
304  orderCalc_( new OrderCalculator() )
305  {
306  }
307 
308  template< class GridPart, class OrderCalculator>
309  inline int LPNorm< GridPart, OrderCalculator>::order(const int spaceOrder) const
310  {
311  return spaceOrder * orderCalc_->operator() (p_);
312  }
313 
314 
315  template< class GridPart, class OrderCalculator>
316  template< class DiscreteFunctionType >
317  inline typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
318  LPNorm< GridPart, OrderCalculator >::norm ( const DiscreteFunctionType &u ) const
319  {
320  typedef typename DiscreteFunctionType::RangeFieldType RangeFieldType;
321  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
322  typedef FieldVector< RealType, 1 > ReturnType ;
323 
324  // calculate integral over each element
325  ReturnType sum = BaseType :: forEach( u, ReturnType(0) );
326 
327  // return result
328  return std::pow ( comm().sum( sum[ 0 ] ), (1.0 / p_) );
329  }
330 
331  template< class GridPart, 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
336  {
337  typedef typename UDiscreteFunctionType::RangeFieldType RangeFieldType;
338  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
339  typedef FieldVector< RealType, 1 > ReturnType ;
340 
341  // calculate integral over each element
342  ReturnType sum = BaseType :: forEach( u, v, ReturnType(0) );
343 
344  // return result
345  return std::pow( comm().sum( sum[ 0 ] ), (1.0/p_) );
346  }
347 
348  template< class GridPart, class OrderCalculator >
349  template< class LocalFunctionType, class ReturnType >
350  inline void
351  LPNorm< GridPart, OrderCalculator >::normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const
352  {
353  IntegratorType integrator( order );
354 
356 
357  integrator.integrateAdd( entity, uLocalp, sum );
358  }
359 
360  template< class GridPart, class OrderCalculator >
361  template< class ULocalFunctionType, class VLocalFunctionType, class ReturnType >
362  inline void
363  LPNorm< GridPart, OrderCalculator >::distanceLocal ( const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum ) const
364  {
366 
367  IntegratorType integrator( order );
368 
369  LocalDistanceType dist( uLocal, vLocal );
371 
372  integrator.integrateAdd( entity, distp, sum );
373  }
374 
375 
376  template< class GridPart, class OrderCalculator >
377  template< class Function >
378  struct LPNorm< GridPart, OrderCalculator >::FunctionMultiplicator
379  {
381 
383  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
384  typedef FieldVector< RealType, 1 > RangeType ;
385 
386  explicit FunctionMultiplicator ( const FunctionType &function, double p )
387  : function_( function ),
388  p_(p)
389  {}
390 
391  template< class Point >
392  void evaluate ( const Point &x, RangeType &ret ) const
393  {
394  typename FunctionType::RangeType phi;
395  function_.evaluate( x, phi );
396  ret = std :: pow ( phi.two_norm(), p_);
397  }
398 
399  private:
400  const FunctionType &function_;
401  double p_;
402  };
403 
404 
405  template< class GridPart, class OrderCalculator >
406  template< class UFunction, class VFunction >
407  struct LPNorm< GridPart, OrderCalculator >::FunctionDistance
408  {
409  typedef UFunction UFunctionType;
410  typedef VFunction VFunctionType;
411 
412  typedef typename UFunctionType::RangeFieldType RangeFieldType;
413  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
414  typedef typename UFunctionType::RangeType RangeType;
415  typedef typename UFunctionType::JacobianRangeType JacobianRangeType;
416 
417  FunctionDistance ( const UFunctionType &u, const VFunctionType &v )
418  : u_( u ), v_( v )
419  {}
420 
421  template< class Point >
422  void evaluate ( const Point &x, RangeType &ret ) const
423  {
424  RangeType phi;
425  u_.evaluate( x, ret );
426  v_.evaluate( x, phi );
427  ret -= phi;
428  }
429 
430  template< class Point >
431  void jacobian ( const Point &x, JacobianRangeType &ret ) const
432  {
433  JacobianRangeType phi;
434  u_.jacobian( x, ret );
435  v_.jacobian( x, phi );
436  ret -= phi;
437  }
438 
439  private:
440  const UFunctionType &u_;
441  const VFunctionType &v_;
442  };
443 
444 
445  // Implementation of WeightedL2Norm
446  // --------------------------------
447 
448  template< class WeightFunction, class OrderCalculator >
450  ::WeightedLPNorm ( const WeightFunctionType &weightFunction, double p )
451  : BaseType( weightFunction.space().gridPart(), p ),
452  weightFunction_( weightFunction ),
453  p_(p)
454  {
455  static_assert( (WeightFunctionSpaceType::dimRange == 1),
456  "Weight function must be scalar." );
457  }
458 
459 
460  template< class WeightFunction, class OrderCalculator >
461  template< class LocalFunctionType, class ReturnType >
462  inline void
463  WeightedLPNorm< WeightFunction, OrderCalculator >::normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const
464  {
465  // !!!! order !!!!
466  IntegratorType integrator( order );
467 
468  LocalWeightFunctionType wfLocal = weightFunction_.localFunction( entity );
469 
470  WeightedFunctionMultiplicator< LocalFunctionType > uLocal2( wfLocal, uLocal, p_ );
471 
472  integrator.integrateAdd( entity, uLocal2, sum );
473  }
474 
475 
476  template< class WeightFunction,class OrderCalculator >
477  template< class ULocalFunctionType, class VLocalFunctionType, class ReturnType >
478  inline void
479  WeightedLPNorm< WeightFunction, OrderCalculator >::distanceLocal ( const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum ) const
480  {
481  typedef typename BaseType::template FunctionDistance< ULocalFunctionType, VLocalFunctionType > LocalDistanceType;
482 
483  // !!!! order !!!!
484  IntegratorType integrator( order );
485 
486  LocalWeightFunctionType wfLocal = weightFunction_.localFunction( entity );
487 
488  LocalDistanceType dist( uLocal, vLocal );
490 
491  integrator.integrateAdd( entity, dist2, sum );
492  }
493 
494 
495  template< class WeightFunction, class OrderCalculator>
496  template< class Function >
497  struct WeightedLPNorm< WeightFunction, OrderCalculator>::WeightedFunctionMultiplicator
498  {
500 
502  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
503  typedef FieldVector< RealType, 1 > RangeType;
504 
506  const FunctionType &function,
507  double p )
508  : weightFunction_( weightFunction ),
509  function_( function ),
510  p_(p)
511  {}
512 
513  template< class Point >
514  void evaluate ( const Point &x, RangeType &ret ) const
515  {
516  WeightType weight;
517  weightFunction_.evaluate( x, weight );
518 
519  typename FunctionType::RangeType phi;
520  function_.evaluate( x, phi );
521  ret = weight[ 0 ] * std::pow ( phi.two_norm(), p_);
522  }
523 
524  private:
525  const LocalWeightFunctionType &weightFunction_;
526  const FunctionType &function_;
527  double p_;
528  };
529 
530  } // namespace Fem
531 
532 } // namespace Dune
533 
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
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
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
Definition: lpnorm.hh:26
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
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
Definition: lpnorm.hh:43
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: lpnorm.hh:150