51void Foam::multiLevelDecomp::createMethodsDict()
73 nTotal = (domains.empty() ? 0 : 1);
75 for (
const label
n : domains)
98 domains.setSize(old.size()+1);
103 domains[i+1] = old[i];
105 nTotal *= domains[0];
107 Info<<
" inferred level0 with " << domains[0]
108 <<
" domains" <<
nl <<
nl;
111 if (!nLevels || nTotal !=
nDomains())
114 <<
"Top level decomposition specifies " <<
nDomains()
115 <<
" domains which is not equal to the product of"
116 <<
" all sub domains " << nTotal
124 const dictionary& subMethodCoeffsDict
129 defaultMethod +
"Coeffs",
134 for (
const label
n : domains)
136 const word levelName(
"level" +
Foam::name(nLevels++));
138 entry* dictptr = methodsDict_.set(levelName, dictionary());
140 dictionary&
dict = dictptr->dict();
141 dict.
add(
"method", defaultMethod);
145 if (subMethodCoeffsDict.size())
147 dict.
add(subMethodCoeffsDict.dictName(), subMethodCoeffsDict);
159 for (
const entry& dEntry : coeffsDict_)
172 const bool addDefaultMethod
175 && !defaultMethod.empty()
178 entry*
e = methodsDict_.add(dEntry);
180 if (addDefaultMethod &&
e &&
e->isDict())
182 e->dict().add(
"method", defaultMethod);
190void Foam::multiLevelDecomp::setMethods()
198 methods_.setSize(methodsDict_.size());
199 for (
const entry& dEntry : methodsDict_)
214 methods_.setSize(nLevels);
219 <<
"Decompose " <<
type() <<
" [" << nDomains() <<
"] in "
220 << nLevels <<
" levels:" <<
endl;
225 Info<<
" level " << i <<
" : " << methods_[i].type()
226 <<
" [" << methods_[i].nDomains() <<
"]" <<
endl;
228 nTotal *= methods_[i].nDomains();
231 if (nTotal != nDomains())
234 <<
"Top level decomposition specifies " << nDomains()
235 <<
" domains which is not equal to the product of"
236 <<
" all sub domains " << nTotal
244void Foam::multiLevelDecomp::subsetGlobalCellCells
246 const label nDomains,
256 const globalIndex globalCells(cellCells.size());
262 subCellCells = UIndirectList<labelList>(cellCells, set);
265 List<Map<label>> compactMap;
266 mapDistribute map(globalCells, subCellCells, compactMap);
267 map.distribute(oldToNew);
269 map.distribute(allDist);
278 const globalIndex globalSubCells(
set.size());
283 cutConnections.setSize(nDomains);
286 forAll(subCellCells, subCelli)
288 labelList& cCells = subCellCells[subCelli];
295 const label nbrCelli = oldToNew[cCells[i]];
298 cutConnections[allDist[cCells[i]]]++;
305 const label celli =
set[subCelli];
306 const label oldNbrCelli = cellCells[celli][i];
308 const label proci = globalCells.whichProcID(oldNbrCelli);
310 cCells[newI++] = globalSubCells.toGlobal(proci, nbrCelli);
313 cCells.setSize(newI);
324 const label currLevel,
325 const label leafOffset,
332 methods_[currLevel].decompose
341 const label nextLevel = currLevel+1;
344 const label nCurrDomains = methods_[currLevel].nDomains();
354 for (label i = 0; i <= currLevel; ++i)
356 sizes *= methods_[i].nDomains();
360 sizes = this->nDomains() / sizes;
364 domainLookup[i] = i * sizes + leafOffset;
370 Info<<
"Distribute at level " << currLevel
371 <<
" to domains" <<
nl
378 const label orig = pointMap[i];
379 finalDecomp[orig] = domainLookup[dist[i]];
382 if (nextLevel < methods_.size())
392 Pout<<
"Decomposition at level " << currLevel <<
" :" <<
endl;
395 for (label domainI = 0; domainI < nCurrDomains; ++domainI)
402 scalarField subWeights(pointWeights, domainPoints);
407 subsetGlobalCellCells
429 for (
const label nConnect : nOutsideConnections)
441 Pout<<
" Domain " << domainI <<
nl
442 <<
" Number of cells = " <<
nPoints <<
nl
443 <<
" Number of inter-domain patches = " <<
nPatches
445 <<
" Number of inter-domain faces = " <<
nFaces <<
nl
458 domainLookup[domainI],
472 const label nNext = methods_[nextLevel].nDomains();
473 const label nTotal = nCurrDomains * nNext;
476 dictionary level0Dict;
477 for (
const entry& dEntry : methodsDict_)
481 level0Dict = dEntry.dict();
485 level0Dict.set(
"numberOfSubdomains", nTotal);
489 Pout<<
"Reference decomposition with " << level0Dict <<
" :"
507 for (label blockI = 0; blockI < nCurrDomains; ++blockI)
513 forAll(pointPoints, pointi)
515 if ((dist[pointi] / nNext) == blockI)
519 const labelList& pPoints = pointPoints[pointi];
523 const label distBlockI = dist[pPoints[i]] / nNext;
524 if (distBlockI != blockI)
526 nOutsideConnections[distBlockI]++;
541 for (
const label nConnect : nOutsideConnections)
552 Pout<<
" Domain " << blockI <<
nl
553 <<
" Number of cells = " <<
nPoints <<
nl
554 <<
" Number of inter-domain patches = "
556 <<
" Number of inter-domain faces = " <<
nFaces
596 if (!meth.parallelAware())
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A packed storage unstructured matrix of objects of type <T> using an offset table for access.
List< SubListType > unpack() const
Return non-compact list of lists.
static void listCombineAllGather(const List< commsStruct > &comms, List< T > &values, const CombineOp &cop, const int tag, const label comm)
After completion all processors have the same data.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
void size(const label n)
Older name for setAddressableSize.
Abstract base class for domain decomposition.
selectionType
Selection type when handling the coefficients dictionary.
static const dictionary & findCoeffsDict(const dictionary &dict, const word &coeffsName, int select=selectionType::DEFAULT)
label nDomains() const noexcept
Number of domains.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
void clear()
Clear the dictionary.
entry * add(entry *entryPtr, bool mergeEntry=false)
Add a new entry.
entry * set(entry *entryPtr)
Assign a new entry, overwriting any existing entry.
Decompose given using consecutive application of decomposers.
virtual bool parallelAware() const
Is method parallel aware?
Mesh consisting of general polyhedral cells.
const string & prefix() const noexcept
Return the stream prefix.
splitCell * master() const
bool decompose() const noexcept
Query the decompose flag (normally off)
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Foam::word regionName(Foam::polyMesh::defaultRegion)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const labelList nFaces(UPstream::listGatherValues< label >(aMesh.nFaces()))
void set(List< bool > &bools, const labelUList &locations)
Set the listed locations (assign 'true').
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
List< label > labelList
A List of labels.
messageStream Info
Information stream (stdout output on master, null elsewhere)
vectorField pointField
pointField is a vectorField.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
List< labelList > labelListList
A List of labelList.
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
static constexpr const zero Zero
Global zero (0)
labelList invert(const label len, const labelUList &map)
Create an inverse one-to-one mapping.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
labelListList invertOneToMany(const label len, const labelUList &map)
Invert one-to-many map. Unmapped elements will be size 0.
labelList findIndices(const ListType &input, typename ListType::const_reference val, label start=0)
Linear search to find all occurrences of given element.
errorManipArg< error, int > exit(error &err, const int errNo=1)
UList< label > labelUList
A UList of labels.
UIndirectList< label > labelUIndList
UIndirectList of labels.
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.