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 =
184 int oldPrecision =
Info().precision(4);
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
203 Info().precision(oldPrecision);
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_);
302 GAMGAgglomeration::typeName
306 const word agglomeratorType
314 "geometricGAMGAgglomerationLibs",
315 lduMeshConstructorTablePtr_
318 auto cstrIter = lduMeshConstructorTablePtr_->cfind(agglomeratorType);
320 if (!cstrIter.found())
323 <<
"Unknown GAMGAgglomeration type "
324 << agglomeratorType <<
".\n"
325 <<
"Valid matrix GAMGAgglomeration types :"
326 << lduMatrixConstructorTablePtr_->sortedToc() <<
endl
327 <<
"Valid geometric GAMGAgglomeration types :"
328 << lduMeshConstructorTablePtr_->sortedToc()
338 GAMGAgglomeration::typeName
356 GAMGAgglomeration::typeName
360 const word agglomeratorType
368 "algebraicGAMGAgglomerationLibs",
369 lduMatrixConstructorTablePtr_
374 !lduMatrixConstructorTablePtr_
375 || !lduMatrixConstructorTablePtr_->found(agglomeratorType)
383 lduMatrixConstructorTablePtr_->cfind(agglomeratorType);
385 return store(cstrIter()(matrix,
controlDict).ptr());
392 GAMGAgglomeration::typeName
406 const word agglomeratorType
414 "geometricGAMGAgglomerationLibs",
415 geometryConstructorTablePtr_
418 auto cstrIter = geometryConstructorTablePtr_->cfind(agglomeratorType);
420 if (!cstrIter.found())
423 <<
"Unknown GAMGAgglomeration type "
424 << agglomeratorType <<
".\n"
425 <<
"Valid geometric GAMGAgglomeration types :"
426 << geometryConstructorTablePtr_->sortedToc()
462 return meshLevels_[i - 1];
475 return meshLevels_.set(i - 1);
487 return meshInterfaces_;
491 return meshLevels_[i - 1].rawInterfaces();
500 meshLevels_.set(i - 1,
nullptr);
502 if (i < nCells_.size())
505 restrictAddressing_.set(i,
nullptr);
507 faceRestrictAddressing_.set(i,
nullptr);
508 faceFlipMap_.set(i,
nullptr);
509 nPatchFaces_.set(i,
nullptr);
510 patchFaceRestrictAddressing_.set(i,
nullptr);
521 return procAgglomMap_[leveli];
530 return agglomProcIDs_[leveli];
536 return procCommunicator_[leveli] != -1;
542 return procCommunicator_[leveli];
551 return procCellOffsets_[leveli];
560 return procFaceMap_[leveli];
569 return procBoundaryMap_[leveli];
578 return procBoundaryFaceMap_[leveli];
591 if (fineAddressing.
size() != restriction.
size())
594 <<
"nCells:" << fineAddressing.
size()
595 <<
" agglom:" << restriction.
size()
612 label own =
lower[facei];
613 label nei =
upper[facei];
615 if (restriction[own] == restriction[nei])
619 if (master[own] < master[nei])
621 master[nei] = master[own];
624 else if (master[own] > master[nei])
626 master[own] = master[nei];
644 forAll(restriction, celli)
646 labelList& masters = coarseToMasters[restriction[celli]];
648 if (!masters.found(master[celli]))
650 masters.
append(master[celli]);
655 if (nNewCoarse > nCoarse)
666 nNewCoarse = nCoarse;
668 forAll(coarseToMasters, coarseI)
670 const labelList& masters = coarseToMasters[coarseI];
672 labelList& newCoarse = coarseToNewCoarse[coarseI];
673 newCoarse.
setSize(masters.size());
674 newCoarse[0] = coarseI;
675 for (label i=1; i<newCoarse.size(); i++)
677 newCoarse[i] = nNewCoarse++;
682 forAll(restriction, celli)
684 const label coarseI = restriction[celli];
686 label index = coarseToMasters[coarseI].find(master[celli]);
687 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_
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 Time & time() const
Return time.
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.
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
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
Return a list of pointers for each patch.
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.
messageStream Info
Information stream (uses stdout - output is on the master only)
virtual const lduAddressing & lduAddr() const
Return ldu addressing.
runTime controlDict().readEntry("adjustTimeStep"
void clearLevel(const label leveli)
Istream and Ostream manipulators taking arguments.
const Type & lookupObject(const word &name, const bool recursive=false) const
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.
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
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.
virtual lduInterfacePtrsList interfaces() const
Return a list of pointers for each patch.
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.
void size(const label n) noexcept
Override size to be inconsistent with allocated storage.
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.
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.
void setSize(const label newSize)
Alias for resize(const label)
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.
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.