42void Foam::faMeshDecomposition::distributeFaces()
44 const word& polyMeshRegionName =
mesh().
name();
46 Info<<
"\nCalculating distribution of finiteArea faces" <<
endl;
50 for (label proci = 0; proci <
nProcs(); ++proci)
64 processorDb.timeName(),
82 ioAddr.rename(
"faceProcAddressing");
85 labelHashSet faceProcAddressingHash(2*fvFaceProcAddressing.size());
95 if (hasGlobalFaceZones_)
98 ioAddr.rename(
"owner");
101 for (
int i = 0; i < ownerSize; ++i)
103 faceProcAddressingHash.insert(fvFaceProcAddressing[i]);
108 faceProcAddressingHash.insert
110 static_cast<labelList&
>(fvFaceProcAddressing)
117 const label index = (
faceLabels()[facei] + 1);
119 if (faceProcAddressingHash.found(index))
121 faceToProc_[facei] = proci;
126 Info<<
"\nFinished decomposition in "
127 << decompositionTime.elapsedCpuTime()
137 const label nProcessors,
142 nProcs_(nProcessors),
144 hasGlobalFaceZones_(false),
146 procFaceLabels_(nProcs_),
147 procMeshEdgesMap_(nProcs_),
148 procNInternalEdges_(nProcs_,
Zero),
149 procPatchEdgeLabels_(nProcs_),
150 procPatchPointAddressing_(nProcs_),
151 procPatchEdgeAddressing_(nProcs_),
152 procEdgeAddressing_(nProcs_),
153 procFaceAddressing_(nProcs_),
154 procBoundaryAddressing_(nProcs_),
155 procPatchSize_(nProcs_),
156 procPatchStartIndex_(nProcs_),
157 procNeighbourProcessors_(nProcs_),
158 procProcessorPatchSize_(nProcs_),
159 procProcessorPatchStartIndex_(nProcs_),
160 globallySharedPoints_(),
161 cyclicParallel_(false)
175 if (params.
found(
"globalFaceZones"))
177 hasGlobalFaceZones_ =
true;
189 Info<<
"\nDistributing faces to processors" <<
endl;
195 forAll(faceToProc_, facei)
197 const label proci = faceToProc_[facei];
199 if (proci < 0 || proci >= nProcs_)
202 <<
"Invalid processor label " << proci
203 <<
" for face " << facei <<
nl
208 ++nLocalFaces[proci];
213 forAll(nLocalFaces, proci)
215 procFaceAddressing_[proci].resize(nLocalFaces[proci]);
216 nLocalFaces[proci] = 0;
220 forAll(faceToProc_, facei)
222 const label proci = faceToProc_[facei];
223 const label localFacei = nLocalFaces[proci];
224 ++nLocalFaces[proci];
226 procFaceAddressing_[proci][localFacei] = facei;
232 for (label procI = 0; procI < nProcs(); procI++)
238 time().caseName()/(
"processor" +
Foam::name(procI))
264 ioAddr.
rename(
"pointProcAddressing");
271 ioAddr.
rename(
"faceProcAddressing");
274 fvFaceProcAddressingHash.
resize(2*fvFaceProcAddressing.
size());
276 forAll(fvFaceProcAddressing, facei)
278 fvFaceProcAddressingHash.
insert
280 fvFaceProcAddressing[facei],
286 const labelList& curProcFaceAddressing = procFaceAddressing_[procI];
288 labelList& curFaceLabels = procFaceLabels_[procI];
292 forAll(curProcFaceAddressing, faceI)
294 curFaceLabels[faceI] =
295 fvFaceProcAddressingHash.
find
297 faceLabels()[curProcFaceAddressing[faceI]] + 1
313 const label
nIntEdges = patch.nInternalEdges();
315 for (label edgei = 0; edgei <
nIntEdges; ++edgei)
317 edgesHash.
insert(patch.edges()[edgei], edgesHash.
size());
323 const label size =
boundary()[patchi].labelList::size();
325 for (label edgei=0; edgei < size; ++edgei)
329 patch.edges()[
boundary()[patchi][edgei]],
340 labelList& curPatchPointAddressing = procPatchPointAddressing_[procI];
341 curPatchPointAddressing.
resize(procMeshPoints.
size(), -1);
343 forAll(procMeshPoints, pointi)
345 curPatchPointAddressing[pointi] =
346 map[fvPointProcAddressing[procMeshPoints[pointi]]];
349 labelList& curPatchEdgeAddressing = procPatchEdgeAddressing_[procI];
350 curPatchEdgeAddressing.
resize(procEdges.
size(), -1);
352 Map<label>& curMap = procMeshEdgesMap_[procI];
358 edge curGlobalEdge(curPatchPointAddressing, procEdges[edgeI]);
359 curPatchEdgeAddressing[edgeI] = edgesHash.
find(curGlobalEdge).val();
362 forAll(curPatchEdgeAddressing, edgeI)
364 curMap.
insert(curPatchEdgeAddressing[edgeI], edgeI);
371 Info <<
"\nDistributing edges to processors" <<
endl;
378 const edgeList& edges = this->edges();
380 const labelList& neighbour = edgeNeighbour();
386 forAll(procEdgeList, procI)
388 for(label i=0; i<procNInternalEdges_[procI]; i++)
390 procEdgeList[procI].
append(procPatchEdgeAddressing_[procI][i]);
407 if (faceToProc_[owner[edgeI]] != faceToProc_[neighbour[edgeI]])
413 bool interProcBouFound =
false;
415 const label ownProc = faceToProc_[owner[edgeI]];
416 const label neiProc = faceToProc_[neighbour[edgeI]];
418 auto curInterProcBdrsOwnIter =
419 interProcBoundaries[ownProc].
cbegin();
421 auto curInterProcBEdgesOwnIter =
422 interProcBEdges[ownProc].
begin();
429 curInterProcBdrsOwnIter.good()
430 && curInterProcBEdgesOwnIter.good();
431 ++curInterProcBdrsOwnIter,
432 ++curInterProcBEdgesOwnIter
435 if (curInterProcBdrsOwnIter() == neiProc)
438 interProcBouFound =
true;
440 bool neighbourFound =
false;
442 curInterProcBEdgesOwnIter().append(edgeI);
444 auto curInterProcBdrsNeiIter =
445 interProcBoundaries[neiProc].
cbegin();
447 auto curInterProcBEdgesNeiIter =
448 interProcBEdges[neiProc].
begin();
455 curInterProcBdrsNeiIter.good()
456 && curInterProcBEdgesNeiIter.good();
457 ++curInterProcBdrsNeiIter,
458 ++curInterProcBEdgesNeiIter
461 if (curInterProcBdrsNeiIter() == ownProc)
464 neighbourFound =
true;
466 curInterProcBEdgesNeiIter().append(edgeI);
469 if (neighbourFound)
break;
472 if (interProcBouFound && !neighbourFound)
475 (
"faDomainDecomposition::decomposeMesh()")
476 <<
"Inconsistency in inter - "
477 <<
"processor boundary lists for processors "
478 << ownProc <<
" and " << neiProc
483 if (interProcBouFound)
break;
486 if (!interProcBouFound)
494 interProcBoundaries[ownProc].
append(neiProc);
498 interProcBoundaries[neiProc].
append(ownProc);
499 interProcBEdges[neiProc]
515 forAll(procPatchSize_, procI)
518 procPatchStartIndex_[procI].setSize(
patches.
size());
524 forAll(procPatchSize_, procI)
526 procPatchSize_[procI][patchI] = 0;
527 procPatchStartIndex_[procI][patchI] =
528 procEdgeList[procI].
size();
541 const label size =
patches[patchI].labelList::size();
545 for(
int eI=0; eI<size; eI++)
547 patchEdgeFaces[eI] = eF[
patches[patchI][eI]][0];
550 forAll(patchEdgeFaces, edgeI)
552 const label curProc = faceToProc_[patchEdgeFaces[edgeI]];
555 procEdgeList[curProc].
append(patchStart + edgeI);
558 procPatchSize_[curProc][patchI]++;
567 const label cycOffset = cPatch.
size()/2;
583 forAll(firstEdgeFaces, edgeI)
587 faceToProc_[firstEdgeFaces[edgeI]]
588 != faceToProc_[secondEdgeFaces[edgeI]]
597 cyclicParallel_ =
true;
599 bool interProcBouFound =
false;
601 const label ownProc =
602 faceToProc_[firstEdgeFaces[edgeI]];
603 const label neiProc =
604 faceToProc_[secondEdgeFaces[edgeI]];
606 auto curInterProcBdrsOwnIter =
607 interProcBoundaries[ownProc].
cbegin();
609 auto curInterProcBEdgesOwnIter =
610 interProcBEdges[ownProc].
begin();
617 curInterProcBdrsOwnIter.good()
618 && curInterProcBEdgesOwnIter.good();
619 ++curInterProcBdrsOwnIter,
620 ++curInterProcBEdgesOwnIter
623 if (curInterProcBdrsOwnIter() == neiProc)
627 interProcBouFound =
true;
629 bool neighbourFound =
false;
631 curInterProcBEdgesOwnIter()
632 .append(patchStart + edgeI);
634 auto curInterProcBdrsNeiIter
635 = interProcBoundaries[neiProc].
cbegin();
637 auto curInterProcBEdgesNeiIter =
638 interProcBEdges[neiProc].
begin();
645 curInterProcBdrsNeiIter.good()
646 && curInterProcBEdgesNeiIter.good();
647 ++curInterProcBdrsNeiIter,
648 ++curInterProcBEdgesNeiIter
651 if (curInterProcBdrsNeiIter() == ownProc)
654 neighbourFound =
true;
656 curInterProcBEdgesNeiIter()
665 if (neighbourFound)
break;
668 if (interProcBouFound && !neighbourFound)
672 "faDomainDecomposition::decomposeMesh()"
673 ) <<
"Inconsistency in inter-processor "
674 <<
"boundary lists for processors "
675 << ownProc <<
" and " << neiProc
676 <<
" in cyclic boundary matching"
681 if (interProcBouFound)
break;
684 if (!interProcBouFound)
692 interProcBoundaries[ownProc].
append(neiProc);
693 interProcBEdges[ownProc]
697 interProcBoundaries[neiProc].
append(ownProc);
698 interProcBEdges[neiProc]
713 label ownProc = faceToProc_[firstEdgeFaces[edgeI]];
716 procEdgeList[ownProc].
append(patchStart + edgeI);
719 procPatchSize_[ownProc][patchI]++;
730 forAll(secondEdgeFaces, edgeI)
734 faceToProc_[firstEdgeFaces[edgeI]]
735 == faceToProc_[secondEdgeFaces[edgeI]]
739 label ownProc = faceToProc_[firstEdgeFaces[edgeI]];
742 procEdgeList[ownProc].
append
743 (patchStart + cycOffset + edgeI);
746 procPatchSize_[ownProc][patchI]++;
754 forAll(procEdgeList, procI)
760 labelList& curProcEdgeAddressing = procEdgeAddressing_[procI];
763 procNeighbourProcessors_[procI];
766 procProcessorPatchSize_[procI];
768 labelList& curProcProcessorPatchStartIndex =
769 procProcessorPatchStartIndex_[procI];
772 label nEdgesOnProcessor = curProcEdges.size();
774 for (
const auto& bedges : interProcBEdges[procI])
776 nEdgesOnProcessor += bedges.size();
779 curProcEdgeAddressing.
setSize(nEdgesOnProcessor);
792 for (
const label procEdgei : curProcEdges)
794 curProcEdgeAddressing[
nEdges] = procEdgei;
802 curProcNeighbourProcessors.
setSize
804 interProcBoundaries[procI].size()
807 curProcProcessorPatchSize.
setSize
809 interProcBoundaries[procI].size()
812 curProcProcessorPatchStartIndex.
setSize
814 interProcBoundaries[procI].size()
819 auto curInterProcBdrsIter =
820 interProcBoundaries[procI].
cbegin();
822 auto curInterProcBEdgesIter =
823 interProcBEdges[procI].
cbegin();
828 curInterProcBdrsIter.good()
829 && curInterProcBEdgesIter.good();
830 ++curInterProcBdrsIter,
831 ++curInterProcBEdgesIter
835 curInterProcBdrsIter();
847 for (
const label edgei : *curInterProcBEdgesIter)
853 if (faceToProc_[owner[edgei]] == procI)
855 curProcEdgeAddressing[
nEdges] = edgei;
861 curProcEdgeAddressing[
nEdges] = edgei;
876 Info <<
"\nCalculating processor boundary addressing" <<
endl;
882 forAll(procPatchSize_, procI)
885 const labelList oldPatchSizes = procPatchSize_[procI];
887 const labelList oldPatchStarts = procPatchStartIndex_[procI];
889 labelList& curPatchSizes = procPatchSize_[procI];
891 labelList& curPatchStarts = procPatchStartIndex_[procI];
893 const labelList& curProcessorPatchSizes =
894 procProcessorPatchSize_[procI];
896 labelList& curBoundaryAddressing = procBoundaryAddressing_[procI];
901 + curProcessorPatchSizes.
size()
906 forAll(oldPatchSizes, patchI)
912 curBoundaryAddressing[
nPatches] = patchI;
914 curPatchSizes[
nPatches] = oldPatchSizes[patchI];
916 curPatchStarts[
nPatches] = oldPatchStarts[patchI];
926 forAll(curProcessorPatchSizes, procPatchI)
928 curBoundaryAddressing[
nPatches] = -1;
956 for (label procI = 0; procI < nProcs(); procI++)
959 const labelList& curEdgeLabels = procEdgeAddressing_[procI];
962 const labelList& curProcessorPatchStarts =
963 procProcessorPatchStartIndex_[procI];
965 const labelList& curProcessorPatchSizes =
966 procProcessorPatchSize_[procI];
971 forAll(curProcessorPatchStarts, patchI)
973 const label curStart = curProcessorPatchStarts[patchI];
974 const label curEnd = curStart + curProcessorPatchSizes[patchI];
978 label edgeI = curStart;
985 const label curE = curEdgeLabels[edgeI];
987 const edge&
e = edges[curE];
991 if (pointsUsage[
e[pointI]] == 0)
994 pointsUsage[
e[pointI]] = patchI + 1;
996 else if (pointsUsage[
e[pointI]] != patchI + 1)
999 gSharedPoints.
insert(
e[pointI]);
1007 globallySharedPoints_ = gSharedPoints.
sortedToc();
1013 for (label procI = 0; procI < nProcs(); procI++)
1017 time().caseName()/(
"processor" +
Foam::name(procI))
1048 const labelList& curEdgeAddressing = procEdgeAddressing_[procI];
1050 const labelList& curPatchStartIndex = procPatchStartIndex_[procI];
1051 const labelList& curPatchSize = procPatchSize_[procI];
1053 const labelList& curProcessorPatchStartIndex =
1054 procProcessorPatchStartIndex_[procI];
1056 const labelList& curProcessorPatchSize =
1057 procProcessorPatchSize_[procI];
1059 labelListList& curPatchEdgeLabels = procPatchEdgeLabels_[procI];
1060 curPatchEdgeLabels.
resize
1063 + curProcessorPatchSize.
size()
1066 forAll(curPatchSize, patchI)
1068 labelList& curEdgeLabels = curPatchEdgeLabels[patchI];
1069 curEdgeLabels.
setSize(curPatchSize[patchI], -1);
1075 int i=curPatchStartIndex[patchI];
1076 i<(curPatchStartIndex[patchI]+curPatchSize[patchI]);
1080 curEdgeLabels[edgeI] =
1081 procMeshEdgesMap_[procI][curEdgeAddressing[i]];
1086 forAll(curProcessorPatchSize, patchI)
1089 curPatchEdgeLabels[curPatchSize.
size() + patchI];
1090 curEdgeLabels.
setSize(curProcessorPatchSize[patchI], -1);
1096 int i=curProcessorPatchStartIndex[patchI];
1097 i<(curProcessorPatchStartIndex[patchI]
1098 +curProcessorPatchSize[patchI]);
1102 curEdgeLabels[edgeI] =
1103 procMeshEdgesMap_[procI][curEdgeAddressing[i]];
1115 Info<<
"\nConstructing processor FA meshes" <<
endl;
1118 Map<label> sharedPointLookup(2*globallySharedPoints_.size());
1120 forAll(globallySharedPoints_, pointi)
1122 sharedPointLookup.
insert(globallySharedPoints_[pointi], pointi);
1125 label maxProcFaces = 0;
1126 label totProcFaces = 0;
1127 label maxProcEdges = 0;
1128 label totProcEdges = 0;
1129 label maxProcPatches = 0;
1130 label totProcPatches = 0;
1133 for (label procI = 0; procI < nProcs(); procI++)
1139 time().caseName()/(
"processor" +
Foam::name(procI))
1165 "boundaryProcAddressing",
1183 const labelList& curBoundaryAddressing =
1184 procBoundaryAddressing_[procI];
1186 const labelList& curPatchSizes = procPatchSize_[procI];
1188 const labelList& curNeighbourProcessors =
1189 procNeighbourProcessors_[procI];
1191 const labelList& curProcessorPatchSizes =
1192 procProcessorPatchSize_[procI];
1195 procPatchEdgeLabels_[procI];
1201 curPatchSizes.
size() + curProcessorPatchSizes.
size()
1206 forAll(curPatchSizes, patchi)
1210 const label neiPolyPatchId =
1211 fvBoundaryProcAddressing.
find
1213 meshPatches[curBoundaryAddressing[patchi]]
1214 .ngbPolyPatchIndex()
1220 meshPatches[curBoundaryAddressing[patchi]].clone
1231 forAll(curProcessorPatchSizes, procPatchI)
1245 curNeighbourProcessors[procPatchI]
1261 Info<<
nl <<
"Processor " << procI;
1272 Info<<
"Number of faces = " << procMesh.
nFaces() <<
nl;
1279 totProcFaces += procMesh.
nFaces();
1280 maxProcFaces =
max(maxProcFaces, procMesh.
nFaces());
1282 label nBoundaryEdges = 0;
1288 const auto* ppp = isA<processorFaPatch>(fap);
1292 const auto& procPatch = *ppp;
1294 Info<<
" Number of edges shared with processor "
1295 << procPatch.neighbProcNo() <<
" = "
1296 << procPatch.size() <<
nl;
1303 nBoundaryEdges += fap.size();
1311 <<
" Number of boundary edges = " << nBoundaryEdges <<
nl;
1333 ioAddr.
rename(
"pointProcAddressing");
1337 ioAddr.
rename(
"edgeProcAddressing");
1341 ioAddr.
rename(
"faceProcAddressing");
1345 ioAddr.
rename(
"boundaryProcAddressing");
1352 <<
"Number of processor edges = " << (totProcEdges/2) <<
nl
1353 <<
"Max number of faces = " << maxProcFaces;
1355 if (maxProcFaces != totProcFaces)
1357 scalar avgValue = scalar(totProcFaces)/nProcs_;
1359 Info<<
" (" << 100.0*(maxProcFaces-avgValue)/avgValue
1360 <<
"% above average " << avgValue <<
')';
1364 Info<<
"Max number of processor patches = " << maxProcPatches;
1367 scalar avgValue = scalar(totProcPatches)/nProcs_;
1369 Info<<
" (" << 100.0*(maxProcPatches-avgValue)/avgValue
1370 <<
"% above average " << avgValue <<
')';
1374 Info<<
"Max number of edges between processors = " << maxProcEdges;
1377 scalar avgValue = scalar(totProcEdges)/nProcs_;
1379 Info<<
" (" << 100.0*(maxProcEdges-avgValue)/avgValue
1380 <<
"% above average " << avgValue <<
')';
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
Non-intrusive singly-linked list.
Map from edge (expressed as its endpoints) to value. For easier forward declaration it is currently i...
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
List< Key > sortedToc() const
The table of contents (the keys) in sorted order.
void resize(const label sz)
Resize the hash table for efficiency.
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
label size() const noexcept
The number of elements in table.
iterator find(const Key &key)
Find and return an iterator set at the hashed entry.
void clear()
Clear all entries from table.
A IOList wrapper for writing external data.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const fileName & caseName() const
Return the Time::caseName()
const word & name() const noexcept
Return the object name.
IOobject(const IOobject &)=default
Copy construct.
virtual void rename(const word &newName)
Rename the object.
const fileName & rootPath() const
Return the Time::rootPath()
static unsigned int defaultPrecision() noexcept
Return the default precision.
Template class for non-intrusive linked lists.
void setSize(const label n)
Alias for resize()
void append(const T &val)
Append an element at the end of the list.
void resize(const label len)
Adjust allocated size of list.
A HashTable to objects of type <T> with a label key.
A list of faces which address into the list of points.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
label nInternalEdges() const
Number of internal edges.
const labelList & meshPoints() const
Return labelList of mesh points in patch.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
const T * set(const label i) const
A List obtained as a section of another List.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
static word timeName(const scalar t, const int precision=precision_)
static word controlDictName
The default control dictionary name (normally "controlDict")
iterator begin() noexcept
Return an iterator to begin traversing the UList.
const_iterator cbegin() const noexcept
Return const_iterator to begin traversing the constant UList.
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
void size(const label n)
Older name for setAddressableSize.
label size() const noexcept
The number of elements in the list.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Finite area boundary mesh.
Automatic faMesh decomposition class.
bool writeDecomposition()
Write decomposition.
void updateParameters(const dictionary ¶ms)
Update flags based on the decomposition model settings.
label nProcs() const noexcept
Number of processor in decomposition.
void decomposeMesh()
Decompose mesh.
Finite area mesh (used for 2-D non-Euclidian finite area method) defined using a patch of faces on a ...
const Time & time() const
Return reference to time.
const faBoundaryMesh & boundary() const noexcept
Return constant reference to boundary mesh.
virtual bool write(const bool valid=true) const
Write mesh.
const polyMesh & mesh() const
Return access to polyMesh.
const uindirectPrimitivePatch & patch() const
Return constant reference to primitive patch.
virtual const objectRegistry & thisDb() const
Return reference to the mesh database.
label nPoints() const noexcept
Number of local mesh points.
label nFaces() const noexcept
Number of patch faces.
void addFaPatches(faPatchList &plist, const bool validBoundary=true)
Add boundary patches. Constructor helper.
const labelList & faceLabels() const noexcept
Return the underlying polyMesh face labels.
static word meshSubDir
The mesh sub-directory name (usually "faMesh")
Finite area patch class. Used for 2-D non-Euclidian finite area method.
virtual label size() const
Patch size is the number of edge labels.
const labelUList & edgeFaces() const
Return edge-face addressing.
A class for handling file names.
const word & name() const
Return reference to name.
label start() const noexcept
The start label of boundary faces in the polyMesh face list.
Mesh consisting of general polyhedral cells.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
virtual bool write(const bool valid=true) const
Write using setting from DB.
A class for handling words, derived from Foam::string.
const polyBoundaryMesh & patches
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const labelList nFaces(UPstream::listGatherValues< label >(aMesh.nFaces()))
const labelList nIntEdges(UPstream::listGatherValues< label >(aMesh.nInternalEdges()))
const labelList nProcEdges(UPstream::listGatherValues< label >(nLocalProcEdges))
const labelList nEdges(UPstream::listGatherValues< label >(aMesh.nEdges()))
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
List< label > labelList
A List of labels.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManip< error > abort(error &err)
IOList< label > labelIOList
Label container classes.
static constexpr const zero Zero
Global zero (0)
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
cpuTimeCxx cpuTime
Selection of preferred clock mechanism for the elapsed cpu time.
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.