83static const scalar defaultMergeTol = 1
e-7;
93 const label masterMeshProcStart,
94 const label masterMeshProcEnd,
96 const label meshToAddProcStart,
97 const label meshToAddProcEnd,
99 const scalar mergeDist
102 if (fullMatch || masterMesh.
nCells() == 0)
127 forAll(masterPatches, patchi)
129 const polyPatch& pp = masterPatches[patchi];
131 if (isA<processorPolyPatch>(pp))
135 label proci=meshToAddProcStart;
136 proci<meshToAddProcEnd;
140 const string toProcString(
"to" +
name(proci));
142 pp.
name().rfind(toProcString)
143 == (pp.
name().size()-toProcString.size())
146 label meshFacei = pp.
start();
149 masterFaces.append(meshFacei++);
157 masterFaces.shrink();
171 forAll(addPatches, patchi)
173 const polyPatch& pp = addPatches[patchi];
175 if (isA<processorPolyPatch>(pp))
177 bool isConnected =
false;
181 label mergedProci=masterMeshProcStart;
182 !isConnected && (mergedProci < masterMeshProcEnd);
188 label proci = meshToAddProcStart;
189 proci < meshToAddProcEnd;
193 const word fromProcString
195 processorPolyPatch::newName(proci, mergedProci)
198 if (pp.
name() == fromProcString)
208 label meshFacei = pp.
start();
211 addFaces.append(meshFacei++);
236 const scalar mergeDist,
245 fvMeshAdder::findSharedPoints
252 Info<<
"mergeSharedPoints : detected " << pointToMaster.size()
253 <<
" points that are to be merged." <<
endl;
262 fvMeshAdder::mergePoints(
mesh, pointToMaster, meshMod);
268 mesh.updateMesh(map());
274 forAll(pointProcAddressing, proci)
276 labelList& constructMap = pointProcAddressing[proci];
280 label oldPointi = constructMap[i];
283 label
newPointi = map().reversePointMap()[oldPointi];
296 <<
"Problem. oldPointi:" << oldPointi
314 for (
const Time& procDb : databases)
325 if (pointsInstance != procDb.timeName())
328 <<
"Your time was specified as " << procDb.timeName()
329 <<
" but there is no polyMesh/points in that time." <<
nl
330 <<
"(points file at " << pointsInstance <<
')' <<
nl
331 <<
"Please rerun with the correct time specified"
332 <<
" (through the -constant, -time or -latestTime "
333 <<
"(at your option)."
337 Info<<
"Reading points from "
339 <<
" for time = " << procDb.timeName()
367void writeDistribution
390 forAll(cellProcAddressing, proci)
392 const labelList& pCells = cellProcAddressing[proci];
396 cellDecomposition.write();
398 Info<<
nl <<
"Wrote decomposition to "
399 << cellDecomposition.objectRelPath()
400 <<
" for use in manual decomposition." <<
endl;
406 const label oldIndex =
runTime.timeIndex();
409 runTime.setTime(0, oldIndex+1);
428 forAll(cellDecomposition, celli)
430 cellDist[celli] = cellDecomposition[celli];
433 cellDist.correctBoundaryConditions();
436 Info<<
nl <<
"Wrote decomposition to "
437 << cellDist.objectRelPath()
438 <<
" (volScalarField) for visualization."
456 /
mesh.time().timeName()
458 / polyMesh::meshSubDir
461 Info<<
nl <<
"Writing merged mesh to "
462 <<
mesh.time().relativePath(outputDir) <<
nl <<
endl;
467 <<
"Failed writing polyMesh."
470 topoSet::removeFiles(
mesh);
476 const label masterInternalFaces,
490 / polyMesh::meshSubDir
497 polyMesh::meshSubDir,
505 Info<<
"Writing addressing : " << outputDir <<
nl;
509 Info<<
" pointProcAddressing" <<
endl;
510 ioAddr.rename(
"pointProcAddressing");
514 Info<<
" faceProcAddressing" <<
endl;
515 ioAddr.rename(
"faceProcAddressing");
516 labelIOList faceProcAddr(ioAddr, faceProcAddressing);
520 forAll(faceProcAddr, procFacei)
522 const label masterFacei = faceProcAddr[procFacei];
527 && masterFacei < masterInternalFaces
533 label procOwn = procMesh.
faceOwner()[procFacei];
534 label masterOwn = masterOwner[masterFacei];
536 if (cellProcAddressing[procOwn] == masterOwn)
539 faceProcAddr[procFacei]++;
544 faceProcAddr[procFacei] = -1 - faceProcAddr[procFacei];
550 faceProcAddr[procFacei]++;
554 faceProcAddr.write();
558 Info<<
" cellProcAddressing" <<
endl;
559 ioAddr.rename(
"cellProcAddressing");
564 Info<<
" boundaryProcAddressing" <<
endl;
565 ioAddr.rename(
"boundaryProcAddressing");
575"Merge individual processor meshes back into one master mesh.\n"
576"Use if the original master mesh has been deleted or the processor meshes\n"
577"have been modified (topology change).\n"
578"This tool will write the resulting mesh to a new time step and construct\n"
579"xxxxProcAddressing files in the processor meshes so reconstructPar can be\n"
580"used to regenerate the fields on the master mesh.\n\n"
581"Not well tested & use at your own risk!\n\n";
585int main(
int argc,
char *argv[])
589 "Reconstruct a mesh using geometric/topological information only"
594 timeSelector::addOptions(
true,
true);
596 argList::noParallel();
598 argList::addVerboseOption(
"Additional verbosity");
599 argList::addBoolOption
602 "Create procAddressing only without overwriting the mesh"
608 "The merge distance relative to the bounding box size (default 1e-7)"
610 argList::addBoolOption
613 "Do (slower) geometric matching on all boundary faces"
615 argList::addBoolOption
618 "Do matching on processor faces only"
620 argList::addBoolOption
623 "Write cell distribution as a labelList - for use with 'manual' "
624 "decomposition method or as a volScalarField for post-processing."
634 const bool fullMatch =
args.
found(
"fullMatch");
635 const bool procMatch =
args.
found(
"procMatch");
636 const bool writeCellDist =
args.
found(
"cellDist");
638 const bool writeAddrOnly =
args.
found(
"addressing-only");
640 const scalar mergeTol =
645 Info<<
"Use geometric matching on all boundary faces." <<
nl <<
endl;
649 Info<<
"Use geometric matching on correct procBoundaries only." <<
nl
650 <<
"This assumes a correct decomposition." <<
endl;
654 Info<<
"Merge assuming correct, fully matched procBoundaries." <<
nl
658 if (fullMatch || procMatch)
660 const scalar writeTol =
661 Foam::pow(10.0, -scalar(IOstream::defaultPrecision()));
663 Info<<
"Merge tolerance : " << mergeTol <<
nl
664 <<
"Write tolerance : " << writeTol <<
endl;
666 if (
runTime.writeFormat() == IOstream::ASCII && mergeTol < writeTol)
669 <<
"Your current settings specify ASCII writing with "
670 << IOstream::defaultPrecision() <<
" digits precision." <<
endl
671 <<
"Your merging tolerance (" << mergeTol <<
")"
672 <<
" is finer than this." <<
endl
673 <<
"Please change your writeFormat to binary"
674 <<
" or increase the writePrecision" <<
endl
675 <<
"or adjust the merge tolerance (-mergeTol)."
689 <<
"No regions specified or detected."
692 else if (
regionNames[0] == polyMesh::defaultRegion)
709 <<
"No processor* directories found"
716 Info<<
"Found " << nProcs <<
" processor directories" <<
endl;
719 Info<<
" Reading database "
728 Time::controlDictName,
740 databases[0].times(),
748 runTime.setTime(timeDirs[timei], timei);
753 databases[proci].setTime(timeDirs[timei], timei);
759 label nMeshChanged = 0;
769 polyMesh::regionName(
regionName)/polyMesh::meshSubDir,
783 hasRegionMesh[regioni] =
true;
802 if (!hasRegionMesh[regioni])
822 label masterInternalFaces;
834 const scalar mergeDist = mergeTol*bb.
mag();
836 Info<<
"Overall mesh bounding box : " << bb <<
nl
837 <<
"Relative tolerance : " << mergeTol <<
nl
838 <<
"Absolute matching distance : " << mergeDist <<
nl
845 for (label proci=0; proci<nProcs; proci++)
877 boundProcAddressing[proci] =
902 renumber(map().addedCellMap(), cellProcAddressing[proci]);
903 renumber(map().addedFaceMap(), faceProcAddressing[proci]);
904 renumber(map().addedPointMap(), pointProcAddressing[proci]);
905 renumber(map().addedPatchMap(), boundProcAddressing[proci]);
907 for (label step=2; step<nProcs*2; step*=2)
909 for (label proci=0; proci<nProcs; proci+=step)
911 label next = proci + step/2;
917 Info<<
"Merging mesh " << proci <<
" with "
942 for (label mergedI=proci; mergedI<next; mergedI++)
947 cellProcAddressing[mergedI]
953 faceProcAddressing[mergedI]
959 pointProcAddressing[mergedI]
966 boundProcAddressing[mergedI]
974 addedI<
min(proci+step, nProcs);
980 map().addedCellMap(),
981 cellProcAddressing[addedI]
986 map().addedFaceMap(),
987 faceProcAddressing[addedI]
992 map().addedPointMap(),
993 pointProcAddressing[addedI]
998 map().addedPatchMap(),
999 boundProcAddressing[addedI]
1003 masterMesh.
set(next,
nullptr);
1007 for (label proci=0; proci<nProcs; proci++)
1009 Info<<
"Reading mesh to add from "
1010 << databases[proci].caseName()
1011 <<
" for time = " << databases[proci].timeName()
1017 mergeSharedPoints(mergeDist,masterMesh[0],pointProcAddressing);
1021 masterOwner = masterMesh[0].
faceOwner();
1027 <<
"Disabled writing of merged mesh (-addressing-only)"
1033 writeMesh(masterMesh[0], cellProcAddressing);
1049 for (label proci=0; proci<nProcs; proci++)
1078 const label nGlobalPatches = polyMeshAdder::procPatchPairs
1090 polyMeshAdder::patchFacePairs
1105 const labelList& procFaces = localBoundaryFace[proci];
1106 remoteFaceStart[proci].
setSize(procFaces.
size(), 0);
1118 patchMap[proci] =
identity(
mesh.boundaryMesh().size());
1121 patchMap[proci].setSize(nGlobalPatches);
1123 pointZoneMap[proci] =
identity(
mesh.pointZones().size());
1124 faceZoneMap[proci] =
identity(
mesh.faceZones().size());
1125 cellZoneMap[proci] =
identity(
mesh.cellZones().size());
1133 const labelList oldFaceOwner(fvMeshes[0].faceOwner());
1146 boundProcAddressing,
1159 const auto& pp = pbm[patchi];
1160 if (!isA<processorPolyPatch>(pp) || pp.size())
1162 oldToNew[patchi] = newi++;
1165 const label nNonProcPatches = newi;
1170 if (oldToNew[patchi] == -1)
1172 oldToNew[patchi] = newi++;
1175 fvMeshTools::reorderPatches
1183 masterMeshPtr = fvMeshes[0];
1187 const fvMesh& masterMesh = masterMeshPtr();
1207 <<
"Disabled writing of merged mesh (-addressing-only)"
1212 Time& masterTime =
const_cast<Time&
>(masterMesh.
time());
1219 writeMesh(masterMesh, cellProcAddressing);
1229 masterTime.
caseName() = oldCaseName;
1237 Info<<
"Reconstructing addressing from processor meshes"
1238 <<
" to the newly reconstructed mesh" <<
nl <<
endl;
1242 Info<<
"Processor " << proci <<
nl
1243 <<
"Read processor mesh: "
1244 << (databases[proci].caseName()/
regionDir) << endl;
1258 masterInternalFaces,
1261 cellProcAddressing[proci],
1262 faceProcAddressing[proci],
1263 pointProcAddressing[proci],
1264 boundProcAddressing[proci]
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
bool set(const Key &key, const T &obj)
Copy assign a new entry, overwriting existing entries.
A primitive field of type <T> with automated input and output.
A IOList wrapper for writing external data.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
void setSize(const label n)
Alias for resize()
A HashTable to objects of type <T> with a label key.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
const T * set(const label i) const
bool processorCase() const noexcept
Return true if this is a processor case.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
const fileName & caseName() const
Return case name.
A List with indirect addressing. Like IndirectList but does not store addressing.
void size(const label n)
Older name for setAddressableSize.
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
label size() const noexcept
The number of elements in the list.
const fileName & rootPath() const noexcept
Return root path.
bool found(const word &optName) const
Return true if the named option is found.
fileName path() const
Return the full path to the (processor local) case.
T getOrDefault(const word &optName, const T &deflt) const
Get a value from the named option if present, or return default.
const fileName & caseName() const noexcept
Return case name (parallel run) or global case (serial run)
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A bounding box defined in terms of min/max extrema points.
void add(const boundBox &bb)
Extend to include the second box.
scalar mag() const
The magnitude of the bounding box span.
A class for handling file names.
virtual label nProcs(const fileName &dir, const fileName &local="") const
Get number of processor directories/results. Used for e.g.
Sums a given list of (at least two or more) fields and outputs the result into a new field,...
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return the top-level database.
const Time & time() const noexcept
Return time registry.
const word & name() const noexcept
The patch name.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Mesh consisting of general polyhedral cells.
static const word & regionName(const word ®ion)
The mesh region name or word::null if polyMesh::defaultRegion.
const fileName & facesInstance() const
Return the current instance directory for faces.
virtual const labelList & faceOwner() const
Return face owner.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
A patch is a list of labels that address the faces in the global face list.
label start() const
Return start label of this patch in the polyMesh face list.
Direct mesh changes based on v1.3 polyTopoChange syntax.
bool isInternalFace(const label faceIndex) const noexcept
Return true if given face label is internal to the mesh.
label nInternalFaces() const noexcept
Number of internal faces.
label nPoints() const noexcept
Number of mesh points.
label nCells() const noexcept
Number of mesh cells.
label nFaces() const noexcept
Number of mesh faces.
A class for managing references or pointers (no reference counting)
virtual bool write(const bool valid=true) const
Write using setting from DB.
A class for handling words, derived from Foam::string.
Info<< "Creating field kinetic energy K\n"<< endl;volScalarField K("K", 0.5 *magSqr(U));if(U.nOldTimes()){ volVectorField *Uold=&U.oldTime();volScalarField *Kold=&K.oldTime();*Kold==0.5 *magSqr(*Uold);while(Uold->nOldTimes()) { Uold=&Uold-> oldTime()
Foam::word regionName(Foam::polyMesh::defaultRegion)
Foam::PtrList< Foam::fvMesh > meshes(regionNames.size())
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const fileOperation & fileHandler()
Get current file handler.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
messageStream Info
Information stream (stdout output on master, null elsewhere)
IntListType renumber(const labelUList &oldToNew, const IntListType &input)
Renumber the values (not the indices) of a list.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
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)
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.
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)
Foam::argList args(argc, argv)
#define forAll(list, i)
Loop across all elements in list.
static const char *const typeName
The type name used in ensight case files.