dune-fem 2.8.0
Loading...
Searching...
No Matches
h1norm.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_H1NORM_HH
2#define DUNE_FEM_H1NORM_HH
3
5
6namespace Dune
7{
8
9 namespace Fem
10 {
11
12 template< class GridPart >
13 class H1Norm
14 : public IntegralBase< GridPart, H1Norm< GridPart> >
15 {
18
19 public:
20 typedef GridPart GridPartType;
21
22 template< class Function >
24
25 protected:
27
29 public:
31
32
34 using BaseType::comm;
35
41 explicit H1Norm ( const GridPartType &gridPart,
42 const unsigned int order = 0,
43 const bool communicate = true );
44
46 template< class DiscreteFunctionType, class PartitionSet >
47 typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
48 norm ( const DiscreteFunctionType &u, const PartitionSet& partitionSet ) const;
49
51 template< class DiscreteFunctionType >
52 typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
53 norm ( const DiscreteFunctionType &u ) const
54 {
55 return norm( u, Partitions::interior );
56 }
57
59 template< class UDiscreteFunctionType, class VDiscreteFunctionType, class PartitionSet >
60 typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
61 distance ( const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const PartitionSet& partitionSet ) const;
62
64 template< class UDiscreteFunctionType, class VDiscreteFunctionType >
65 typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
66 distance ( const UDiscreteFunctionType &u, const VDiscreteFunctionType &v ) const
67 {
68 return distance( u, v, Partitions::interior );
69 }
70
71 template< class LocalFunctionType, class ReturnType >
72 void normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const;
73
74 template< class ULocalFunctionType, class VLocalFunctionType, class ReturnType >
75 void distanceLocal ( const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum ) const;
76
77 ThisType operator= ( const ThisType& ) = delete;
78
79 protected:
80 const unsigned int order_;
81 const bool communicate_;
82 };
83
84
85
86 // H1Norm::FunctionJacobianSquare
87 // ------------------------------
88
89 template< class GridPart >
90 template< class Function >
91 struct H1Norm< GridPart >::FunctionJacobianSquare
92 {
94
96 typedef typename Dune::FieldTraits< typename FunctionType::RangeFieldType >::real_type RealType;
97 typedef FieldVector< RealType, 1 > RangeType;
98
99 public:
100 explicit FunctionJacobianSquare ( const FunctionType &function )
101 : function_( function )
102 {}
103
104 template< class Point >
105 void evaluate ( const Point &x, RangeType &ret ) const
106 {
107 const int dimRange = FunctionType::RangeType::dimension;
108
109 typename FunctionType::RangeType phi;
110 function_.evaluate( x, phi );
111 ret[ 0 ] = phi.two_norm2();
112
114 function_.jacobian( x, grad );
115 for( int i = 0; i < dimRange; ++i )
116 ret[ 0 ] += grad[ i ].two_norm2();
117 }
118
119 private:
120 const FunctionType &function_;
121 };
122
123
124
125 // Implementation of H1 Norm
126 // -------------------------
127
128 template< class GridPart >
129 inline H1Norm< GridPart >::H1Norm ( const GridPartType &gridPart, const unsigned int order, const bool communicate )
130 : BaseType( gridPart ),
131 order_( order ),
132 communicate_( communicate )
133 {}
134
135
136 template< class GridPart >
137 template< class DiscreteFunctionType, class PartitionSet >
138 inline typename Dune::FieldTraits< typename DiscreteFunctionType::RangeFieldType >::real_type
139 H1Norm< GridPart >::norm ( const DiscreteFunctionType &u, const PartitionSet& partitionSet ) const
140 {
141 typedef typename DiscreteFunctionType::RangeFieldType RangeFieldType;
142 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
143 typedef FieldVector< RealType, 1 > ReturnType ;
144
145 ReturnType sum = BaseType :: forEach( u, ReturnType( 0 ), partitionSet, order_ );
146
147 // communicate_ indicates global norm
148 if( communicate_ )
149 {
150 sum = comm().sum( sum[ 0 ] );
151 }
152
153 // return result, e.g. sqrt of calculated sum
154 return sqrt( sum[ 0 ] );
155 }
156
157 template< class GridPart >
158 template< class UDiscreteFunctionType, class VDiscreteFunctionType, class PartitionSet >
159 inline typename Dune::FieldTraits< typename UDiscreteFunctionType::RangeFieldType >::real_type
160 H1Norm< GridPart >::distance ( const UDiscreteFunctionType &u,
161 const VDiscreteFunctionType &v,
162 const PartitionSet& partitionSet ) const
163 {
164 typedef typename UDiscreteFunctionType::RangeFieldType RangeFieldType;
165 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
166 typedef FieldVector< RealType, 1 > ReturnType ;
167
168 ReturnType sum = BaseType :: forEach( u, v, ReturnType( 0 ), partitionSet, order_ );
169
170 // communicate_ indicates global norm
171 if( communicate_ )
172 {
173 sum = comm().sum( sum[ 0 ] );
174 }
175
176 // return result, e.g. sqrt of calculated sum
177 return sqrt( sum[ 0 ] );
178 }
179
180 template< class GridPart >
181 template< class ULocalFunctionType, class VLocalFunctionType, class ReturnType >
182 inline void
183 H1Norm< GridPart >::distanceLocal ( const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum ) const
184 {
185 typedef typename L2Norm< GridPart >::template FunctionDistance< ULocalFunctionType, VLocalFunctionType > LocalDistanceType;
186
187 IntegratorType integrator( order );
188
189 LocalDistanceType dist( uLocal, vLocal );
191
192 integrator.integrateAdd( entity, dist2, sum );
193 }
194
195
196 template< class GridPart >
197 template< class LocalFunctionType, class ReturnType >
198 inline void
199 H1Norm< GridPart >::normLocal ( const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum ) const
200 {
201 // evaluate norm locally
202 IntegratorType integrator( order );
203
205
206 integrator.integrateAdd( entity, uLocal2, sum );
207 }
208
209 } // namespace Fem
210
211} // namespace Dune
212
213#endif // #ifndef DUNE_FEM_H1NORM_HH
Definition: bindguard.hh:11
static double sqrt(const Double &v)
Definition: double.hh:886
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
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian type
Definition: common/function.hh:70
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
Definition: h1norm.hh:15
Dune::FieldTraits< typenameDiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u, const PartitionSet &partitionSet) const
|| u ||_H1 on given set of entities (partition set)
Definition: h1norm.hh:139
H1Norm(const GridPartType &gridPart, const unsigned int order=0, const bool communicate=true)
constructor
Definition: h1norm.hh:129
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: h1norm.hh:183
Dune::FieldTraits< typenameDiscreteFunctionType::RangeFieldType >::real_type norm(const DiscreteFunctionType &u) const
|| u ||_H1 on interior partition entities
Definition: h1norm.hh:53
BaseType::EntityType EntityType
Definition: h1norm.hh:26
Dune::FieldTraits< typenameUDiscreteFunctionType::RangeFieldType >::real_type distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const PartitionSet &partitionSet) const
|| u - v ||_H1 on given set of entities (partition set)
Definition: h1norm.hh:160
Integrator< QuadratureType > IntegratorType
Definition: h1norm.hh:30
Dune::FieldTraits< typenameUDiscreteFunctionType::RangeFieldType >::real_type distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v) const
|| u - v ||_H1 on interior partition entities
Definition: h1norm.hh:66
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: h1norm.hh:199
const bool communicate_
Definition: h1norm.hh:81
const unsigned int order_
Definition: h1norm.hh:80
GridPart GridPartType
Definition: h1norm.hh:20
ThisType operator=(const ThisType &)=delete
CachingQuadrature< GridPartType, 0 > QuadratureType
Definition: h1norm.hh:28
void evaluate(const Point &x, RangeType &ret) const
Definition: h1norm.hh:105
Dune::FieldTraits< typenameFunctionType::RangeFieldType >::real_type RealType
Definition: h1norm.hh:96
FunctionJacobianSquare(const FunctionType &function)
Definition: h1norm.hh:100
FunctionType::RangeFieldType RangeFieldType
Definition: h1norm.hh:95
FieldVector< RealType, 1 > RangeType
Definition: h1norm.hh:97
Function FunctionType
Definition: h1norm.hh:93
Definition: l2norm.hh:20
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