89 label diff = neighbour[facei] - owner[facei];
103 const bool calculateIntersect,
109 scalar& sumSqrIntersect
117 label own = owner[facei];
118 label nei = neighbour[facei];
122 cellBandwidth[nei] =
max(cellBandwidth[nei],
diff);
125 bandwidth =
max(cellBandwidth);
129 forAll(cellBandwidth, celli)
131 profile += 1.0*cellBandwidth[celli];
134 sumSqrIntersect = 0.0;
135 if (calculateIntersect)
139 for (
label colI = celli-cellBandwidth[celli]; colI <= celli; colI++)
141 nIntersect[colI] += 1.0;
166 forAll(cellOrder, newCelli)
168 label oldCelli = cellOrder[newCelli];
170 const cell& cFaces = mesh.
cells()[oldCelli];
177 label facei = cFaces[i];
183 if (nbrCelli == newCelli)
185 nbrCelli = reverseCellOrder[mesh.
faceOwner()[facei]];
188 if (newCelli < nbrCelli)
212 if (nbr[index] != -1)
214 oldToNewFace[cFaces[index]] = newFacei++;
220 for (
label facei = newFacei; facei < mesh.
nFaces(); facei++)
222 oldToNewFace[facei] = facei;
227 forAll(oldToNewFace, facei)
229 if (oldToNewFace[facei] == -1)
232 <<
"Did not determine new position" <<
" for face " << facei
256 label prevRegion = -1;
258 forAll(cellOrder, newCelli)
260 label oldCelli = cellOrder[newCelli];
262 if (cellToRegion[oldCelli] != prevRegion)
264 prevRegion = cellToRegion[oldCelli];
267 const cell& cFaces = mesh.
cells()[oldCelli];
273 label facei = cFaces[i];
279 if (nbrCelli == newCelli)
281 nbrCelli = reverseCellOrder[mesh.
faceOwner()[facei]];
284 if (cellToRegion[oldCelli] != cellToRegion[cellOrder[nbrCelli]])
289 else if (newCelli < nbrCelli)
313 oldToNewFace[cFaces[nbr.indices()[i]]] = newFacei++;
319 label nRegions =
max(cellToRegion)+1;
329 if (ownRegion != neiRegion)
332 min(ownRegion, neiRegion)*nRegions
333 +
max(ownRegion, neiRegion);
342 label key = sortKey[i];
354 oldToNewFace[sortKey.indices()[i]] = newFacei++;
359 for (
label facei = newFacei; facei < mesh.
nFaces(); facei++)
361 oldToNewFace[facei] = facei;
366 forAll(oldToNewFace, facei)
368 if (oldToNewFace[facei] == -1)
371 <<
"Did not determine new position"
372 <<
" for face " << facei
414 forAll(newNeighbour, facei)
416 label own = newOwner[facei];
417 label nei = newNeighbour[facei];
421 newFaces[facei].flip();
422 Swap(newOwner[facei], newNeighbour[facei]);
423 flipFaceFlux.insert(facei);
443 NullObjectMove<pointField>(),
467 move(reverseFaceOrder),
468 move(reverseCellOrder),
473 move(oldPatchNMeshPoints),
488 Info<<
"Determining cell order:" <<
endl;
492 label nRegions =
max(cellToRegion)+1;
498 forAll(regionToCells, regionI)
500 Info<<
" region " << regionI <<
" starts at " << celli <<
endl;
508 subsetter.setLargeCellSubset(cellToRegion, regionI);
510 const fvMesh& subMesh = subsetter.subMesh();
521 const labelList& cellMap = subsetter.cellMap();
525 cellOrder[celli++] = cellMap[subCellOrder[i]];
536 int main(
int argc,
char *argv[])
540 "Renumber mesh to minimise bandwidth"
550 "calculate the rms of the frontwidth"
555 "do not update fields"
574 scalar sumSqrIntersect;
598 <<
"Before renumbering :" <<
nl
599 <<
" band : " << band <<
nl
600 <<
" profile : " << profile <<
nl;
604 Info<<
" rms frontwidth : " << rmsFrontwidth <<
nl;
609 bool sortCoupledFaceCells =
false;
610 bool writeMaps =
false;
611 bool orderPoints =
false;
613 bool renumberSets =
true;
621 renumberDictPtr.
reset
625 const dictionary& renumberDict = renumberDictPtr();
631 "sortCoupledFaceCells",
634 if (sortCoupledFaceCells)
636 Info<<
"Sorting cells on coupled boundaries to be last." <<
nl
643 Info<<
"Ordering cells into regions of size " << blockSize
644 <<
" (using decomposition);"
645 <<
" ordering faces into region-internal and region-external."
648 if (blockSize < 0 || blockSize >= mesh.
nCells())
651 <<
"Block size " << blockSize
652 <<
" should be positive integer"
653 <<
" and less than the number of cells in the mesh."
661 Info<<
"Ordering points into internal and boundary points." <<
nl
665 renumberDict.
lookup(
"writeMaps") >> writeMaps;
668 Info<<
"Writing renumber maps (new to old) to polyMesh." <<
nl
676 Info<<
"Using default renumberMethod." <<
nl <<
endl;
681 Info<<
"Selecting renumberMethod " << renumberPtr().type() <<
nl <<
endl;
690 "cellProcAddressing",
703 "faceProcAddressing",
715 "pointProcAddressing",
747 Info<<
"nBlocks = " << nBlocks <<
endl;
750 dictionary decomposeDict(renumberDictPtr().subDict(
"blockCoeffs"));
751 decomposeDict.set(
"numberOfSubdomains", nBlocks);
763 decomposePtr().decompose
780 Info<<
nl <<
"Written decomposition as volScalarField::Internal to "
781 <<
"cellProc for use in postprocessing."
784 cellOrder = regionRenumber(renumberPtr(), mesh, cellToRegion);
787 faceOrder = getRegionFaceOrder
797 cellOrder = renumberPtr().renumber
803 if (sortCoupledFaceCells)
812 if (pbm[
patchi].coupled())
825 if (pbm[
patchi].coupled())
830 label celli = faceCells[i];
832 if (reverseCellOrder[celli] != -1)
834 bndCells[nBndCells] = celli;
835 bndCellMap[nBndCells++] = reverseCellOrder[celli];
836 reverseCellOrder[celli] = -1;
842 bndCellMap.setSize(nBndCells);
855 newReverseCellOrder[origCelli] = --sortedI;
858 Info<<
"Ordered all " << nBndCells <<
" cells with a coupled face"
859 <<
" to the end of the cell list, starting at " << sortedI
864 forAll(cellOrder, newCelli)
866 label origCelli = cellOrder[newCelli];
867 if (newReverseCellOrder[origCelli] == -1)
869 newReverseCellOrder[origCelli] = sortedI++;
879 faceOrder = getFaceOrder
912 pointOrderMap().pointMap()
917 pointOrderMap().reversePointMap(),
918 const_cast<labelList&
>(map().reversePointMap())
929 cellProcAddressing.headerOk()
930 && cellProcAddressing.size() == mesh.
nCells()
933 Info<<
"Renumbering processor cell decomposition map "
934 << cellProcAddressing.name() <<
endl;
943 faceProcAddressing.headerOk()
944 && faceProcAddressing.size() == mesh.
nFaces()
947 Info<<
"Renumbering processor face decomposition map "
948 << faceProcAddressing.name() <<
endl;
959 label facei = iter.key();
960 label masterFacei = faceProcAddressing[facei];
962 faceProcAddressing[facei] = -masterFacei;
964 if (masterFacei == 0)
973 pointProcAddressing.headerOk()
974 && pointProcAddressing.size() == mesh.
nPoints()
977 Info<<
"Renumbering processor point decomposition map "
978 << pointProcAddressing.name() <<
endl;
990 scalar sumSqrIntersect;
1012 Info<<
"After renumbering :" <<
nl
1013 <<
" band : " << band <<
nl
1014 <<
" profile : " << profile <<
nl;
1019 Info<<
" rms frontwidth : " << rmsFrontwidth <<
nl;
1064 <<
" total : " << nTotPoints <<
nl
1065 <<
" internal: " << nTotIntPoints <<
nl
1066 <<
" boundary: " << nTotPoints-nTotIntPoints <<
nl
1068 <<
" total : " << nTotEdges <<
nl
1069 <<
" internal: " << nTotIntEdges <<
nl
1070 <<
" internal using 0 boundary points: "
1071 << nTotInt0Edges <<
nl
1072 <<
" internal using 1 boundary points: "
1073 << nTotInt1Edges-nTotInt0Edges <<
nl
1074 <<
" internal using 2 boundary points: "
1075 << nTotIntEdges-nTotInt1Edges <<
nl
1076 <<
" boundary: " << nTotEdges-nTotIntEdges <<
nl
1090 if (cellProcAddressing.headerOk())
1093 if (cellProcAddressing.size() == mesh.
nCells())
1095 cellProcAddressing.write();
1100 const fileName fName(cellProcAddressing.filePath());
1103 Info<<
"Deleting inconsistent processor cell decomposition"
1104 <<
" map " << fName <<
endl;
1110 if (faceProcAddressing.headerOk())
1113 if (faceProcAddressing.size() == mesh.
nFaces())
1115 faceProcAddressing.write();
1119 const fileName fName(faceProcAddressing.filePath());
1122 Info<<
"Deleting inconsistent processor face decomposition"
1123 <<
" map " << fName <<
endl;
1129 if (pointProcAddressing.headerOk())
1132 if (pointProcAddressing.size() == mesh.
nPoints())
1134 pointProcAddressing.write();
1138 const fileName fName(pointProcAddressing.filePath());
1141 Info<<
"Deleting inconsistent processor point decomposition"
1142 <<
" map " << fName <<
endl;
1165 Info<<
nl <<
"Written current cellID and origCellID as volScalarField"
1166 <<
" for use in postprocessing."
1228 Info<<
"Renumbering cellSets:" <<
endl;
1233 cs.topoChange(map());
1244 Info<<
"Renumbering faceSets:" <<
endl;
1249 fs.topoChange(map());
1260 Info<<
"Renumbering pointSets:" <<
endl;
1265 ps.topoChange(map());
Field reading functions for post-processing utilities.
#define forAll(list, i)
Loop across all elements in list.
#define forAllReverse(list, i)
Reverse loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Cuthill-McKee renumbering.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
List of IOobjects with searching and retrieving facilities.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
void size(const label)
Override size to be inconsistent with allocated storage.
void setSize(const label)
Reset size of List.
A list that is sorted upon construction or when explicitly requested with the sort() method.
A List with indirect addressing.
static bool & parRun()
Is this a parallel run?
label size() const
Return the number of elements in the UPtrList.
static void addNote(const string &)
Add extra notes for the usage information.
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
bool optionFound(const word &opt) const
Return true if the named option is found.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
void reset(T *=nullptr)
If object pointer already set, delete object and set to given.
A collection of cell labels.
A cell is defined as a list of faces with extra functionality.
static autoPtr< decompositionMethod > NewDecomposer(const dictionary &decompositionDict)
Return a reference to the selected decomposition method.
A list of keyword definitions, which are a keyword followed by any number of values (e....
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T, if not found return the given default.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
const word & name() const
Return const reference to name.
A class for handling file names.
Post-processing mesh subset tool. Given the original mesh and the list of selected cells,...
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return the top-level database.
virtual bool write(const bool write=true) const
Write mesh using IO settings from time.
virtual void topoChange(const polyTopoChangeMap &map)
Update mesh corresponding to the given map.
label nTotalCells() const
Return total number of cells in decomposed mesh.
Mesh consisting of general polyhedral cells.
const fileName & facesInstance() const
Return the current instance directory for faces.
virtual const faceList & faces() const
Return raw faces.
void resetPrimitives(pointField &&points, faceList &&faces, labelList &&owner, labelList &&neighbour, const labelList &patchSizes, const labelList &patchStarts, const bool validBoundary=true)
Reset mesh primitive data. Assumes all patch info correct.
virtual const labelList & faceOwner() const
Return face owner.
const globalMeshData & globalData() const
Return parallel info.
const fileName & pointsInstance() const
Return the current instance directory for points.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const labelList & faceNeighbour() const
Return face neighbour.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
void setInstance(const fileName &)
Set the instance for mesh files.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Direct mesh changes based on v1.3 polyTopoChange syntax.
Cell-face mesh analysis engine.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
label nInternal1Edges() const
Internal edges using 0 or 1 boundary point.
const vectorField & cellCentres() const
virtual const labelList & faceOwner() const =0
Face face-owner addressing.
label nInternalFaces() const
label nInternal0Edges() const
Internal edges (i.e. not on boundary face) using.
label nInternalEdges() const
Internal edges using 0,1 or 2 boundary points.
virtual const labelList & faceNeighbour() const =0
Face face-neighbour addressing.
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
label nInternalPoints() const
Points not on boundary.
const cellList & cells() const
virtual bool write(const bool write=true) const
Write using setting from DB.
Abstract base class for renumbering.
static autoPtr< renumberMethod > New(const dictionary &renumberDict)
Return a reference to the selected renumbering method.
virtual labelList renumber(const pointField &) const
Return the order in which cells need to be visited, i.e.
static void addOptions(const bool constant=true, const bool withZero=false)
Add the options handled by timeSelector to argList::validOptions.
static instantList select0(Time &runTime, const argList &args)
Return the set of times selected based on the argList options.
A class for handling words, derived from string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
int main(int argc, char *argv[])
gmvFile<< "tracers "<< particles.size()<< nl;{ pointField positions(particles.size());label particlei=0;forAllConstIter(Cloud< passiveParticle >, particles, iter) { positions[particlei++]=iter().position(mesh);} for(i=0;i< pTraits< point >::nComponents;i++) { forAll(positions, particlei) { gmvFile<< component(positions[particlei], i)<< ' ';} gmvFile<< nl;}}forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const fvPatchList & patches
Info<< "Calculating turbulent flame speed field St\n"<< endl;volScalarField St(IOobject("St", runTime.name(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), flameWrinkling->Xi() *Su);multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
errorManipArg< error, int > exit(error &err, const int errNo=1)
List< label > labelList
A List of labels.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
bool rm(const fileName &)
Remove a file, returning true if successful otherwise false.
IOdictionary systemDict(const word &dictName, const argList &args, const objectRegistry &ob, const word ®ionName=polyMesh::defaultRegion)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
word name(const bool)
Return a word representation of a bool.
labelList invert(const label len, const labelUList &)
Invert one-to-one map. Unmapped elements will be -1.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
int order(const scalar s)
errorManip< error > abort(error &err)
const dimensionSet dimless
ListType renumber(const labelUList &oldToNew, const ListType &)
Renumber the values (not the indices) of a list.
labelListList invertOneToMany(const label len, const labelUList &)
Invert one-to-many map. Unmapped elements will be size 0.
List< scalar > scalarList
A List of scalars.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
dimensionedScalar sqrt(const dimensionedScalar &ds)
ListType reorder(const label size, const typename ListType::value_type &defaultValue, const labelUList &oldToNew, const ListType &lst)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
void inplaceRenumber(const labelUList &oldToNew, ListType &)
Inplace renumber the values of a list.
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
IOList< label > labelIOList
Label container classes.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
static const label labelMax
labelList identityMap(const label len)
Create identity map (map[i] == i) of given length.
void sortedOrder(const UList< T > &, labelList &order)
Generate the (stable) sort order for the list.
Foam::argList args(argc, argv)