1 #ifndef DUNE_FEM_INTEGRAL_HH
2 #define DUNE_FEM_INTEGRAL_HH
6 #include <dune/grid/common/rangegenerators.hh>
13 #include <dune/common/bartonnackmanifcheck.hh>
24 template<
class Gr
idPart,
class NormImplementation >
39 typedef typename GridPartType::template Codim< 0 >::EntityType
EntityType;
41 template <
bool uDiscrete,
bool vDiscrete>
44 template <
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class ReturnType,
class PartitionSet>
45 static ReturnType
forEach (
const ThisType &norm,
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
46 const ReturnType &initialValue,
47 const PartitionSet& partitionSet,
50 static_assert( uDiscrete && vDiscrete,
"Distance can only be calculated between GridFunctions" );
55 for(
const EntityType &entity : elements( norm.gridPart_, partitionSet ) )
57 uLocal.bind( entity );
58 vLocal.bind( entity );
59 const unsigned int uOrder = uLocal.order();
60 const unsigned int vOrder = vLocal.order();
61 const unsigned int orderLocal = (order == 0 ? 2*
std::max( uOrder, vOrder ) : order);
71 template <
bool vDiscrete>
75 class VDiscreteFunctionType,
80 const F& f,
const VDiscreteFunctionType& v,
81 const ReturnType& initialValue,
82 const PartitionSet& partitionSet,
83 const unsigned int order )
86 GridFunction u(
"Integral::adapter" , f , v.space().gridPart(), v.space().order() );
89 forEach( norm, u, v, initialValue, partitionSet, order );
94 template <
bool uDiscrete>
97 template <
class UDiscreteFunctionType,
103 const UDiscreteFunctionType& u,
105 const ReturnType& initialValue,
106 const PartitionSet& partitionSet,
107 const unsigned int order )
110 forEach( norm, f, u, initialValue, partitionSet, order );
114 template<
class DiscreteFunctionType,
class ReturnType,
class PartitionSet >
115 ReturnType
forEach (
const DiscreteFunctionType &u,
const ReturnType &initialValue,
116 const PartitionSet& partitionSet,
117 unsigned int order = 0 )
const
119 static_assert( (std::is_base_of<Fem::HasLocalFunction, DiscreteFunctionType>::value),
120 "Norm only implemented for quantities implementing a local function!" );
124 for(
const EntityType &entity : elements( gridPart_, partitionSet ) )
126 uLocal.bind( entity );
127 const unsigned int orderLocal = (order == 0 ? 2*uLocal.order() : order);
134 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class ReturnType,
class PartitionSet >
135 ReturnType
forEach (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
136 const ReturnType &initialValue,
const PartitionSet& partitionSet,
137 unsigned int order = 0 )
const
139 enum { uDiscrete = std::is_convertible<UDiscreteFunctionType, HasLocalFunction>::value };
140 enum { vDiscrete = std::is_convertible<VDiscreteFunctionType, HasLocalFunction>::value };
145 forEach( *
this, u, v, initialValue, partitionSet, order );
154 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
155 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &
sum )
const
157 CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
asImp().
distanceLocal( entity, order, uLocal, vLocal,
sum ) );
160 template<
class LocalFunctionType,
class ReturnType >
163 CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
asImp().
normLocal( entity, order, uLocal,
sum ) );
168 const typename GridPartType::CollectiveCommunicationType&
comm ()
const
176 bool commMax = communicate;
177 assert( communicate ==
comm().
max( commMax ) );
189 template<
class Gr
idPart >
203 template<
class UFunction,
class VFunction >
218 const unsigned int order = 0,
219 const bool communicate =
true );
223 template<
class DiscreteFunctionType,
class PartitionSet >
224 typename DiscreteFunctionType::RangeType
225 norm (
const DiscreteFunctionType &u,
const PartitionSet& partitionSet )
const;
228 template<
class DiscreteFunctionType >
229 typename DiscreteFunctionType::RangeType
230 norm (
const DiscreteFunctionType &u )
const
232 return norm( u, Partitions::interior );
236 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class PartitionSet >
237 typename UDiscreteFunctionType::RangeType
238 distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
const PartitionSet& partitionSet )
const;
241 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType >
242 typename UDiscreteFunctionType::RangeType
243 distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v )
const
245 return distance( u, v, Partitions::interior );
248 template<
class LocalFunctionType,
class ReturnType >
249 void normLocal (
const EntityType &entity,
unsigned int order,
const LocalFunctionType &uLocal, ReturnType &
sum )
const;
251 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
252 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &
sum )
const;
260 template<
class Gr
idPart >
264 communicate_(
BaseType::checkCommunicateFlag( communicate ) )
268 template<
class Gr
idPart >
269 template<
class DiscreteFunctionType,
class PartitionSet >
270 inline typename DiscreteFunctionType::RangeType
273 typedef typename DiscreteFunctionType::RangeType RangeType;
274 typedef RangeType ReturnType ;
289 template<
class Gr
idPart >
290 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class PartitionSet >
291 inline typename UDiscreteFunctionType::RangeType
293 ::distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
const PartitionSet& partitionSet )
const
295 typedef typename UDiscreteFunctionType::RangeType RangeType;
296 typedef RangeType ReturnType ;
310 template<
class Gr
idPart >
311 template<
class LocalFunctionType,
class ReturnType >
320 template<
class Gr
idPart >
321 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
329 LocalDistanceType dist( uLocal, vLocal );
334 template<
class Gr
idPart >
335 template<
class UFunction,
class VFunction >
349 template<
class Po
int >
353 u_.evaluate( x, ret );
354 v_.evaluate( x, phi );
358 template<
class Po
int >
362 u_.jacobian( x, ret );
363 v_.jacobian( x, phi );
368 const UFunctionType &u_;
369 const VFunctionType &v_;
Definition: bindguard.hh:11
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:517
static double max(const Double &v, const double p)
Definition: double.hh:398
static void forEach(IndexRange< T, sz > range, F &&f)
Definition: hybrid.hh:129
static constexpr T max(T a)
Definition: utility.hh:77
static constexpr std::decay_t< T > sum(T a)
Definition: utility.hh:33
BasicGridFunctionAdapter provides local functions for a Function.
Definition: gridfunctionadapter.hh:79
Definition: bartonnackmaninterface.hh:17
const NormImplementation & asImp() const
Definition: bartonnackmaninterface.hh:37
Definition: domainintegral.hh:28
GridPart GridPartType
Definition: domainintegral.hh:34
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: domainintegral.hh:155
IntegralBase(const GridPartType &gridPart)
Definition: domainintegral.hh:149
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: domainintegral.hh:161
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: domainintegral.hh:39
bool checkCommunicateFlag(bool communicate) const
Definition: domainintegral.hh:173
ReturnType forEach(const DiscreteFunctionType &u, const ReturnType &initialValue, const PartitionSet &partitionSet, unsigned int order=0) const
Definition: domainintegral.hh:115
const GridPartType & gridPart() const
Definition: domainintegral.hh:166
const GridPartType::CollectiveCommunicationType & comm() const
Definition: domainintegral.hh:168
ReturnType forEach(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const ReturnType &initialValue, const PartitionSet &partitionSet, unsigned int order=0) const
Definition: domainintegral.hh:135
Definition: domainintegral.hh:43
static ReturnType forEach(const ThisType &norm, const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const ReturnType &initialValue, const PartitionSet &partitionSet, unsigned int order)
Definition: domainintegral.hh:45
static ReturnType forEach(const ThisType &norm, const F &f, const VDiscreteFunctionType &v, const ReturnType &initialValue, const PartitionSet &partitionSet, const unsigned int order)
Definition: domainintegral.hh:79
static ReturnType forEach(const ThisType &norm, const UDiscreteFunctionType &u, const F &f, const ReturnType &initialValue, const PartitionSet &partitionSet, const unsigned int order)
Definition: domainintegral.hh:102
Definition: domainintegral.hh:191
DiscreteFunctionType::RangeType norm(const DiscreteFunctionType &u) const
|| u ||_L1 on interior partition entities
Definition: domainintegral.hh:230
CachingQuadrature< GridPartType, 0 > QuadratureType
Definition: domainintegral.hh:207
GridPart GridPartType
Definition: domainintegral.hh:196
const bool communicate_
Definition: domainintegral.hh:210
UDiscreteFunctionType::RangeType distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v) const
|| u - v ||_L2 on interior partition entities
Definition: domainintegral.hh:243
BaseType::EntityType EntityType
Definition: domainintegral.hh:204
DiscreteFunctionType::RangeType norm(const DiscreteFunctionType &u, const PartitionSet &partitionSet) const
|| u ||_L1 on given set of entities (partition set)
Definition: domainintegral.hh:271
Integral(const GridPartType &gridPart, const unsigned int order=0, const bool communicate=true)
constructor
Definition: domainintegral.hh:261
UDiscreteFunctionType::RangeType distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const PartitionSet &partitionSet) const
|| u - v ||_L2 on given set of entities (partition set)
Definition: domainintegral.hh:293
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: domainintegral.hh:313
const unsigned int order_
Definition: domainintegral.hh:209
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: domainintegral.hh:323
Definition: domainintegral.hh:337
UFunctionType::RangeType RangeType
Definition: domainintegral.hh:342
void evaluate(const Point &x, RangeType &ret) const
Definition: domainintegral.hh:350
UFunction UFunctionType
Definition: domainintegral.hh:338
FunctionDistance(const UFunctionType &u, const VFunctionType &v)
Definition: domainintegral.hh:345
UFunctionType::RangeFieldType RangeFieldType
Definition: domainintegral.hh:341
VFunction VFunctionType
Definition: domainintegral.hh:339
void jacobian(const Point &x, JacobianRangeType &ret) const
Definition: domainintegral.hh:359
UFunctionType::JacobianRangeType JacobianRangeType
Definition: domainintegral.hh:343
BaseType::EntityType EntityType
Definition: lpnorm.hh:143
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