renumberMesh: Now supports sets

Resolves bug-report http://www.openfoam.org/mantisbt/view.php?id=1377
This commit is contained in:
Henry Weller 2015-10-21 16:49:57 +01:00
parent e64a8929a6
commit 135414cc16
2 changed files with 195 additions and 73 deletions

View file

@ -3,7 +3,7 @@
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\/ M anispulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
This file is part of OpenFOAM. This file is part of OpenFOAM.
@ -46,11 +46,15 @@ Description
#include "zeroGradientFvPatchFields.H" #include "zeroGradientFvPatchFields.H"
#include "CuthillMcKeeRenumber.H" #include "CuthillMcKeeRenumber.H"
#include "fvMeshSubset.H" #include "fvMeshSubset.H"
#include "cellSet.H"
#include "faceSet.H"
#include "pointSet.H"
#ifdef FOAM_USE_ZOLTAN #ifdef FOAM_USE_ZOLTAN
#include "zoltanRenumber.H" #include "zoltanRenumber.H"
#endif #endif
using namespace Foam; using namespace Foam;
@ -163,7 +167,7 @@ void getBand
labelList getFaceOrder labelList getFaceOrder
( (
const primitiveMesh& mesh, const primitiveMesh& mesh,
const labelList& cellOrder // new to old cell const labelList& cellOrder // New to old cell
) )
{ {
labelList reverseCellOrder(invert(cellOrder.size(), cellOrder)); labelList reverseCellOrder(invert(cellOrder.size(), cellOrder));
@ -244,8 +248,7 @@ labelList getFaceOrder
( (
"getFaceOrder" "getFaceOrder"
"(const primitiveMesh&, const labelList&, const labelList&)" "(const primitiveMesh&, const labelList&, const labelList&)"
) << "Did not determine new position" ) << "Did not determine new position" << " for face " << faceI
<< " for face " << faceI
<< abort(FatalError); << abort(FatalError);
} }
} }
@ -259,12 +262,10 @@ labelList getFaceOrder
labelList getRegionFaceOrder labelList getRegionFaceOrder
( (
const primitiveMesh& mesh, const primitiveMesh& mesh,
const labelList& cellOrder, // new to old cell const labelList& cellOrder, // New to old cell
const labelList& cellToRegion // old cell to region const labelList& cellToRegion // Old cell to region
) )
{ {
//Pout<< "Determining face order:" << endl;
labelList reverseCellOrder(invert(cellOrder.size(), cellOrder)); labelList reverseCellOrder(invert(cellOrder.size(), cellOrder));
labelList oldToNewFace(mesh.nFaces(), -1); labelList oldToNewFace(mesh.nFaces(), -1);
@ -280,8 +281,6 @@ labelList getRegionFaceOrder
if (cellToRegion[oldCellI] != prevRegion) if (cellToRegion[oldCellI] != prevRegion)
{ {
prevRegion = cellToRegion[oldCellI]; prevRegion = cellToRegion[oldCellI];
//Pout<< " region " << prevRegion << " internal faces start at "
// << newFaceI << endl;
} }
const cell& cFaces = mesh.cells()[oldCellI]; const cell& cFaces = mesh.cells()[oldCellI];
@ -368,9 +367,6 @@ labelList getRegionFaceOrder
if (prevKey != key) if (prevKey != key)
{ {
//Pout<< " faces inbetween region " << key/nRegions
// << " and " << key%nRegions
// << " start at " << newFaceI << endl;
prevKey = key; prevKey = key;
} }
@ -399,7 +395,6 @@ labelList getRegionFaceOrder
<< abort(FatalError); << abort(FatalError);
} }
} }
//Pout<< endl;
return invert(mesh.nFaces(), oldToNewFace); return invert(mesh.nFaces(), oldToNewFace);
} }
@ -407,7 +402,7 @@ labelList getRegionFaceOrder
// cellOrder: old cell for every new cell // cellOrder: old cell for every new cell
// faceOrder: old face for every new face. Ordering of boundary faces not // faceOrder: old face for every new face. Ordering of boundary faces not
// changed. // changed.
autoPtr<mapPolyMesh> reorderMesh autoPtr<mapPolyMesh> reorderMesh
( (
polyMesh& mesh, polyMesh& mesh,
@ -528,7 +523,7 @@ autoPtr<mapPolyMesh> reorderMesh
( (
new mapPolyMesh new mapPolyMesh
( (
mesh, //const polyMesh& mesh, mesh, // const polyMesh& mesh,
mesh.nPoints(), // nOldPoints, mesh.nPoints(), // nOldPoints,
mesh.nFaces(), // nOldFaces, mesh.nFaces(), // nOldFaces,
mesh.nCells(), // nOldCells, mesh.nCells(), // nOldCells,
@ -634,17 +629,20 @@ int main(int argc, char *argv[])
"calculate the rms of the frontwidth" "calculate the rms of the frontwidth"
); );
// Force linker to include zoltan symbols. This section is only needed since
// Zoltan is a static library
#ifdef FOAM_USE_ZOLTAN
Info<< "renumberMesh built with zoltan support." << nl << endl;
(void)zoltanRenumber::typeName;
#endif
#include "setRootCase.H" #include "setRootCase.H"
#include "createTime.H" #include "createTime.H"
runTime.functionObjects().off(); runTime.functionObjects().off();
// Force linker to include zoltan symbols. This section is only needed since
// Zoltan is a static library
#ifdef FOAM_USE_ZOLTAN
Info<< "renumberMesh built with zoltan support." << nl << endl;
(void)zoltanRenumber::typeName;
#endif
// Get times list // Get times list
instantList Times = runTime.times(); instantList Times = runTime.times();
@ -682,24 +680,26 @@ int main(int argc, char *argv[])
( (
sumSqrIntersect, sumSqrIntersect,
sumOp<scalar>() sumOp<scalar>()
) )/mesh.globalData().nTotalCells()
/mesh.globalData().nTotalCells()
); );
Info<< "Mesh size: " << mesh.globalData().nTotalCells() << nl Info<< "Mesh size: " << mesh.globalData().nTotalCells() << nl
<< "Before renumbering :" << nl << "Before renumbering :" << nl
<< " band : " << band << nl << " band : " << band << nl
<< " profile : " << profile << nl; << " profile : " << profile << nl;
if (doFrontWidth) if (doFrontWidth)
{ {
Info<< " rms frontwidth : " << rmsFrontwidth << nl; Info<< " rms frontwidth : " << rmsFrontwidth << nl;
} }
Info<< endl; Info<< endl;
bool sortCoupledFaceCells = false; bool sortCoupledFaceCells = false;
bool writeMaps = false; bool writeMaps = false;
bool orderPoints = false; bool orderPoints = false;
label blockSize = 0; label blockSize = 0;
bool renumberSets = true;
// Construct renumberMethod // Construct renumberMethod
autoPtr<IOdictionary> renumberDictPtr; autoPtr<IOdictionary> renumberDictPtr;
@ -717,7 +717,6 @@ int main(int argc, char *argv[])
renumberPtr = renumberMethod::New(renumberDict); renumberPtr = renumberMethod::New(renumberDict);
sortCoupledFaceCells = renumberDict.lookupOrDefault sortCoupledFaceCells = renumberDict.lookupOrDefault
( (
"sortCoupledFaceCells", "sortCoupledFaceCells",
@ -760,6 +759,8 @@ int main(int argc, char *argv[])
Info<< "Writing renumber maps (new to old) to polyMesh." << nl Info<< "Writing renumber maps (new to old) to polyMesh." << nl
<< endl; << endl;
} }
renumberSets = renumberDict.lookupOrDefault("renumberSets", true);
} }
else else
{ {
@ -827,6 +828,7 @@ int main(int argc, char *argv[])
// Read objects in time directory // Read objects in time directory
IOobjectList objects(mesh, runTime.timeName()); IOobjectList objects(mesh, runTime.timeName());
// Read vol fields. // Read vol fields.
PtrList<volScalarField> vsFlds; PtrList<volScalarField> vsFlds;
@ -844,6 +846,7 @@ int main(int argc, char *argv[])
PtrList<volTensorField> vtFlds; PtrList<volTensorField> vtFlds;
ReadFields(mesh, objects, vtFlds); ReadFields(mesh, objects, vtFlds);
// Read surface fields. // Read surface fields.
PtrList<surfaceScalarField> ssFlds; PtrList<surfaceScalarField> ssFlds;
@ -861,6 +864,25 @@ int main(int argc, char *argv[])
PtrList<surfaceTensorField> stFlds; PtrList<surfaceTensorField> stFlds;
ReadFields(mesh, objects, stFlds); ReadFields(mesh, objects, stFlds);
// Read point fields.
PtrList<pointScalarField> psFlds;
ReadFields(pointMesh::New(mesh), objects, psFlds);
PtrList<pointVectorField> pvFlds;
ReadFields(pointMesh::New(mesh), objects, pvFlds);
PtrList<pointSphericalTensorField> pstFlds;
ReadFields(pointMesh::New(mesh), objects, pstFlds);
PtrList<pointSymmTensorField> psymtFlds;
ReadFields(pointMesh::New(mesh), objects, psymtFlds);
PtrList<pointTensorField> ptFlds;
ReadFields(pointMesh::New(mesh), objects, ptFlds);
Info<< endl; Info<< endl;
// From renumbering: // From renumbering:
@ -873,7 +895,7 @@ int main(int argc, char *argv[])
// Renumbering in two phases. Should be done in one so mapping of // Renumbering in two phases. Should be done in one so mapping of
// fields is done correctly! // fields is done correctly!
label nBlocks = mesh.nCells() / blockSize; label nBlocks = mesh.nCells()/blockSize;
Info<< "nBlocks = " << nBlocks << endl; Info<< "nBlocks = " << nBlocks << endl;
// Read decompositionMethod dictionary // Read decompositionMethod dictionary
@ -1011,7 +1033,7 @@ int main(int argc, char *argv[])
faceOrder = getFaceOrder faceOrder = getFaceOrder
( (
mesh, mesh,
cellOrder // new to old cell cellOrder // New to old cell
); );
} }
@ -1032,17 +1054,19 @@ int main(int argc, char *argv[])
autoPtr<mapPolyMesh> pointOrderMap = meshMod.changeMesh autoPtr<mapPolyMesh> pointOrderMap = meshMod.changeMesh
( (
mesh, mesh,
false, //inflate false, // inflate
true, //syncParallel true, // syncParallel
false, //orderCells false, // orderCells
orderPoints //orderPoints orderPoints // orderPoints
); );
// Combine point reordering into map. // Combine point reordering into map.
const_cast<labelList&>(map().pointMap()) = UIndirectList<label> const_cast<labelList&>(map().pointMap()) = UIndirectList<label>
( (
map().pointMap(), map().pointMap(),
pointOrderMap().pointMap() pointOrderMap().pointMap()
)(); )();
inplaceRenumber inplaceRenumber
( (
pointOrderMap().reversePointMap(), pointOrderMap().reversePointMap(),
@ -1055,7 +1079,6 @@ int main(int argc, char *argv[])
mesh.updateMesh(map); mesh.updateMesh(map);
// Update proc maps // Update proc maps
if (cellProcAddressing.headerOk())
if if
( (
cellProcAddressing.headerOk() cellProcAddressing.headerOk()
@ -1070,7 +1093,6 @@ int main(int argc, char *argv[])
UIndirectList<label>(cellProcAddressing, map().cellMap()) UIndirectList<label>(cellProcAddressing, map().cellMap())
); );
} }
if (faceProcAddressing.headerOk())
if if
( (
faceProcAddressing.headerOk() faceProcAddressing.headerOk()
@ -1101,7 +1123,6 @@ int main(int argc, char *argv[])
} }
} }
} }
if (pointProcAddressing.headerOk())
if if
( (
pointProcAddressing.headerOk() pointProcAddressing.headerOk()
@ -1147,17 +1168,19 @@ int main(int argc, char *argv[])
( (
sumSqrIntersect, sumSqrIntersect,
sumOp<scalar>() sumOp<scalar>()
) )/mesh.globalData().nTotalCells()
/ mesh.globalData().nTotalCells()
); );
Info<< "After renumbering :" << nl Info<< "After renumbering :" << nl
<< " band : " << band << nl << " band : " << band << nl
<< " profile : " << profile << nl; << " profile : " << profile << nl;
if (doFrontWidth) if (doFrontWidth)
{ {
Info<< " rms frontwidth : " << rmsFrontwidth << nl; Info<< " rms frontwidth : " << rmsFrontwidth << nl;
} }
Info<< endl; Info<< endl;
} }
@ -1225,47 +1248,83 @@ int main(int argc, char *argv[])
Info<< "Writing mesh to " << mesh.facesInstance() << endl; Info<< "Writing mesh to " << mesh.facesInstance() << endl;
mesh.write(); mesh.write();
if (cellProcAddressing.headerOk()) if (cellProcAddressing.headerOk())
if
(
cellProcAddressing.headerOk()
&& cellProcAddressing.size() == mesh.nCells()
)
{ {
cellProcAddressing.instance() = mesh.facesInstance(); cellProcAddressing.instance() = mesh.facesInstance();
cellProcAddressing.write(); if (cellProcAddressing.size() == mesh.nCells())
} {
if (faceProcAddressing.headerOk()) cellProcAddressing.write();
if }
( else
faceProcAddressing.headerOk() {
&& faceProcAddressing.size() == mesh.nFaces() // procAddressing file no longer valid. Delete it.
) const fileName fName(cellProcAddressing.filePath());
{ if (fName.size())
faceProcAddressing.instance() = mesh.facesInstance(); {
faceProcAddressing.write(); Info<< "Deleting inconsistent processor cell decomposition"
} << " map " << fName << endl;
if (pointProcAddressing.headerOk()) rm(fName);
if }
( }
pointProcAddressing.headerOk()
&& pointProcAddressing.size() == mesh.nPoints()
)
{
pointProcAddressing.instance() = mesh.facesInstance();
pointProcAddressing.write();
}
if (boundaryProcAddressing.headerOk())
if
(
boundaryProcAddressing.headerOk()
&& boundaryProcAddressing.size() == mesh.boundaryMesh().size()
)
{
boundaryProcAddressing.instance() = mesh.facesInstance();
boundaryProcAddressing.write();
} }
if (faceProcAddressing.headerOk())
{
faceProcAddressing.instance() = mesh.facesInstance();
if (faceProcAddressing.size() == mesh.nFaces())
{
faceProcAddressing.write();
}
else
{
const fileName fName(faceProcAddressing.filePath());
if (fName.size())
{
Info<< "Deleting inconsistent processor face decomposition"
<< " map " << fName << endl;
rm(fName);
}
}
}
if (pointProcAddressing.headerOk())
{
pointProcAddressing.instance() = mesh.facesInstance();
if (pointProcAddressing.size() == mesh.nPoints())
{
pointProcAddressing.write();
}
else
{
const fileName fName(pointProcAddressing.filePath());
if (fName.size())
{
Info<< "Deleting inconsistent processor point decomposition"
<< " map " << fName << endl;
rm(fName);
}
}
}
if (boundaryProcAddressing.headerOk())
{
boundaryProcAddressing.instance() = mesh.facesInstance();
if (boundaryProcAddressing.size() == mesh.boundaryMesh().size())
{
boundaryProcAddressing.write();
}
else
{
const fileName fName(boundaryProcAddressing.filePath());
if (fName.size())
{
Info<< "Deleting inconsistent processor patch decomposition"
<< " map " << fName << endl;
rm(fName);
}
}
}
if (writeMaps) if (writeMaps)
{ {
@ -1276,17 +1335,18 @@ int main(int argc, char *argv[])
"origCellID", "origCellID",
map().cellMap() map().cellMap()
)().write(); )().write();
createScalarField createScalarField
( (
mesh, mesh,
"cellID", "cellID",
identity(mesh.nCells()) identity(mesh.nCells())
)().write(); )().write();
Info<< nl << "Written current cellID and origCellID as volScalarField" Info<< nl << "Written current cellID and origCellID as volScalarField"
<< " for use in postprocessing." << " for use in postprocessing."
<< nl << endl; << nl << endl;
labelIOList labelIOList
( (
IOobject IOobject
@ -1301,6 +1361,7 @@ int main(int argc, char *argv[])
), ),
map().cellMap() map().cellMap()
).write(); ).write();
labelIOList labelIOList
( (
IOobject IOobject
@ -1315,6 +1376,7 @@ int main(int argc, char *argv[])
), ),
map().faceMap() map().faceMap()
).write(); ).write();
labelIOList labelIOList
( (
IOobject IOobject
@ -1331,6 +1393,64 @@ int main(int argc, char *argv[])
).write(); ).write();
} }
// Renumber sets if required
if (renumberSets)
{
Info<< endl;
// Read sets
IOobjectList objects(mesh, mesh.facesInstance(), "polyMesh/sets");
{
IOobjectList cSets(objects.lookupClass(cellSet::typeName));
if (cSets.size())
{
Info<< "Renumbering cellSets:" << endl;
forAllConstIter(IOobjectList, cSets, iter)
{
cellSet cs(*iter());
Info<< " " << cs.name() << endl;
cs.updateMesh(map());
cs.instance() = mesh.facesInstance();
cs.write();
}
}
}
{
IOobjectList fSets(objects.lookupClass(faceSet::typeName));
if (fSets.size())
{
Info<< "Renumbering faceSets:" << endl;
forAllConstIter(IOobjectList, fSets, iter)
{
faceSet fs(*iter());
Info<< " " << fs.name() << endl;
fs.updateMesh(map());
fs.instance() = mesh.facesInstance();
fs.write();
}
}
}
{
IOobjectList pSets(objects.lookupClass(pointSet::typeName));
if (pSets.size())
{
Info<< "Renumbering pointSets:" << endl;
forAllConstIter(IOobjectList, pSets, iter)
{
pointSet ps(*iter());
Info<< " " << ps.name() << endl;
ps.updateMesh(map());
ps.instance() = mesh.facesInstance();
ps.write();
}
}
}
}
Info<< "\nEnd\n" << endl; Info<< "\nEnd\n" << endl;
return 0; return 0;

View file

@ -34,6 +34,8 @@ sortCoupledFaceCells false;
// Optional entry: sort points into internal and boundary points // Optional entry: sort points into internal and boundary points
//orderPoints false; //orderPoints false;
// Optional: suppress renumbering cellSets,faceSets,pointSets
//renumberSets false;
method CuthillMcKee; method CuthillMcKee;