dune-fem  2.8-git
l2norm.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_L2NORM_HH
2 #define DUNE_FEM_L2NORM_HH
3 
6 
8 
9 namespace Dune
10 {
11 
12  namespace Fem
13  {
14 
15  // L2Norm
16  // ------
17 
18  template< class GridPart >
19  class L2Norm : public IntegralBase< GridPart, L2Norm< GridPart > >
20  {
23 
24  public:
25  typedef GridPart GridPartType;
26 
27  using BaseType :: gridPart ;
28  using BaseType :: comm ;
29 
30  template< class Function >
31  struct FunctionSquare;
32 
33  template< class UFunction, class VFunction >
34  struct FunctionDistance;
35 
36  protected:
37  typedef typename BaseType::EntityType EntityType;
40 
41  const unsigned int order_;
42  const bool communicate_;
43  public:
49  explicit L2Norm ( const GridPartType &gridPart,
50  const unsigned int order = 0,
51  const bool communicate = true );
52 
54  template< class DiscreteFunctionType, class PartitionSet >
55  typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
56  norm ( const DiscreteFunctionType &u, const PartitionSet& partitionSet ) const;
57 
59  template< class DiscreteFunctionType >
60  typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
61  norm ( const DiscreteFunctionType &u ) const
62  {
63  return norm( u, Partitions::interior );
64  }
65 
67  template< class UDiscreteFunctionType, class VDiscreteFunctionType, class PartitionSet >
68  typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
69  distance ( const UDiscreteFunctionType &u,
70  const VDiscreteFunctionType &v,
71  const PartitionSet& partitionSet ) const;
72 
74  template< class UDiscreteFunctionType, class VDiscreteFunctionType >
75  typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
76  distance ( const UDiscreteFunctionType &u, const VDiscreteFunctionType &v ) const
77  {
78  return distance( u, v, Partitions::interior );
79  }
80 
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;
83 
84  template< class LocalFunctionType, class ReturnType >
85  void normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const;
86  };
87 
88 
89  // Implementation of L2Norm
90  // ------------------------
91 
92  template< class GridPart >
93  inline L2Norm< GridPart >::L2Norm ( const GridPartType &gridPart, const unsigned int order, const bool communicate )
94  : BaseType( gridPart ),
95  order_( order ),
96  communicate_( BaseType::checkCommunicateFlag( communicate ) )
97  {
98  }
99 
100 
101  template< class GridPart >
102  template< class DiscreteFunctionType, class PartitionSet >
103  typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
104  L2Norm< GridPart >::norm ( const DiscreteFunctionType &u, const PartitionSet& partitionSet ) const
105  {
106  typedef typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type RealType;
107  typedef FieldVector< RealType, 1 > ReturnType ;
108 
109  // calculate integral over each element
110  ReturnType sum = BaseType :: forEach( u, ReturnType(0), partitionSet, order_ );
111 
112  // communicate_ indicates global norm
113  if( communicate_ )
114  {
115  sum[ 0 ] = comm().sum( sum[ 0 ] );
116  }
117 
118  // return result, e.g. sqrt of calculated sum
119  return std::sqrt( sum[ 0 ] );
120  }
121 
122 
123  template< class GridPart >
124  template< class UDiscreteFunctionType, class VDiscreteFunctionType, class PartitionSet >
125  typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
127  ::distance ( const UDiscreteFunctionType &u,
128  const VDiscreteFunctionType &v,
129  const PartitionSet& partitionSet ) const
130  {
131  typedef typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type RealType;
132  typedef FieldVector< RealType, 1 > ReturnType ;
133 
134  // calculate integral over each element
135  ReturnType sum = BaseType :: forEach( u, v, ReturnType(0), partitionSet, order_ );
136 
137  // communicate_ indicates global norm
138  if( communicate_ )
139  {
140  sum[ 0 ] = comm().sum( sum[ 0 ] );
141  }
142 
143  // return result, e.g. sqrt of calculated sum
144  return std::sqrt( sum[ 0 ] );
145  }
146 
147  template< class GridPart >
148  template< class LocalFunctionType, class ReturnType >
149  inline void L2Norm< GridPart >::normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const
150  {
151  // evaluate norm locally
152  IntegratorType integrator( order );
153 
154  FunctionSquare< LocalFunctionType > uLocal2( uLocal );
155 
156  integrator.integrateAdd( entity, uLocal2, sum );
157  }
158 
159  template< class GridPart >
160  template< class ULocalFunctionType, class VLocalFunctionType, class ReturnType >
161  inline void
162  L2Norm< GridPart >::distanceLocal ( const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum ) const
163  {
164  // evaluate norm locally
165  IntegratorType integrator( order );
166 
168 
169  LocalDistanceType dist( uLocal, vLocal );
171 
172  integrator.integrateAdd( entity, dist2, sum );
173  }
174 
175 
176  template< class GridPart >
177  template< class Function >
178  struct L2Norm< GridPart >::FunctionSquare
179  {
181 
183  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
184  typedef FieldVector< RealType, 1 > RangeType;
185 
186  explicit FunctionSquare ( const FunctionType &function )
187  : function_( function )
188  {}
189 
190  template< class Point >
191  void evaluate ( const Point &x, RangeType &ret ) const
192  {
193  typename FunctionType::RangeType phi;
194  function_.evaluate( x, phi );
195  ret = phi.two_norm2();
196  }
197 
198  private:
199  const FunctionType &function_;
200  };
201 
202 
203  template< class GridPart >
204  template< class UFunction, class VFunction >
205  struct L2Norm< GridPart >::FunctionDistance
206  {
207  typedef UFunction UFunctionType;
208  typedef VFunction VFunctionType;
209 
210  typedef typename UFunctionType::RangeFieldType RangeFieldType;
211  typedef typename UFunctionType::RangeType RangeType;
212  typedef typename UFunctionType::JacobianRangeType JacobianRangeType;
213 
215  : u_( u ), v_( v )
216  {}
217 
218  template< class Point >
219  void evaluate ( const Point &x, RangeType &ret ) const
220  {
221  RangeType phi;
222  u_.evaluate( x, ret );
223  v_.evaluate( x, phi );
224  ret -= phi;
225  }
226 
227  template< class Point >
228  void jacobian ( const Point &x, JacobianRangeType &ret ) const
229  {
230  JacobianRangeType phi;
231  u_.jacobian( x, ret );
232  v_.jacobian( x, phi );
233  ret -= phi;
234  }
235 
236  private:
237  const UFunctionType &u_;
238  const VFunctionType &v_;
239  };
240 
241 
242 
243  // WeightedL2Norm
244  // --------------
245 
246  template< class WeightFunction >
248  : public L2Norm< typename WeightFunction::DiscreteFunctionSpaceType::GridPartType >
249  {
252 
253  public:
254  typedef WeightFunction WeightFunctionType;
255 
256  typedef typename WeightFunctionType::DiscreteFunctionSpaceType WeightFunctionSpaceType;
257  typedef typename WeightFunctionSpaceType::GridPartType GridPartType;
258 
259  protected:
260  template< class Function >
261  struct WeightedFunctionSquare;
262 
264  typedef typename WeightFunctionType::RangeType WeightType;
265 
268 
269  using BaseType::gridPart;
270  using BaseType::comm;
271 
272  public:
273 
274  using BaseType::norm;
275  using BaseType::distance;
276 
277  explicit WeightedL2Norm ( const WeightFunctionType &weightFunction, const unsigned int order = 0, const bool communicate = true );
278 
279  template< class LocalFunctionType, class ReturnType >
280  void normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const;
281 
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;
284 
285  private:
286  LocalWeightFunctionType wfLocal_;
287  };
288 
289 
290 
291 
292  // Implementation of WeightedL2Norm
293  // --------------------------------
294 
295  template< class WeightFunction >
297  ::WeightedL2Norm ( const WeightFunctionType &weightFunction, const unsigned int order, const bool communicate )
298  : BaseType( weightFunction.space().gridPart(), order, communicate ),
299  wfLocal_( weightFunction )
300  {
301  static_assert( (WeightFunctionSpaceType::dimRange == 1),
302  "Wight function must be scalar." );
303  }
304 
305 
306  template< class WeightFunction >
307  template< class LocalFunctionType, class ReturnType >
308  inline void
309  WeightedL2Norm< WeightFunction >::normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const
310  {
311  // !!!! order !!!!
312  IntegratorType integrator( order );
313 
314  wfLocal_.bind( entity );
315  WeightedFunctionSquare< LocalFunctionType > uLocal2( wfLocal_, uLocal );
316  integrator.integrateAdd( entity, uLocal2, sum );
317  wfLocal_.unbind();
318  }
319 
320 
321  template< class WeightFunction >
322  template< class ULocalFunctionType, class VLocalFunctionType, class ReturnType >
323  inline void
324  WeightedL2Norm< WeightFunction >::distanceLocal ( const EntityType& entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType& sum ) const
325  {
326  typedef typename BaseType::template FunctionDistance< ULocalFunctionType, VLocalFunctionType > LocalDistanceType;
327 
328  // !!!! order !!!!
329  IntegratorType integrator( order );
330 
331  wfLocal_.bind( entity );
332  LocalDistanceType dist( uLocal, vLocal );
333  WeightedFunctionSquare< LocalDistanceType > dist2( wfLocal_, dist );
334 
335  integrator.integrateAdd( entity, dist2, sum );
336  wfLocal_.unbind();
337  }
338 
339 
340  template< class WeightFunction >
341  template< class Function >
342  struct WeightedL2Norm< WeightFunction >::WeightedFunctionSquare
343  {
345 
347  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
348  typedef FieldVector< RealType, 1 > RangeType;
349 
351  const FunctionType &function )
352  : weightFunction_( weightFunction ),
353  function_( function )
354  {}
355 
356  template< class Point >
357  void evaluate ( const Point &x, RangeType &ret ) const
358  {
359  WeightType weight;
360  weightFunction_.evaluate( x, weight );
361 
362  typename FunctionType::RangeType phi;
363  function_.evaluate( x, phi );
364  ret = weight[ 0 ] * (phi * phi);
365  }
366 
367  private:
368  const LocalWeightFunctionType &weightFunction_;
369  const FunctionType &function_;
370  };
371 
372  } // end namespace Fem
373 
374 } // end namespace Dune
375 
376 #endif // #ifndef DUNE_FEM_L2NORM_HH
Definition: bindguard.hh:11
static double sqrt(const Double &v)
Definition: double.hh:886
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:517
static void forEach(IndexRange< T, sz > range, F &&f)
Definition: hybrid.hh:129
static constexpr std::decay_t< T > sum(T a)
Definition: utility.hh:33
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:28
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: domainintegral.hh:39
const GridPartType & gridPart() const
Definition: domainintegral.hh:166
const GridPartType::CollectiveCommunicationType & comm() const
Definition: domainintegral.hh:168
Definition: l2norm.hh:20
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: l2norm.hh:149
GridPart GridPartType
Definition: l2norm.hh:25
Dune::FieldTraits< typename UDiscreteFunctionType::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
const unsigned int order_
Definition: l2norm.hh:41
Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v) const
|| u - v ||_L2 on interior partition entities
Definition: l2norm.hh:76
Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u, const PartitionSet &partitionSet) const
|| u ||_L2 on given set of entities (partition set)
Definition: l2norm.hh:104
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: l2norm.hh:162
Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u) const
|| u ||_L2 on interior partition entities
Definition: l2norm.hh:61
Integrator< QuadratureType > IntegratorType
Definition: l2norm.hh:39
BaseType::EntityType EntityType
Definition: l2norm.hh:34
const bool communicate_
Definition: l2norm.hh:42
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
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: l2norm.hh:183
FieldVector< RealType, 1 > RangeType
Definition: l2norm.hh:184
void evaluate(const Point &x, RangeType &ret) const
Definition: l2norm.hh:191
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:261
BaseType::EntityType EntityType
Definition: l2norm.hh:266
WeightedL2Norm(const WeightFunctionType &weightFunction, const unsigned int order=0, const bool communicate=true)
Definition: l2norm.hh:297
WeightFunctionSpaceType::GridPartType GridPartType
Definition: l2norm.hh:257
WeightFunctionType::RangeType WeightType
Definition: l2norm.hh:264
WeightedFunctionSquare(const LocalWeightFunctionType &weightFunction, const FunctionType &function)
Definition: l2norm.hh:350
Function FunctionType
Definition: l2norm.hh:344
FieldVector< RealType, 1 > RangeType
Definition: l2norm.hh:348
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: l2norm.hh:347
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