dune-fem  2.8-git
l1norm.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_L1NORM_HH
2 #define DUNE_FEM_L1NORM_HH
3 
5 
7 
8 namespace Dune
9 {
10 
11  namespace Fem
12  {
13 
14  // L1Norm
15  // ------
16 
17  template< class GridPart >
18  class L1Norm : public IntegralBase< GridPart, L1Norm< GridPart > >
19  {
22 
23  public:
24  typedef GridPart GridPartType;
25 
26  using BaseType :: gridPart ;
27  using BaseType :: comm ;
28 
29  protected:
30  template< class Function >
31  struct FunctionAbs;
32 
33  template< class UFunction, class VFunction >
34  struct FunctionDistance;
35 
36  typedef typename BaseType::EntityType EntityType;
38 
39  const unsigned int order_;
40  const bool communicate_;
41  public:
47  explicit L1Norm ( const GridPartType &gridPart,
48  const unsigned int order = 0,
49  const bool communicate = true );
50 
51 
53  template< class DiscreteFunctionType, class PartitionSet >
54  typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
55  norm ( const DiscreteFunctionType &u, const PartitionSet& partitionSet ) const;
56 
58  template< class DiscreteFunctionType >
59  typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
60  norm ( const DiscreteFunctionType &u ) const
61  {
62  return norm( u, Partitions::interior );
63  }
64 
66  template< class UDiscreteFunctionType, class VDiscreteFunctionType, class PartitionSet >
67  typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
68  distance ( const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const PartitionSet& partitionSet ) const;
69 
71  template< class UDiscreteFunctionType, class VDiscreteFunctionType >
72  typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
73  distance ( const UDiscreteFunctionType &u, const VDiscreteFunctionType &v ) const
74  {
75  return distance( u, v, Partitions::interior );
76  }
77 
78  template< class LocalFunctionType, class ReturnType >
79  void normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const;
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 
85 
86 
87  // Implementation of L1Norm
88  // ------------------------
89 
90  template< class GridPart >
91  inline L1Norm< GridPart >::L1Norm ( const GridPartType &gridPart, const unsigned int order, const bool communicate )
92  : BaseType( gridPart ),
93  order_( order ),
94  communicate_( BaseType::checkCommunicateFlag( communicate ) )
95  {}
96 
97 
98  template< class GridPart >
99  template< class DiscreteFunctionType, class PartitionSet >
100  inline typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
101  L1Norm< GridPart >::norm ( const DiscreteFunctionType &u, const PartitionSet& partitionSet ) const
102  {
103  typedef typename DiscreteFunctionType::RangeFieldType RangeFieldType;
104  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
105  typedef FieldVector< RealType, 1 > ReturnType ;
106 
107  // calculate integral over each element
108  ReturnType sum = BaseType :: forEach( u, ReturnType(0), partitionSet, order_ );
109 
110  // communicate_ indicates global norm
111  if( communicate_ )
112  {
113  sum[ 0 ] = comm().sum( sum[ 0 ] );
114  }
115 
116  return sum[ 0 ];
117  }
118 
119 
120  template< class GridPart >
121  template< class UDiscreteFunctionType, class VDiscreteFunctionType, class PartitionSet >
122  inline typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
124  ::distance ( const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const PartitionSet& partitionSet ) const
125  {
126  typedef typename UDiscreteFunctionType::RangeFieldType RangeFieldType;
127  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
128  typedef FieldVector< RealType, 1 > ReturnType ;
129 
130  // calculate integral over each element
131  ReturnType sum = BaseType :: forEach( u, v, ReturnType(0), partitionSet, order_ );
132 
133  // communicate_ indicates global norm
134  if( communicate_ )
135  {
136  sum[ 0 ] = comm().sum( sum[ 0 ] );
137  }
138 
139  return sum[ 0 ];
140  }
141 
142  template< class GridPart >
143  template< class LocalFunctionType, class ReturnType >
144  inline void
145  L1Norm< GridPart >::normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const
146  {
147  Integrator< QuadratureType > integrator( order );
148 
149  FunctionAbs< LocalFunctionType > uLocalAbs( uLocal );
150 
151  integrator.integrateAdd( entity, uLocalAbs, sum );
152  }
153 
154  template< class GridPart >
155  template< class ULocalFunctionType, class VLocalFunctionType, class ReturnType >
156  inline void
157  L1Norm< GridPart >::distanceLocal ( const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum ) const
158  {
159  Integrator< QuadratureType > integrator( order );
160 
162 
163  LocalDistanceType dist( uLocal, vLocal );
164  FunctionAbs< LocalDistanceType > distAbs( dist );
165 
166  integrator.integrateAdd( entity, distAbs, sum );
167  }
168 
169 
170  template< class GridPart >
171  template< class Function >
172  struct L1Norm< GridPart >::FunctionAbs
173  {
175 
177  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
178  typedef FieldVector< RealType, 1 > RangeType;
179 
180  explicit FunctionAbs ( const FunctionType &function )
181  : function_( function )
182  {}
183 
184  template< class Point >
185  void evaluate ( const Point &x, RangeType &ret ) const
186  {
187  typename FunctionType::RangeType phi;
188  function_.evaluate( x, phi );
189  ret = phi.one_norm();
190  }
191 
192  private:
193  const FunctionType &function_;
194  };
195 
196 
197  template< class GridPart >
198  template< class UFunction, class VFunction >
199  struct L1Norm< GridPart >::FunctionDistance
200  {
201  typedef UFunction UFunctionType;
202  typedef VFunction VFunctionType;
203 
204  typedef typename UFunctionType::RangeFieldType RangeFieldType;
205  typedef typename UFunctionType::RangeType RangeType;
206  typedef typename UFunctionType::JacobianRangeType JacobianRangeType;
207 
209  : u_( u ), v_( v )
210  {}
211 
212  template< class Point >
213  void evaluate ( const Point &x, RangeType &ret ) const
214  {
215  RangeType phi;
216  u_.evaluate( x, ret );
217  v_.evaluate( x, phi );
218  ret -= phi;
219  }
220 
221  template< class Point >
222  void jacobian ( const Point &x, JacobianRangeType &ret ) const
223  {
224  JacobianRangeType phi;
225  u_.jacobian( x, ret );
226  v_.jacobian( x, phi );
227  ret -= phi;
228  }
229 
230  private:
231  const UFunctionType &u_;
232  const VFunctionType &v_;
233  };
234 
235  } // namespace Fem
236 
237 } // namespace Dune
238 
239 #endif // #ifndef DUNE_FEM_L1NORM_HH
Definition: bindguard.hh:11
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: l1norm.hh:19
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: l1norm.hh:157
BaseType::EntityType EntityType
Definition: l1norm.hh:34
Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u) const
|| u ||_L1 on interior partition entities
Definition: l1norm.hh:60
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: l1norm.hh:124
GridPart GridPartType
Definition: l1norm.hh:24
Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v) const
|| u - v ||_L2 on interior partition entities
Definition: l1norm.hh:73
Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u, const PartitionSet &partitionSet) const
|| u ||_L1 on given set of entities (partition set)
Definition: l1norm.hh:101
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: l1norm.hh:145
const bool communicate_
Definition: l1norm.hh:40
L1Norm(const GridPartType &gridPart, const unsigned int order=0, const bool communicate=true)
constructor
Definition: l1norm.hh:91
CachingQuadrature< GridPartType, 0 > QuadratureType
Definition: l1norm.hh:37
const unsigned int order_
Definition: l1norm.hh:39
Definition: l1norm.hh:173
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: l1norm.hh:177
FunctionType::RangeFieldType RangeFieldType
Definition: l1norm.hh:176
FieldVector< RealType, 1 > RangeType
Definition: l1norm.hh:178
FunctionAbs(const FunctionType &function)
Definition: l1norm.hh:180
void evaluate(const Point &x, RangeType &ret) const
Definition: l1norm.hh:185
Function FunctionType
Definition: l1norm.hh:174
Definition: l1norm.hh:200
FunctionDistance(const UFunctionType &u, const VFunctionType &v)
Definition: l1norm.hh:208
VFunction VFunctionType
Definition: l1norm.hh:202
UFunction UFunctionType
Definition: l1norm.hh:201
void jacobian(const Point &x, JacobianRangeType &ret) const
Definition: l1norm.hh:222
UFunctionType::JacobianRangeType JacobianRangeType
Definition: l1norm.hh:206
void evaluate(const Point &x, RangeType &ret) const
Definition: l1norm.hh:213
UFunctionType::RangeFieldType RangeFieldType
Definition: l1norm.hh:204
UFunctionType::RangeType RangeType
Definition: l1norm.hh:205
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