1#ifndef DUNE_FEM_INTEGRAL_HH
2#define DUNE_FEM_INTEGRAL_HH
7#include <dune/grid/common/rangegenerators.hh>
14#include <dune/common/bartonnackmanifcheck.hh>
29 template<
class Gr
idPart,
class NormImplementation >
44 typedef typename GridPartType::template Codim< 0 >::EntityType
EntityType;
46 template <
bool uDiscrete,
bool vDiscrete>
50 template <
class IteratorRange,
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class ReturnType>
51 static ReturnType forEachLocal (
const ThisType &norm,
const IteratorRange& iterators,
52 const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
53 const ReturnType &initialValue,
56 static_assert( uDiscrete && vDiscrete,
"Distance can only be calculated between GridFunctions" );
64 auto uGuard =
bindGuard( uLocal, entity );
65 auto vGuard =
bindGuard( vLocal, entity );
67 const unsigned int uOrder = uLocal.order();
68 const unsigned int vOrder = vLocal.order();
69 const unsigned int orderLocal = (order == 0 ? 2*
std::max( uOrder, vOrder ) : order);
70 norm.
distanceLocal( entity, orderLocal, uLocal, vLocal, sum );
77 template <
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class ReturnType,
class PartitionSet>
78 static ReturnType
forEach (
const ThisType &norm,
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
79 const ReturnType &initialValue,
80 const PartitionSet& partitionSet,
85 return forEachLocal( norm, elements( norm.gridPart_, partitionSet ), u, v, initialValue, order );
87 std::vector< ReturnType > sums( nThreads, ReturnType(0) );
90 iterators( norm.gridPart_ );
92 auto doIntegrate = [ &norm, &iterators, &sums, &u, &v, &initialValue, &order ] ()
94 sums[
MPIManager::thread() ] = forEachLocal( norm, iterators, u, v, initialValue, order );
104 return forEachLocal( norm, elements( norm.gridPart_, partitionSet ), u, v, initialValue, order );
107 return std::accumulate( sums.begin(), sums.end(), ReturnType(0) );
112 template <
bool vDiscrete>
116 class VDiscreteFunctionType,
121 const F& f,
const VDiscreteFunctionType& v,
122 const ReturnType& initialValue,
123 const PartitionSet& partitionSet,
124 const unsigned int order )
127 GridFunction u(
"Integral::adapter" , f , v.space().gridPart(), v.space().order() );
130 forEach( norm, u, v, initialValue, partitionSet, order );
135 template <
bool uDiscrete>
138 template <
class UDiscreteFunctionType,
144 const UDiscreteFunctionType& u,
146 const ReturnType& initialValue,
147 const PartitionSet& partitionSet,
148 const unsigned int order )
151 forEach( norm, f, u, initialValue, partitionSet, order );
155 template <
class IteratorRange,
class UDiscreteFunctionType,
class ReturnType>
157 const UDiscreteFunctionType &u,
158 const ReturnType &initialValue,
159 unsigned int order )
const
166 auto uGuard =
bindGuard( uLocal, entity );
168 const unsigned int uOrder = uLocal.order();
169 const unsigned int orderLocal = (order == 0 ? 2*uOrder : order);
170 normLocal( entity, orderLocal, uLocal, sum );
176 template<
class DiscreteFunctionType,
class ReturnType,
class PartitionSet >
177 ReturnType
forEach (
const DiscreteFunctionType &u,
const ReturnType &initialValue,
178 const PartitionSet& partitionSet,
179 unsigned int order = 0 )
const
181 static_assert( (std::is_base_of<Fem::HasLocalFunction, DiscreteFunctionType>::value),
182 "Norm only implemented for quantities implementing a local function!" );
185 return forEachLocal( elements( gridPart_, partitionSet ), u, initialValue, order );
187 std::vector< ReturnType > sums( nThreads, ReturnType(0) );
190 iterators( gridPart_ );
192 auto doIntegrate = [
this, &iterators, &sums, &u, &initialValue, &order ] ()
204 return forEachLocal( elements( gridPart_, partitionSet ), u, initialValue, order );
207 return std::accumulate( sums.begin(), sums.end(), ReturnType(0) );
210 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class ReturnType,
class PartitionSet >
211 ReturnType
forEach (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
212 const ReturnType &initialValue,
const PartitionSet& partitionSet,
213 unsigned int order = 0 )
const
215 enum { uDiscrete = std::is_convertible<UDiscreteFunctionType, HasLocalFunction>::value };
216 enum { vDiscrete = std::is_convertible<VDiscreteFunctionType, HasLocalFunction>::value };
221 forEach( *
this, u, v, initialValue, partitionSet, order );
230 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
231 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &sum )
const
233 CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
asImp().
distanceLocal( entity, order, uLocal, vLocal, sum ) );
236 template<
class LocalFunctionType,
class ReturnType >
237 void normLocal (
const EntityType &entity,
unsigned int order,
const LocalFunctionType &uLocal, ReturnType &sum )
const
239 CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
asImp().
normLocal( entity, order, uLocal, sum ) );
244 const typename GridPartType::CollectiveCommunicationType&
comm ()
const
252 bool commMax = communicate;
253 assert( communicate ==
comm().
max( commMax ) );
266 template<
class Gr
idPart >
280 template<
class UFunction,
class VFunction >
295 const unsigned int order = 0,
296 const bool communicate =
true );
300 template<
class DiscreteFunctionType,
class PartitionSet >
301 typename DiscreteFunctionType::RangeType
302 norm (
const DiscreteFunctionType &u,
const PartitionSet& partitionSet )
const;
305 template<
class DiscreteFunctionType >
306 typename DiscreteFunctionType::RangeType
307 norm (
const DiscreteFunctionType &u )
const
309 return norm( u, Partitions::interior );
313 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class PartitionSet >
314 typename UDiscreteFunctionType::RangeType
315 distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
const PartitionSet& partitionSet )
const;
318 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType >
319 typename UDiscreteFunctionType::RangeType
320 distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v )
const
322 return distance( u, v, Partitions::interior );
325 template<
class LocalFunctionType,
class ReturnType >
326 void normLocal (
const EntityType &entity,
unsigned int order,
const LocalFunctionType &uLocal, ReturnType &sum )
const;
328 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
329 void distanceLocal (
const EntityType &entity,
unsigned int order,
const ULocalFunctionType &uLocal,
const VLocalFunctionType &vLocal, ReturnType &sum )
const;
337 template<
class Gr
idPart >
341 communicate_(
BaseType::checkCommunicateFlag( communicate ) )
345 template<
class Gr
idPart >
346 template<
class DiscreteFunctionType,
class PartitionSet >
347 inline typename DiscreteFunctionType::RangeType
350 typedef typename DiscreteFunctionType::RangeType RangeType;
351 typedef RangeType ReturnType ;
354 ReturnType sum = BaseType :: forEach( u, ReturnType(0), partitionSet, order_ );
359 sum = comm().sum( sum );
366 template<
class Gr
idPart >
367 template<
class UDiscreteFunctionType,
class VDiscreteFunctionType,
class PartitionSet >
368 inline typename UDiscreteFunctionType::RangeType
370 ::distance (
const UDiscreteFunctionType &u,
const VDiscreteFunctionType &v,
const PartitionSet& partitionSet )
const
372 typedef typename UDiscreteFunctionType::RangeType RangeType;
373 typedef RangeType ReturnType ;
376 ReturnType sum = BaseType :: forEach( u, v, ReturnType(0), partitionSet, order_ );
381 sum = comm().sum( sum );
387 template<
class Gr
idPart >
388 template<
class LocalFunctionType,
class ReturnType >
397 template<
class Gr
idPart >
398 template<
class ULocalFunctionType,
class VLocalFunctionType,
class ReturnType >
406 LocalDistanceType dist( uLocal, vLocal );
411 template<
class Gr
idPart >
412 template<
class UFunction,
class VFunction >
426 template<
class Po
int >
430 u_.evaluate( x, ret );
431 v_.evaluate( x, phi );
435 template<
class Po
int >
439 u_.jacobian( x, ret );
440 v_.jacobian( x, phi );
445 const UFunctionType &u_;
446 const VFunctionType &v_;
double max(const Dune::Fem::Double &v, const double p)
Definition: double.hh:965
Definition: bindguard.hh:11
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:604
static double max(const Double &v, const double p)
Definition: double.hh:398
static auto bindGuard(Object &object, Args &&... args) -> std::enable_if_t< isBindable< Object, Args... >::value, BindGuard< Object > >
Definition: bindguard.hh:67
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:33
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: domainintegral.hh:44
ReturnType forEachLocal(const IteratorRange &iterators, const UDiscreteFunctionType &u, const ReturnType &initialValue, unsigned int order) const
Definition: domainintegral.hh:156
GridPart GridPartType
Definition: domainintegral.hh:39
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: domainintegral.hh:231
IntegralBase(const GridPartType &gridPart)
Definition: domainintegral.hh:225
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: domainintegral.hh:237
const GridPartType::CollectiveCommunicationType & comm() const
Definition: domainintegral.hh:244
bool checkCommunicateFlag(bool communicate) const
Definition: domainintegral.hh:249
ReturnType forEach(const DiscreteFunctionType &u, const ReturnType &initialValue, const PartitionSet &partitionSet, unsigned int order=0) const
Definition: domainintegral.hh:177
ReturnType forEach(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const ReturnType &initialValue, const PartitionSet &partitionSet, unsigned int order=0) const
Definition: domainintegral.hh:211
const GridPartType & gridPart() const
Definition: domainintegral.hh:242
Definition: domainintegral.hh:48
static ReturnType forEach(const ThisType &norm, const UDiscreteFunctionType &u, const VDiscreteFunctionType &v, const ReturnType &initialValue, const PartitionSet &partitionSet, unsigned int order)
Definition: domainintegral.hh:78
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:120
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:143
Definition: domainintegral.hh:268
DiscreteFunctionType::RangeType norm(const DiscreteFunctionType &u) const
|| u ||_L1 on interior partition entities
Definition: domainintegral.hh:307
CachingQuadrature< GridPartType, 0 > QuadratureType
Definition: domainintegral.hh:284
GridPart GridPartType
Definition: domainintegral.hh:273
const bool communicate_
Definition: domainintegral.hh:287
UDiscreteFunctionType::RangeType distance(const UDiscreteFunctionType &u, const VDiscreteFunctionType &v) const
|| u - v ||_L2 on interior partition entities
Definition: domainintegral.hh:320
BaseType::EntityType EntityType
Definition: domainintegral.hh:283
DiscreteFunctionType::RangeType norm(const DiscreteFunctionType &u, const PartitionSet &partitionSet) const
|| u ||_L1 on given set of entities (partition set)
Definition: domainintegral.hh:348
Integral(const GridPartType &gridPart, const unsigned int order=0, const bool communicate=true)
constructor
Definition: domainintegral.hh:338
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:370
void normLocal(const EntityType &entity, unsigned int order, const LocalFunctionType &uLocal, ReturnType &sum) const
Definition: domainintegral.hh:390
const unsigned int order_
Definition: domainintegral.hh:286
void distanceLocal(const EntityType &entity, unsigned int order, const ULocalFunctionType &uLocal, const VLocalFunctionType &vLocal, ReturnType &sum) const
Definition: domainintegral.hh:400
Definition: domainintegral.hh:414
UFunctionType::RangeType RangeType
Definition: domainintegral.hh:419
void evaluate(const Point &x, RangeType &ret) const
Definition: domainintegral.hh:427
UFunction UFunctionType
Definition: domainintegral.hh:415
FunctionDistance(const UFunctionType &u, const VFunctionType &v)
Definition: domainintegral.hh:422
UFunctionType::RangeFieldType RangeFieldType
Definition: domainintegral.hh:418
VFunction VFunctionType
Definition: domainintegral.hh:416
void jacobian(const Point &x, JacobianRangeType &ret) const
Definition: domainintegral.hh:436
UFunctionType::JacobianRangeType JacobianRangeType
Definition: domainintegral.hh:420
BaseType::EntityType EntityType
Definition: lpnorm.hh:143
Exception thrown when a code segment that is supposed to be only accessed in single thread mode is ac...
Definition: mpimanager.hh:43
static void run(F &&f, Args &&... args)
run functor f with given arguments args in threaded mode
Definition: mpimanager.hh:440
static int thread()
return thread number
Definition: mpimanager.hh:424
static int numThreads()
return number of current threads
Definition: mpimanager.hh:421
Thread iterators.
Definition: threaditerator.hh:22
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