50 nCells_.setSize(nCreatedLevels);
51 restrictAddressing_.setSize(nCreatedLevels);
52 nFaces_.setSize(nCreatedLevels);
53 faceRestrictAddressing_.setSize(nCreatedLevels);
54 faceFlipMap_.setSize(nCreatedLevels);
55 nPatchFaces_.setSize(nCreatedLevels);
56 patchFaceRestrictAddressing_.setSize(nCreatedLevels);
57 meshLevels_.setSize(nCreatedLevels);
60 procCommunicator_.setSize(nCreatedLevels + 1);
61 if (processorAgglomerate())
63 procAgglomMap_.setSize(nCreatedLevels);
64 agglomProcIDs_.setSize(nCreatedLevels);
65 procCellOffsets_.setSize(nCreatedLevels);
66 procFaceMap_.setSize(nCreatedLevels);
67 procBoundaryMap_.setSize(nCreatedLevels);
68 procBoundaryFaceMap_.setSize(nCreatedLevels);
70 procAgglomeratorPtr_().agglomerate();
76 if (processorAgglomerate() && debug)
78 Info<<
"GAMGAgglomeration:" <<
nl 79 <<
" local agglomerator : " <<
type() <<
nl;
80 if (processorAgglomerate())
82 Info<<
" processor agglomerator : " 83 << procAgglomeratorPtr_().type() << nl
88 <<
setw(20) <<
"nFaces/nCells" 89 <<
setw(20) <<
"nInterfaces" 90 <<
setw(20) <<
"nIntFaces/nCells" 91 <<
setw(12) <<
"profile" 94 <<
setw(8) <<
"nProcs" 110 <<
setw(8) <<
"-----" 111 <<
setw(8) <<
"------" 129 for (
label levelI = 0; levelI <= size(); levelI++)
133 scalar faceCellRatio = 0;
134 label nInterfaces = 0;
137 scalar profile = 0.0;
139 if (hasMeshLevel(levelI))
143 const lduMesh& fineMesh = meshLevel(levelI);
152 if (interfaces.set(i))
155 nIntFaces += interfaces[i].faceCells().size();
158 ratio = scalar(nIntFaces)/nCells;
168 scalar maxFaceCellRatio =
170 scalar totFaceCellRatio =
181 int oldPrecision =
Info().precision(4);
184 <<
setw(8) << totNprocs
186 <<
setw(8) << totNCells/totNprocs
187 <<
setw(8) << maxNCells
189 <<
setw(8) << totFaceCellRatio/totNprocs
190 <<
setw(8) << maxFaceCellRatio
192 <<
setw(8) << scalar(totNInt)/totNprocs
193 <<
setw(8) << maxNInt
195 <<
setw(8) << totRatio/totNprocs
196 <<
setw(8) << maxRatio
197 <<
setw(12) << totProfile/totNprocs
200 Info().precision(oldPrecision);
209 const label nCoarseCells
213 bool contAgg = nCoarseCells >= nCellsInCoarsestLevel_;
221 Foam::GAMGAgglomeration::GAMGAgglomeration
231 nCellsInCoarsestLevel_
240 && controlDict.
found(
"processorAgglomerator")
244 controlDict.
lookup(
"processorAgglomerator"),
252 restrictAddressing_(maxLevels_),
254 faceRestrictAddressing_(maxLevels_),
255 faceFlipMap_(maxLevels_),
256 nPatchFaces_(maxLevels_),
257 patchFaceRestrictAddressing_(maxLevels_),
259 meshLevels_(maxLevels_)
261 procCommunicator_.setSize(maxLevels_ + 1, -1);
262 if (processorAgglomerate())
264 procAgglomMap_.setSize(maxLevels_);
265 agglomProcIDs_.setSize(maxLevels_);
266 procCellOffsets_.setSize(maxLevels_);
267 procFaceMap_.setSize(maxLevels_);
268 procBoundaryMap_.setSize(maxLevels_);
269 procBoundaryFaceMap_.setSize(maxLevels_);
284 GAMGAgglomeration::typeName
288 const word agglomeratorType
296 "geometricGAMGAgglomerationLibs",
297 lduMeshConstructorTablePtr_
300 lduMeshConstructorTable::iterator cstrIter =
301 lduMeshConstructorTablePtr_->
find(agglomeratorType);
303 if (cstrIter == lduMeshConstructorTablePtr_->end())
307 "GAMGAgglomeration::New" 308 "(const lduMesh& mesh, const dictionary& controlDict)" 309 ) <<
"Unknown GAMGAgglomeration type " 310 << agglomeratorType <<
".\n" 311 <<
"Valid matrix GAMGAgglomeration types are :" 312 << lduMatrixConstructorTablePtr_->sortedToc() <<
endl 313 <<
"Valid geometric GAMGAgglomeration types are :" 314 << lduMeshConstructorTablePtr_->sortedToc()
318 return store(cstrIter()(mesh, controlDict).ptr());
324 GAMGAgglomeration::typeName
342 GAMGAgglomeration::typeName
346 const word agglomeratorType
354 "algebraicGAMGAgglomerationLibs",
355 lduMatrixConstructorTablePtr_
360 !lduMatrixConstructorTablePtr_
361 || !lduMatrixConstructorTablePtr_->found(agglomeratorType)
364 return New(mesh, controlDict);
368 lduMatrixConstructorTable::iterator cstrIter =
369 lduMatrixConstructorTablePtr_->find(agglomeratorType);
371 return store(cstrIter()(matrix, controlDict).ptr());
378 GAMGAgglomeration::typeName
392 const word agglomeratorType
400 "geometricGAMGAgglomerationLibs",
401 geometryConstructorTablePtr_
404 geometryConstructorTable::iterator cstrIter =
405 geometryConstructorTablePtr_->
find(agglomeratorType);
407 if (cstrIter == geometryConstructorTablePtr_->end())
411 "GAMGAgglomeration::New" 412 "(const lduMesh& mesh, const scalarField&" 413 ", const vectorField&, const dictionary& controlDict)" 414 ) <<
"Unknown GAMGAgglomeration type " 415 << agglomeratorType <<
".\n" 416 <<
"Valid geometric GAMGAgglomeration types are :" 417 << geometryConstructorTablePtr_->sortedToc()
453 return meshLevels_[i - 1];
466 return meshLevels_.set(i - 1);
478 return meshInterfaces_;
482 return meshLevels_[i - 1].rawInterfaces();
491 meshLevels_.set(i - 1, NULL);
493 if (i < nCells_.size())
496 restrictAddressing_.set(i, NULL);
498 faceRestrictAddressing_.set(i, NULL);
499 faceFlipMap_.set(i, NULL);
500 nPatchFaces_.set(i, NULL);
501 patchFaceRestrictAddressing_.set(i, NULL);
512 return procAgglomMap_[leveli];
521 return agglomProcIDs_[leveli];
527 return procCommunicator_[leveli] != -1;
533 return procCommunicator_[leveli];
542 return procCellOffsets_[leveli];
551 return procFaceMap_[leveli];
560 return procBoundaryMap_[leveli];
569 return procBoundaryFaceMap_[leveli];
582 if (fineAddressing.
size() != restrict.
size())
586 "checkRestriction(..)" 587 ) <<
"nCells:" << fineAddressing.
size()
588 <<
" agglom:" << restrict.
size()
605 label own = lower[faceI];
606 label nei = upper[faceI];
608 if (restrict[own] == restrict[nei])
612 if (master[own] < master[nei])
614 master[nei] = master[own];
617 else if (master[own] > master[nei])
619 master[own] = master[nei];
639 labelList& masters = coarseToMasters[restrict[cellI]];
641 if (
findIndex(masters, master[cellI]) == -1)
643 masters.
append(master[cellI]);
648 if (nNewCoarse > nCoarse)
659 nNewCoarse = nCoarse;
661 forAll(coarseToMasters, coarseI)
663 const labelList& masters = coarseToMasters[coarseI];
665 labelList& newCoarse = coarseToNewCoarse[coarseI];
667 newCoarse[0] = coarseI;
668 for (
label i = 1; i < newCoarse.
size(); i++)
670 newCoarse[i] = nNewCoarse++;
677 label coarseI = restrict[cellI];
680 newRestrict[cellI] = coarseToNewCoarse[coarseI][index];
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
bool continueAgglomerating(const label nCoarseCells) const
Check the need for further agglomeration.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
label size() const
Return number of equations.
virtual const labelUList & upperAddr() const =0
Return upper addressing.
const labelListList & boundaryMap(const label fineLeveli) const
Mapping from processor to procMesh boundary.
void reduce(T &Value, const BinaryOp &bop) const
Helper: reduce with current communicator.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
virtual label comm() const =0
Return communicator used for parallel communication.
const labelList & cellOffsets(const label fineLeveli) const
Mapping from processor to procMesh cells.
Geometric agglomerated algebraic multigrid agglomeration class.
bool foundObject(const word &name) const
Is the named Type found?
static autoPtr< GAMGProcAgglomeration > New(const word &type, GAMGAgglomeration &agglom, const dictionary &controlDict)
Return the selected agglomerator.
virtual lduInterfacePtrsList interfaces() const =0
Return a list of pointers for each patch.
static const GAMGAgglomeration & New(const lduMesh &mesh, const dictionary &controlDict)
Return the selected geometric agglomerator.
A class for handling words, derived from string.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
void size(const label)
Override size to be inconsistent with allocated storage.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
static bool checkRestriction(labelList &newRestrict, label &nNewCoarse, const lduAddressing &fineAddressing, const labelUList &restrict, const label nCoarse)
Given restriction determines if coarse cells are connected.
A list of keyword definitions, which are a keyword followed by any number of values (e...
virtual const objectRegistry & thisDb() const
Return the object registry.
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
const labelList & procAgglomMap(const label fineLeveli) const
Mapping from processor to agglomerated processor (global, all.
virtual const lduAddressing & lduAddr() const =0
Return ldu addressing.
const lduMesh & meshLevel(const label leveli) const
Return LDU mesh of given level.
Tuple2< label, scalar > band() const
Calculate bandwidth and profile of addressing.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
const labelListList & faceMap(const label fineLeveli) const
Mapping from processor to procMesh face.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
void setSize(const label)
Reset size of List.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const labelListListList & boundaryFaceMap(const label fineLeveli) const
Mapping from processor to procMesh boundary face.
bool hasMeshLevel(const label leveli) const
Do we have mesh for given level?
lduMatrix is a general matrix class in which the coefficients are stored as three arrays...
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
label size() const
Return the number of elements in the UList.
virtual const labelUList & lowerAddr() const =0
Return lower addressing.
Omanip< int > setw(const int i)
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
errorManip< error > abort(error &err)
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Istream and Ostream manipulators taking arguments.
void compactLevels(const label nCreatedLevels)
Shrink the number of levels to that specified.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
const labelList & agglomProcIDs(const label fineLeveli) const
Set of processors to agglomerate. Element 0 is the.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
void clearLevel(const label leveli)
runTime controlDict().lookup("adjustTimeStep") >> adjustTimeStep
bool hasProcMesh(const label fineLeveli) const
Check that level has combined mesh.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
const lduInterfacePtrsList & interfaceLevel(const label leveli) const
Return LDU interface addressing of given level.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
const lduMesh & mesh() const
Return the LDU mesh from which the addressing is obtained.
void append(const T &)
Append an element at the end of the list.
const Time & time() const
Return time.
The class contains the addressing required by the lduMatrix: upper, lower and losort.
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
~GAMGAgglomeration()
Destructor.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
label procCommunicator(const label fineLeveli) const
Communicator for current level or -1.
defineTypeNameAndDebug(combustionModel, 0)