Go to the documentation of this file.
81 Info<<
"GAMGAgglomeration:" <<
nl
82 <<
" local agglomerator : " <<
type() <<
nl;
85 Info<<
" processor agglomerator : "
91 <<
setw(20) <<
"nFaces/nCells"
92 <<
setw(20) <<
"nInterfaces"
93 <<
setw(20) <<
"nIntFaces/nCells"
94 <<
setw(12) <<
"profile"
97 <<
setw(8) <<
"nProcs"
113 <<
setw(8) <<
"-----"
114 <<
setw(8) <<
"------"
132 for (label levelI = 0; levelI <=
size(); levelI++)
136 scalar faceCellRatio = 0;
137 label nInterfaces = 0;
140 scalar profile = 0.0;
155 if (interfaces.set(i))
158 nIntFaces += interfaces[i].faceCells().size();
161 ratio = scalar(nIntFaces)/
nCells;
171 scalar maxFaceCellRatio =
173 scalar totFaceCellRatio =
187 <<
setw(8) << totNprocs
189 <<
setw(8) << totNCells/totNprocs
190 <<
setw(8) << maxNCells
192 <<
setw(8) << totFaceCellRatio/totNprocs
193 <<
setw(8) << maxFaceCellRatio
195 <<
setw(8) << scalar(totNInt)/totNprocs
196 <<
setw(8) << maxNInt
198 <<
setw(8) << totRatio/totNprocs
199 <<
setw(8) << maxRatio
200 <<
setw(12) << totProfile/totNprocs
212 const label nFineCells,
213 const label nCoarseCells
224 return nTotalCoarseCells < nTotalFineCells;
241 nCellsInCoarsestLevel_
243 controlDict.getOrDefault<label>(
"nCellsInCoarsestLevel", 10)
262 restrictAddressing_(maxLevels_),
264 faceRestrictAddressing_(maxLevels_),
265 faceFlipMap_(maxLevels_),
266 nPatchFaces_(maxLevels_),
267 patchFaceRestrictAddressing_(maxLevels_),
269 meshLevels_(maxLevels_)
273 nCellsInCoarsestLevel_ =
279 procCommunicator_.setSize(maxLevels_ + 1, -1);
280 if (processorAgglomerate())
282 procAgglomMap_.setSize(maxLevels_);
283 agglomProcIDs_.setSize(maxLevels_);
284 procCellOffsets_.setSize(maxLevels_);
285 procFaceMap_.setSize(maxLevels_);
286 procBoundaryMap_.setSize(maxLevels_);
287 procBoundaryFaceMap_.setSize(maxLevels_);
301 GAMGAgglomeration::typeName
310 const word agglomeratorType
318 "geometricGAMGAgglomerationLibs",
319 lduMeshConstructorTablePtr_
322 auto* ctorPtr = lduMeshConstructorTable(agglomeratorType);
327 <<
"Unknown GAMGAgglomeration type "
328 << agglomeratorType <<
".\n"
329 <<
"Valid matrix GAMGAgglomeration types :"
330 << lduMatrixConstructorTablePtr_->sortedToc() <<
endl
331 <<
"Valid geometric GAMGAgglomeration types :"
332 << lduMeshConstructorTablePtr_->sortedToc()
352 GAMGAgglomeration::typeName
361 const word agglomeratorType
369 "algebraicGAMGAgglomerationLibs",
370 lduMatrixConstructorTablePtr_
373 auto* ctorPtr = lduMatrixConstructorTable(agglomeratorType);
399 GAMGAgglomeration::typeName
408 const word agglomeratorType
416 "geometricGAMGAgglomerationLibs",
417 geometryConstructorTablePtr_
420 auto* ctorPtr = geometryConstructorTable(agglomeratorType);
425 <<
"Unknown GAMGAgglomeration type "
426 << agglomeratorType <<
".\n"
427 <<
"Valid geometric GAMGAgglomeration types :"
428 << geometryConstructorTablePtr_->sortedToc()
465 return meshLevels_[i - 1];
478 return meshLevels_.set(i - 1);
490 return meshInterfaces_;
494 return meshLevels_[i - 1].rawInterfaces();
503 meshLevels_.set(i - 1,
nullptr);
505 if (i < nCells_.size())
508 restrictAddressing_.set(i,
nullptr);
510 faceRestrictAddressing_.set(i,
nullptr);
511 faceFlipMap_.set(i,
nullptr);
512 nPatchFaces_.set(i,
nullptr);
513 patchFaceRestrictAddressing_.set(i,
nullptr);
524 return procAgglomMap_[leveli];
533 return agglomProcIDs_[leveli];
539 return procCommunicator_[leveli] != -1;
545 return procCommunicator_[leveli];
554 return procCellOffsets_[leveli];
563 return procFaceMap_[leveli];
572 return procBoundaryMap_[leveli];
581 return procBoundaryFaceMap_[leveli];
594 if (fineAddressing.
size() != restriction.
size())
597 <<
"nCells:" << fineAddressing.
size()
598 <<
" agglom:" << restriction.
size()
615 const label own =
lower[facei];
616 const label nei =
upper[facei];
618 if (restriction[own] == restriction[nei])
622 if (master[own] < master[nei])
624 master[nei] = master[own];
627 else if (master[own] > master[nei])
629 master[own] = master[nei];
647 forAll(restriction, celli)
649 labelList& masters = coarseToMasters[restriction[celli]];
651 if (!masters.found(master[celli]))
653 masters.
append(master[celli]);
658 if (nNewCoarse > nCoarse)
669 nNewCoarse = nCoarse;
671 forAll(coarseToMasters, coarseI)
673 const labelList& masters = coarseToMasters[coarseI];
675 labelList& newCoarse = coarseToNewCoarse[coarseI];
676 newCoarse.
setSize(masters.size());
677 newCoarse[0] = coarseI;
678 for (label i=1; i<newCoarse.size(); i++)
680 newCoarse[i] = nNewCoarse++;
685 forAll(restriction, celli)
687 const label coarseI = restriction[celli];
689 const label index = coarseToMasters[coarseI].find(master[celli]);
690 newRestrict[celli] = coarseToNewCoarse[coarseI][index];
int debug
Static debugging option.
The class contains the addressing required by the lduMatrix: upper, lower and losort.
PtrList< labelListList > procFaceMap_
Mapping from processor to procMeshLevel face.
autoPtr< GAMGProcAgglomeration > procAgglomeratorPtr_
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Geometric agglomerated algebraic multigrid agglomeration class.
virtual const objectRegistry & thisDb() const
Return the object registry - resolve conflict polyMesh/lduMesh.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
labelList nCells_
The number of cells in each level.
lduMatrix is a general matrix class in which the coefficients are stored as three arrays,...
const lduMesh & meshLevel(const label leveli) const
Return LDU mesh of given level.
const labelListList & faceMap(const label fineLeveli) const
Mapping from processor to procMesh face.
void append(const T &val)
Append an element at the end of the list.
Ostream & endl(Ostream &os)
Add newline and flush stream.
PtrList< labelList > faceRestrictAddressing_
Face restriction addressing array.
bool hasMeshLevel(const label leveli) const
Do we have mesh for given level?
PtrList< lduPrimitiveMesh > meshLevels_
Hierarchy of mesh addressing.
label size() const
Return number of equations.
virtual const labelUList & upperAddr() const =0
Return upper addressing.
PtrList< labelField > restrictAddressing_
Cell restriction addressing array.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
const labelList & cellOffsets(const label fineLeveli) const
Mapping from processor to procMesh cells.
#define forAll(list, i)
Loop across all elements in list.
virtual lduInterfacePtrsList interfaces() const =0
const labelList & agglomProcIDs(const label fineLeveli) const
virtual label comm() const
Return communicator used for parallel communication.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
const labelListList & boundaryMap(const label fineLeveli) const
Mapping from processor to procMesh boundary.
string lower(const std::string &s)
Return string copy transformed with std::tolower on each character.
GAMGAgglomeration(const GAMGAgglomeration &)=delete
No copy construct.
OSstream & stream(OSstream *alternative=nullptr)
virtual lduInterfacePtrsList interfaces() const
Return a list of pointers for each patch.
messageStream Info
Information stream (stdout output on master, null elsewhere)
virtual const lduAddressing & lduAddr() const
Return ldu addressing.
void setSize(const label n)
Alias for resize()
runTime controlDict().readEntry("adjustTimeStep"
void clearLevel(const label leveli)
Istream and Ostream manipulators taking arguments.
static autoPtr< GAMGProcAgglomeration > New(const word &type, GAMGAgglomeration &agglom, const dictionary &controlDict)
Return the selected agglomerator.
static bool checkRestriction(labelList &newRestrict, label &nNewCoarse, const lduAddressing &fineAddressing, const labelUList &restriction, const label nCoarse)
Given restriction determines if coarse cells are connected.
PtrList< labelList > nPatchFaces_
The number of (coarse) patch faces in each level.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
dlLibraryTable & libs() const
Mutable access to the loaded dynamic libraries.
PtrList< labelList > procAgglomMap_
Per level, per processor the processor it agglomerates into.
errorManip< error > abort(error &err)
Omanip< int > setw(const int i)
static const GAMGAgglomeration & New(const lduMesh &mesh, const dictionary &controlDict)
Return the selected geometric agglomerator.
bool open(bool verbose=true)
errorManipArg< error, int > exit(error &err, const int errNo=1)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
virtual int precision() const
Get precision of output field.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const labelList & procAgglomMap(const label fineLeveli) const
const Type * cfindObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
labelList nFaces_
The number of (coarse) faces in each level.
const lduInterfacePtrsList & interfaceLevel(const label leveli) const
Return LDU interface addressing of given level.
virtual const labelUList & lowerAddr() const =0
Return lower addressing.
PtrList< labelListList > procBoundaryMap_
Mapping from processor to procMeshLevel boundary.
void compactLevels(const label nCreatedLevels)
Shrink the number of levels to that specified.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
const lduMesh & mesh() const
Return the LDU mesh from which the addressing is obtained.
const labelListListList & boundaryFaceMap(const label fineLeveli) const
Mapping from processor to procMesh boundary face.
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
bool hasProcMesh(const label fineLeveli) const
Check that level has combined mesh.
void size(const label n)
Older name for setAddressableSize.
label procCommunicator(const label fineLeveli) const
Communicator for current level or -1.
string upper(const std::string &s)
Return string copy transformed with std::toupper on each character.
label nCells(const label leveli) const
Return number of coarse cells (before processor agglomeration)
PtrList< boolList > faceFlipMap_
Face flip: for faces mapped to internal faces stores whether.
PtrList< labelListListList > procBoundaryFaceMap_
Mapping from processor to procMeshLevel boundary face.
defineTypeNameAndDebug(combustionModel, 0)
PtrList< labelListList > patchFaceRestrictAddressing_
Patch-local face restriction addressing array.
PtrList< labelList > procCellOffsets_
Mapping from processor to procMeshLevel cells.
virtual const lduAddressing & lduAddr() const =0
Return ldu addressing.
bool processorAgglomerate() const
Whether to agglomerate across processors.
const Time & time() const noexcept
Return time registry.
PtrList< labelList > agglomProcIDs_
Per level the set of processors to agglomerate. Element 0 is.
static label nProcs(const label communicator=worldComm)
Number of processes in parallel run, and 1 for serial run.
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Tuple2< label, scalar > band() const
Calculate bandwidth and profile of addressing.
~GAMGAgglomeration()
Destructor.
bool continueAgglomerating(const label nCells, const label nCoarseCells) const
Check the need for further agglomeration.
labelList procCommunicator_
Communicator for given level.