Peano
Loading...
Searching...
No Matches
Loop.h File Reference
#include <bitset>
#include "config.h"
#include "peano4/utils/Globals.h"
#include "tarch/compiler/CompilerSpecificSettings.h"
#include "tarch/accelerator/accelerator.h"
#include "tarch/la/Vector.h"
#include <oneapi/tbb/blocked_range.h>
#include <oneapi/tbb/blocked_range2d.h>
#include <oneapi/tbb/blocked_range3d.h>
#include <oneapi/tbb/parallel_for.h>
#include "tbb/Loop.h"
Include dependency graph for Loop.h:

Go to the source code of this file.

Namespaces

namespace  peano4
 
namespace  peano4::utils
 

Macros

#define CPUGPUMethod
 
#define dfor(counter, max)
 d-dimensional Loop
 
#define parallelFor(counter, max)    parallelForWithSchedulerInstructions(counter, max, ::peano4::utils::LoopPlacement::SpreadOut)
 
#define simtFor(counter, max)    simtForWithSchedulerInstructions(counter, max, ::peano4::utils::LoopPlacement::SpreadOut)
 
#define parallelDfor(counter, max)    parallelDforWithSchedulerInstructions(counter, max, ::peano4::utils::LoopPlacement::SpreadOut)
 
#define dfor4(counter)
 Shortcut For dfor(counter,4)
 
#define dfor3(counter)
 Shortcut For dfor(counter,3)
 
#define dfor5(counter)
 
#define dfor7(counter)
 
#define dfor9(counter)
 
#define d2for(counter, max)
 If Dimensions is not set to two, we might nevertheless need two-dimensional loops.
 
#define d2for2(counter)
 If Dimensions is not set to two, we might nevertheless need two-dimensional loops.
 
#define d3for2(counter)
 If Dimensions is not set to three, we might nevertheless need three-dimensional loops.
 
#define d3for(counter, max)
 If Dimensions is not set to three, we might nevertheless need two-dimensional loops.
 
#define d2for3(counter)
 If Dimensions is not set to two, we might nevertheless need two-dimensional loops.
 
#define d3for3(counter)
 If Dimensions is not set to three, we might nevertheless need three-dimensional loops.
 
#define dfor2(counter)
 Shortcut For dfor(counter,2)
 
#define enddforx   }
 I prefer to use this macro for dforx instead of a closing bracket as many syntax parser fail otherwise.
 
#define dfore(counter, max, dim, value)
 This is an exclusive d-dimensional for loop.
 
#define zfor(counter, max, direction)
 This is a d-dimensional z-loop.
 
#define zfor3(counter, direction)   zfor(counter, 3, direction)
 
#define endzfor
 

Typedefs

typedef std::bitset< Dimensions > peano4::utils::LoopDirection
 Is used by the z-loop.
 

Enumerations

enum class  peano4::utils::LoopPlacement { peano4::utils::Serial , peano4::utils::Nested , peano4::utils::SpreadOut }
 Guide loop-level parallelism. More...
 

Functions

CPUGPUMethod void peano4::utils::dInc (tarch::la::Vector< Dimensions, int > &counter, int max)
 d-dimensional counterpart of increment operator
 
void peano4::utils::dDec (tarch::la::Vector< Dimensions, int > &counter, int max)
 d-dimensional counterpart of decrement operator
 
void peano4::utils::dInc (tarch::la::Vector< Dimensions, int > &counter, const tarch::la::Vector< Dimensions, int > &max)
 d-dimensional counterpart of increment, where individual vector components have different max values
 
void peano4::utils::dIncByVector (tarch::la::Vector< Dimensions, int > &counter, int max, int increment)
 Perform a d-dimensional increment by value increment: The first component of the counter is incremented by increment.
 
void peano4::utils::dIncByScalar (tarch::la::Vector< Dimensions, int > &counter, int max, int increment)
 Perform a scalar increment of a vector: The operation equals a sequence of increment calls to dInc().
 
void peano4::utils::dInc (tarch::la::Vector< Dimensions, int > &counter, int max, int doNotExamine)
 Same operation as dInc(tarch::la::Vector<Dimensions,int>,int), but now one dimension is not taken into consideration.
 
void peano4::utils::dInc (tarch::la::Vector< Dimensions, int > &counter, int max, LoopDirection &direction)
 Operation similar to dInc, but is given a direction bitset that identifies whether the counters has to be incremented or decremented.
 
CPUGPUMethod int peano4::utils::dCmp (const tarch::la::Vector< Dimensions, int > &counter, int max)
 Element-wise comparison for the for d-dimensional for loops.
 
int peano4::utils::dCmp (const tarch::la::Vector< Dimensions, int > &counter, const tarch::la::Vector< Dimensions, int > &max)
 Element-wise comparison for the for d-dimensional for loops.
 
bool peano4::utils::dCmpLinearOrder (const tarch::la::Vector< Dimensions, int > &counter, const tarch::la::Vector< Dimensions, int > &max)
 Compares two vectors with regards to their linearised value.
 
CPUGPUMethod int peano4::utils::dLinearised (const tarch::la::Vector< Dimensions, int > &counter, int max)
 Map d-dimensional vector onto integer index.
 
CPUGPUMethod int peano4::utils::dLinearised (const std::bitset< Dimensions > &counter)
 
int peano4::utils::d2Linearised (const tarch::la::Vector< 2, int > &counter, int max)
 Special 2d variant of dLinearised that works also if you compile with other dimensions.
 
int peano4::utils::d3Linearised (const tarch::la::Vector< 3, int > &counter, int max)
 Special 3d variant of dLinearised that works also if you compile with other dimensions.
 
int peano4::utils::dLinearisedWithoutLookup (const tarch::la::Vector< Dimensions, int > &counter, int max)
 Linearisation not Optimised.
 
tarch::la::Vector< Dimensions, intpeano4::utils::dDelinearised (int value, int max)
 Counterpart of dLinearised().
 
tarch::la::Vector< Dimensions, intpeano4::utils::dDelinearisedWithoutLookup (int value, int max)
 Delinearization not optimised.
 
void peano4::utils::setupLookupTableForDLinearised ()
 
void peano4::utils::setupLookupTableForDDelinearised ()
 
CPUGPUMethod tarch::la::Vector< Dimensions, intpeano4::utils::dStartVector ()
 Construct start vector (0,0,....) for d-dimensional loop.
 
tarch::la::Vector< Dimensions, intpeano4::utils::dStartVector (int dim, int value)
 
tarch::la::Vector< Dimensions, intpeano4::utils::dStartVector (int max, const LoopDirection &direction)
 Creates a start vector.
 
std::string toString (peano4::utils::LoopPlacement)
 

Macro Definition Documentation

◆ CPUGPUMethod

#define CPUGPUMethod

Definition at line 35 of file Loop.h.

◆ d2for

#define d2for ( counter,
max )
Value:
for (int counter##Scalar = 0; counter##Scalar < tarch::la::aPowI(2, max); counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = 2 - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= max; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}
int aPowI(int i, int a)
Computes the i-th power of a in integer arithmetic.
Simple vector class.
Definition Vector.h:134

If Dimensions is not set to two, we might nevertheless need two-dimensional loops.

So this is the corresponding macro. It is way slower than dfor if you compile with -DDimensions=2.

Please use this macro with an enddforx macro closing your scope rather than brackets.

Please note that counterScalar is already a linearised version of your counter.

Please note that you need a specialised linearisation function (depending on d explicitly) to work with 2d index vectors within such a loop. Do not just use dLinearised, but use the d2Linearised or d3Linearised variant instead.

Definition at line 723 of file Loop.h.

Referenced by toolbox::finiteelements::BSplinesStencilFactory::getElementWiseAssemblyMatrix().

◆ d2for2

#define d2for2 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < tarch::la::aPowI(2, 2); counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = 2 - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 2; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

If Dimensions is not set to two, we might nevertheless need two-dimensional loops.

So this is the corresponding macro.

Please use enddforx to close a loop started with this macro.

Definition at line 744 of file Loop.h.

◆ d2for3

#define d2for3 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < tarch::la::aPowI(2, 3); counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = 2 - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 3; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

If Dimensions is not set to two, we might nevertheless need two-dimensional loops.

So this is the corresponding macro.

Please use enddforx to close a loop started with this macro.

Definition at line 817 of file Loop.h.

◆ d3for

#define d3for ( counter,
max )
Value:
for (int counter##Scalar = 0; counter##Scalar < tarch::la::aPowI(3, max); counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = 3 - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= max; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

If Dimensions is not set to three, we might nevertheless need two-dimensional loops.

So this is the corresponding macro. It is way slower than dfor if you compile with -DDimensions=2.

Please use this macro with an enddforx macro closing your scope rather than brackets.

Please note that counterScalar is already a linearised version of your counter.

Please note that you need a specialised linearisation function (depending on d explicitly) to work with 2d index vectors within such a loop. Do not just use dLinearised, but use the d2Linearised or d3Linearised variant instead.

Definition at line 795 of file Loop.h.

Referenced by toolbox::finiteelements::BSplinesStencilFactory::getElementWiseAssemblyMatrix().

◆ d3for2

#define d3for2 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < tarch::la::aPowI(3, 2); counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = 3 - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 2; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

If Dimensions is not set to three, we might nevertheless need three-dimensional loops.

So this is the corresponding macro.

Please use enddforx to close a loop started with this macro.

Definition at line 765 of file Loop.h.

◆ d3for3

#define d3for3 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < tarch::la::aPowI(3, 3); counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = 3 - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 3; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

If Dimensions is not set to three, we might nevertheless need three-dimensional loops.

So this is the corresponding macro.

Please use enddforx to close a loop started with this macro.

Definition at line 838 of file Loop.h.

Referenced by toolbox::finiteelements::exchangeCoordinates().

◆ dfor

#define dfor ( counter,
max )
Value:
peano4::utils::dInc(counter, max))
CPUGPUMethod tarch::la::Vector< Dimensions, int > dStartVector()
Construct start vector (0,0,....) for d-dimensional loop.
Definition Loop.cpp:316
CPUGPUMethod int dCmp(const tarch::la::Vector< Dimensions, int > &counter, int max)
Element-wise comparison for the for d-dimensional for loops.
Definition Loop.cpp:292

d-dimensional Loop

Very often one needs a d-dimensional for loop. A d-dimensional for loop is something like

for (x(0)=0; x(0)<N; x(0)++)
for (x(1)=0; x(1)<N; x(1)++)
for (x(2)=0; x(2)<N; x(2)++)

with d nested for loops. Thus, one has to code such loops for every d manually. This macro offers a d-independent alternative, just write

dfor (x,N) {
...
}
#define dfor(counter, max)
d-dimensional Loop
Definition Loop.h:313

The precompiler extracts this macro and within the loop body, you are able to use the integer tinyvector x.

Here is an example:

dfor(a,2) {
std::cout << a << ",";
}
const float const float const float struct part *restrict struct part *restrict const float a
Definition hydro_iact.h:55

results in [0,0], [1,0], [0,1], [1,1] if Dimensions equals 2. If DIMENSION equals 3 the same construct gives you [0,0,0], [1,0,0], [0,1,0], [1,1,0], [0,0,1], [1,0,1], [0,1,1], [1,1,1].

There are some short-cuts for frequently used iteration ranges, i.e. popular N. Please note that this routine works for quadratic iteration ranges, i.e. if you make N a scalar. It also works if max has the type tarch::la::Vector<Dimensions,int> and hosts different entries per dimension.

Definition at line 313 of file Loop.h.

Referenced by exahype2::dg::cellIntegral_patchwise_in_situ_GaussLegendre(), exahype2::dg::cellIntegral_patchwise_in_situ_GaussLegendre_functors(), toolbox::blockstructured::computeGradient(), toolbox::blockstructured::computeGradientAndReturnMaxDifference(), toolbox::particles::createEquallySpacedParticles(), toolbox::particles::createParticlesAlignedWithGlobalCartesianMesh(), toolbox::blockstructured::internal::createPiecewiseConstantInterpolationMatrix(), exahype2::dg::evaluatePolynomial(), toolbox::finiteelements::getPoissonMatrixWithJumpingCoefficient(), exahype2::dg::integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre(), exahype2::dg::integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre(), toolbox::blockstructured::interpolateCellDataAssociatedToVolumesIntoOverlappingCell_linear(), toolbox::blockstructured::interpolateCellDataAssociatedToVolumesIntoOverlappingCell_secondOrder(), toolbox::blockstructured::interpolateHaloLayer_AoS_matrix(), exahype2::dg::multiplyWithInvertedMassMatrix_GaussLegendre(), exahype2::fv::plotPatch(), toolbox::blockstructured::internal::projectCells_AoS(), exahype2::dg::reduceMaxEigenvalue_patchwise_functors(), exahype2::fd::reduceMaxEigenvalue_patchwise_functors(), toolbox::blockstructured::restrictCell_AoS_inject(), toolbox::blockstructured::restrictCellIntoOverlappingCell_inject(), toolbox::blockstructured::restrictCellIntoOverlappingCell_inject_and_average(), toolbox::blockstructured::restrictInnerHalfOfHaloLayer_AoS_matrix(), exahype2::dg::tests::CellIntegralTest::runEulerOrder2OnStationarySetup(), peano4::utils::setupLookupTableForDLinearised(), exahype2::fv::rusanov::tests::CopyPatchTest::testCopyPatch(), peano4::utils::tests::ParallelDForTest::testParallelDFor(), exahype2::fv::rusanov::internal::timeStepWithRusanovBatchedFunctors(), exahype2::fv::rusanov::internal::timeStepWithRusanovBatchedStateless(), exahype2::fv::rusanov::internal::timeStepWithRusanovPatchwiseFunctors(), exahype2::fv::rusanov::internal::timeStepWithRusanovPatchwiseStateless(), and exahype2::fv::validatePatch().

◆ dfor2

#define dfor2 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < TwoPowerD; counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = Dimensions - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 2; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}
#define TwoPowerD
Definition Globals.h:19

Shortcut For dfor(counter,2)

The usage of this optimised shortcut differs from dfor: You have to replace both the dfor and the opening bracket by this macro, i.e.

dfor(counter,2) {

becomes

dfor2(counter)
#define dfor2(counter)
Shortcut For dfor(counter,2)
Definition Loop.h:911

You usually use this macro with

#pragma unroll(TwoPowerD)

or

#pragma omp parallel for schedule(static)

If you work with this specialised version of dfor on a variable k, two counter variables are available within the loop's scope. The variable k itself with type tarch::la::Vector<Dimensions,int>. Furthermore, there's always a variable kScalar giving you k's value linearised.

Definition at line 911 of file Loop.h.

Referenced by peano4::datamanagement::VertexMarker::areAdjacentCellsLocal(), peano4::grid::Spacetree::areAllVerticesNonHanging(), peano4::grid::Spacetree::areAllVerticesRefined(), peano4::grid::Spacetree::areAllVerticesUnrefined(), peano4::grid::GridTraversalEventGenerator::createGenericCellTraversalEvent(), peano4::grid::Spacetree::descend(), peano4::grid::Spacetree::evaluateGridControlEvents(), toolbox::finiteelements::extractElementStencil(), toolbox::particles::getAdjacentCellsOwningParticle(), peano4::grid::Spacetree::getAdjacentRanksForNewVertex(), toolbox::finiteelements::getElementWiseAssemblyMatrix(), toolbox::finiteelements::getElementWiseAssemblyMatrix(), toolbox::finiteelements::getElementWiseAssemblyMatrix(), toolbox::finiteelements::getElementWiseAssemblyMatrix(), peano4::parallel::Node::getPeriodicBoundaryNumber(), toolbox::finiteelements::getPoissonMatrixWithJumpingCoefficient(), peano4::grid::GridTraversalEventGenerator::getTreeOwningSpacetreeNode(), peano4::grid::Spacetree::incrementNumberOfAdjacentRefinedLocalCells(), peano4::grid::Spacetree::isCellSplitCandidate(), peano4::grid::GridTraversalEventGenerator::isSpacetreeNodeLocal(), peano4::grid::isSpacetreeNodeRefined(), peano4::grid::Spacetree::markVerticesAroundForkedCell(), peano4::grid::TraversalVTKPlotter::plotCell(), exahype2::fv::internal::projectValueOntoParticle_piecewiseLinear(), peano4::grid::Spacetree::receiveAndMergeGridVertexAtHorizontalBoundary(), toolbox::finiteelements::reconstructStencilFragments(), toolbox::finiteelements::reconstructUniformStencilFragments(), peano4::grid::tests::GridTraversalEventGeneratorTest::testCreateLeaveCellTraversalEvent1(), peano4::grid::Spacetree::traverse(), peano4::grid::Spacetree::updateVertexBeforeStore(), and peano4::grid::Spacetree::updateVertexRanksWithinCell().

◆ dfor3

#define dfor3 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < ThreePowerD; counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = Dimensions - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 3; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}
#define ThreePowerD
Definition Globals.h:24

Shortcut For dfor(counter,3)

The usage of this optimised shortcut differs from dfor: You have to replace both the dfor and the opening bracket by this macro, i.e.

dfor(counter,3) {

becomes

dfor3(counter)
#define dfor3(counter)
Shortcut For dfor(counter,3)
Definition Loop.h:647

You usually use this macro with

#pragma unroll(ThreePowerD)

or

#pragma omp parallel for schedule(static)

If you work with this specialised version of dfor on a variable k, two counter variables are available within the loop's scope. The variable k itself with type tarch::la::Vector<Dimensions,int>. Furthermore, there's always a variable kScalar giving you k's value linearised.

Definition at line 647 of file Loop.h.

Referenced by toolbox::blockstructured::internal::createPiecewiseConstantInterpolationMatrix(), peano4::grid::Spacetree::descend(), and peano4::parallel::Node::getOutputStacksForPeriodicBoundaryExchange().

◆ dfor4

#define dfor4 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < FOUR_POWER_D; counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = Dimensions - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 4; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

Shortcut For dfor(counter,4)

The usage of this optimised shortcut differs from dfor: You have to replace both the dfor and the opening bracket by this macro, i.e.

dfor(counter,4) {

becomes

dfor4(counter)
#define dfor4(counter)
Shortcut For dfor(counter,4)
Definition Loop.h:601

You usually use this macro with

#pragma unroll(FOUR_POWER_D)

or

#pragma omp parallel for schedule(static)

If you work with this specialised version of dfor on a variable k, two counter variables are available within the loop's scope. The variable k itself with type tarch::la::Vector<Dimensions,int>. Furthermore, there's always a variable kScalar giving you k's value linearised.

Definition at line 601 of file Loop.h.

◆ dfor5

#define dfor5 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < FIVE_POWER_D; counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = Dimensions - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 5; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

Definition at line 663 of file Loop.h.

◆ dfor7

#define dfor7 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < SEVEN_POWER_D; counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = Dimensions - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 7; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

Definition at line 678 of file Loop.h.

◆ dfor9

#define dfor9 ( counter)
Value:
for (int counter##Scalar = 0; counter##Scalar < NINE_POWER_D; counter##Scalar++) { \
{ \
int copy##counter##Scalar = counter##Scalar; \
for (int counter##ddd = Dimensions - 1; counter##ddd >= 0; counter##ddd--) { \
int counter##aPowI = 1; \
for (int counter##jjj = 0; counter##jjj < counter##ddd; counter##jjj++) { \
counter##aPowI *= 9; \
} \
counter(counter##ddd) = copy##counter##Scalar / counter##aPowI; \
copy##counter##Scalar -= counter(counter##ddd) * counter##aPowI; \
} \
}

Definition at line 694 of file Loop.h.

◆ dfore

#define dfore ( counter,
max,
dim,
value )
Value:
peano4::utils::dCmp(counter, max); \
peano4::utils::dInc(counter, max, dim))

This is an exclusive d-dimensional for loop.

Exclusive means, there is one dimension that is not manipulated during the for loop. This dimension (entry of the counter) is specified by dim and has the value value throughout the for-loop.

Definition at line 942 of file Loop.h.

Referenced by exahype2::fd::applyBoundaryConditions(), exahype2::fv::applyBoundaryConditions(), exahype2::dg::applyBoundaryConditions(), exahype2::fd::applySommerfeldConditions(), toolbox::blockstructured::internal::clearHalfOfHaloLayerAoS(), toolbox::blockstructured::clearHaloLayerAoS(), exahype2::fv::copyHalfOfHalo(), toolbox::blockstructured::extrapolatePatchSolutionAndProjectExtrapolatedHaloOntoFaces(), peano4::grid::GridTraversalEventGenerator::getAdjacentRanksOfFace(), toolbox::blockstructured::interpolateHaloLayer_AoS_tensor_product(), exahype2::dg::interpolateRiemannSolution(), exahype2::fv::mapInnerNeighbourVoxelAlongBoundayOntoAuxiliaryVariable(), toolbox::finiteelements::preprocessBoundaryStencil(), toolbox::blockstructured::internal::projectInterpolatedFineCellsOnHaloLayer_AoS(), toolbox::blockstructured::projectPatchHaloOntoFaces(), toolbox::blockstructured::projectPatchSolutionOntoFaces(), exahype2::dg::restrictAndAccumulateProjectedFacePolynomial(), toolbox::blockstructured::restrictInnerHalfOfHaloLayer_AoS_averaging(), toolbox::blockstructured::restrictInnerHalfOfHaloLayer_AoS_inject(), toolbox::blockstructured::restrictInnerHalfOfHaloLayer_AoS_tensor_product(), exahype2::dg::fluxaverage::solveRiemannProblem_pointwise_in_situ(), exahype2::dg::rusanov::solveRiemannProblem_pointwise_in_situ(), exahype2::dg::laxfriedrichs::solveRiemannProblem_pointwise_in_situ(), and exahype2::dg::average::solveRiemannProblem_pointwise_in_situ().

◆ enddforx

#define enddforx   }

I prefer to use this macro for dforx instead of a closing bracket as many syntax parser fail otherwise.

Definition at line 933 of file Loop.h.

Referenced by peano4::datamanagement::VertexMarker::areAdjacentCellsLocal(), peano4::grid::Spacetree::areAllVerticesNonHanging(), peano4::grid::Spacetree::areAllVerticesRefined(), peano4::grid::Spacetree::areAllVerticesUnrefined(), peano4::grid::GridTraversalEventGenerator::createGenericCellTraversalEvent(), toolbox::blockstructured::internal::createPiecewiseConstantInterpolationMatrix(), peano4::grid::Spacetree::descend(), peano4::grid::Spacetree::evaluateGridControlEvents(), toolbox::finiteelements::exchangeCoordinates(), toolbox::finiteelements::extractElementStencil(), toolbox::particles::getAdjacentCellsOwningParticle(), peano4::grid::Spacetree::getAdjacentRanksForNewVertex(), toolbox::finiteelements::getElementWiseAssemblyMatrix(), toolbox::finiteelements::getElementWiseAssemblyMatrix(), toolbox::finiteelements::getElementWiseAssemblyMatrix(), toolbox::finiteelements::getElementWiseAssemblyMatrix(), toolbox::finiteelements::BSplinesStencilFactory::getElementWiseAssemblyMatrix(), toolbox::finiteelements::BSplinesStencilFactory::getElementWiseAssemblyMatrix(), peano4::parallel::Node::getOutputStacksForPeriodicBoundaryExchange(), peano4::parallel::Node::getPeriodicBoundaryNumber(), toolbox::finiteelements::getPoissonMatrixWithJumpingCoefficient(), peano4::grid::GridTraversalEventGenerator::getTreeOwningSpacetreeNode(), peano4::grid::Spacetree::incrementNumberOfAdjacentRefinedLocalCells(), peano4::grid::Spacetree::isCellSplitCandidate(), peano4::grid::GridTraversalEventGenerator::isSpacetreeNodeLocal(), peano4::grid::isSpacetreeNodeRefined(), peano4::grid::Spacetree::markVerticesAroundForkedCell(), peano4::grid::TraversalVTKPlotter::plotCell(), exahype2::fv::internal::projectValueOntoParticle_piecewiseLinear(), peano4::grid::Spacetree::receiveAndMergeGridVertexAtHorizontalBoundary(), toolbox::finiteelements::reconstructStencilFragments(), toolbox::finiteelements::reconstructUniformStencilFragments(), peano4::grid::Spacetree::traverse(), peano4::grid::Spacetree::updateVertexBeforeStore(), and peano4::grid::Spacetree::updateVertexRanksWithinCell().

◆ endzfor

#define endzfor
Value:
} \
}

Definition at line 1049 of file Loop.h.

Referenced by peano4::grid::Spacetree::descend().

◆ parallelDfor

#define parallelDfor ( counter,
max )    parallelDforWithSchedulerInstructions(counter, max, ::peano4::utils::LoopPlacement::SpreadOut)

Definition at line 452 of file Loop.h.

Referenced by peano4::utils::tests::ParallelDForTest::testParallelDFor().

◆ parallelFor

#define parallelFor ( counter,
max )    parallelForWithSchedulerInstructions(counter, max, ::peano4::utils::LoopPlacement::SpreadOut)

Definition at line 446 of file Loop.h.

Referenced by runBenchmarks().

◆ simtFor

#define simtFor ( counter,
max )    simtForWithSchedulerInstructions(counter, max, ::peano4::utils::LoopPlacement::SpreadOut)

Definition at line 449 of file Loop.h.

◆ zfor

#define zfor ( counter,
max,
direction )
Value:
{ \
peano4::utils::dCmp(counter, max); \
peano4::utils::dInc(counter, max, direction)) {

This is a d-dimensional z-loop.

A z-loop is a d-dimensional loop the counter direction changes every time an inner loop direction has changed. So this is the loop corresponding to a Peano curve. The for loop is passed a counter name, the number of steps to perform in each direction and a direction flag that identifies the initial direction. Note that this argument has to be a real variable, it might not be a constant. The direction flag array identifies for each direction, whether the initial loop goes along the axis or not. The type of direction is LoopDirection.

Here are some examples for two dimensions:

LoopDirection d(3); // equals {true,true} and identifies the standard
// Peano Leitmotiv
zfor( a, 3, d ) {
std::cout << a;
}
#define zfor(counter, max, direction)
This is a d-dimensional z-loop.
Definition Loop.h:979

yields in [0,0],[1,0],[2,0],[2,1],[1,1],[0,1],[0,2],[1,2],[2,2].

LoopDirection d(1); // equals {true, false} and specifies a Peano curve
// from the left top to right bottom
zfor( a, 3, d ) {
std::cout << a;
}

yields in [0,2],[1,2],[2,2],[2,1],[1,1],[0,1],[0,0],[1,0],[2,0].

Definition at line 979 of file Loop.h.

◆ zfor3

#define zfor3 ( counter,
direction )   zfor(counter, 3, direction)

Definition at line 1044 of file Loop.h.

Referenced by peano4::grid::Spacetree::descend().

Function Documentation

◆ toString()

std::string toString ( peano4::utils::LoopPlacement arg)