meshCheck: New library for mesh checking functions

Used by the checkMesh utility and functionObject
This commit is contained in:
Henry Weller
2023-10-19 15:00:28 +01:00
parent 738f205ac8
commit 715dcd97d0
31 changed files with 107 additions and 84 deletions

View File

@ -1,6 +1,3 @@
checkTools.C
checkTopology.C
checkGeometry.C
checkMeshQuality.C
checkMesh.C

View File

@ -1,14 +1,14 @@
EXE_INC = \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/meshCheck/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/surfMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude
EXE_LIBS = \
-lmeshTools \
-lgenericPatches \
-lmeshCheck \
-lfiniteVolume \
-lsampling \
-lsurfMesh \
-ldynamicMesh
-ldynamicMesh \
-lgenericPatches

View File

@ -1,78 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Description
Routines for checking mesh geometry
SourceFiles
checkGeometry.C
\*---------------------------------------------------------------------------*/
#ifndef checkGeometry_H
#define checkGeometry_H
#include "label.H"
#include "HashSet.H"
#include "labelVector.H"
#include "setWriter.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class polyMesh;
class wedgePolyPatch;
class surfaceWriter;
//- Find wedge with opposite orientation. Note: does not actually check
// that it is opposite, only that it has opposite normal and same axis.
label findOppositeWedge(const polyMesh&, const wedgePolyPatch&);
//- Check wedge orientation
bool checkWedges
(
const polyMesh&,
const bool report,
const Vector<label>&,
labelHashSet*
);
//- Check 0th vertex on coupled faces
bool checkCoupledPoints(const polyMesh&, const bool report, labelHashSet*);
//- Check the geometry
label checkGeometry
(
const polyMesh& mesh,
const bool allGeometry,
const autoPtr<surfaceWriter>&,
const autoPtr<setWriter>&
);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -196,13 +196,13 @@ int main(int argc, char *argv[])
// Reconstruct globalMeshData
mesh.globalData();
printMeshStats(mesh, allTopology);
meshTools::printMeshStats(mesh, allTopology);
label nFailedChecks = 0;
if (!noTopology)
{
nFailedChecks += checkTopology
nFailedChecks += meshCheck::checkTopology
(
mesh,
allTopology,
@ -211,7 +211,7 @@ int main(int argc, char *argv[])
);
}
nFailedChecks += checkGeometry
nFailedChecks += meshCheck::checkGeometry
(
mesh,
allGeometry,
@ -242,7 +242,7 @@ int main(int argc, char *argv[])
{
Info<< "Time = " << runTime.userTimeName() << nl << endl;
label nFailedChecks = checkGeometry
label nFailedChecks = meshCheck::checkGeometry
(
mesh,
allGeometry,

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -58,7 +58,7 @@ Foam::label Foam::checkMeshQuality
faces.write();
if (writer.valid())
{
mergeAndWrite(writer(), faces);
meshTools::mergeAndWrite(writer(), faces);
}
}
}

View File

@ -1,467 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2015-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "checkTools.H"
#include "polyMesh.H"
#include "globalMeshData.H"
#include "hexMatcher.H"
#include "wedgeMatcher.H"
#include "prismMatcher.H"
#include "pyrMatcher.H"
#include "tetWedgeMatcher.H"
#include "tetMatcher.H"
#include "IOmanip.H"
#include "pointSet.H"
#include "faceSet.H"
#include "cellSet.H"
#include "Time.H"
#include "surfaceWriter.H"
#include "syncTools.H"
#include "globalIndex.H"
#include "PatchTools.H"
#include "writeFile.H"
#include "coordSet.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
void Foam::printMeshStats(const polyMesh& mesh, const bool allTopology)
{
Info<< "Mesh stats" << nl
<< " points: "
<< returnReduce(mesh.points().size(), sumOp<label>()) << nl;
label nInternalPoints = returnReduce
(
mesh.nInternalPoints(),
sumOp<label>()
);
if (nInternalPoints != -Pstream::nProcs())
{
Info<< " internal points: " << nInternalPoints << nl;
if (returnReduce(mesh.nInternalPoints(), minOp<label>()) == -1)
{
WarningInFunction
<< "Some processors have their points sorted into internal"
<< " and external and some do not." << endl
<< "This can cause problems later on." << endl;
}
}
if (allTopology && nInternalPoints != -Pstream::nProcs())
{
label nEdges = returnReduce(mesh.nEdges(), sumOp<label>());
label nInternalEdges = returnReduce
(
mesh.nInternalEdges(),
sumOp<label>()
);
label nInternal1Edges = returnReduce
(
mesh.nInternal1Edges(),
sumOp<label>()
);
label nInternal0Edges = returnReduce
(
mesh.nInternal0Edges(),
sumOp<label>()
);
Info<< " edges: " << nEdges << nl
<< " internal edges: " << nInternalEdges << nl
<< " internal edges using one boundary point: "
<< nInternal1Edges-nInternal0Edges << nl
<< " internal edges using two boundary points: "
<< nInternalEdges-nInternal1Edges << nl;
}
label nFaces = returnReduce(mesh.faces().size(), sumOp<label>());
label nIntFaces = returnReduce(mesh.faceNeighbour().size(), sumOp<label>());
label nCells = returnReduce(mesh.cells().size(), sumOp<label>());
Info<< " faces: " << nFaces << nl
<< " internal faces: " << nIntFaces << nl
<< " cells: " << nCells << nl
<< " faces per cell: "
<< scalar(nFaces + nIntFaces)/max(1, nCells) << nl
<< " boundary patches: " << mesh.boundaryMesh().size() << nl
<< " point zones: " << mesh.pointZones().size() << nl
<< " face zones: " << mesh.faceZones().size() << nl
<< " cell zones: " << mesh.cellZones().size() << nl
<< endl;
// Construct shape recognisers
hexMatcher hex;
prismMatcher prism;
wedgeMatcher wedge;
pyrMatcher pyr;
tetWedgeMatcher tetWedge;
tetMatcher tet;
// Counters for different cell types
label nHex = 0;
label nWedge = 0;
label nPrism = 0;
label nPyr = 0;
label nTet = 0;
label nTetWedge = 0;
label nUnknown = 0;
Map<label> polyhedralFaces;
for (label celli = 0; celli < mesh.nCells(); celli++)
{
if (hex.isA(mesh, celli))
{
nHex++;
}
else if (tet.isA(mesh, celli))
{
nTet++;
}
else if (pyr.isA(mesh, celli))
{
nPyr++;
}
else if (prism.isA(mesh, celli))
{
nPrism++;
}
else if (wedge.isA(mesh, celli))
{
nWedge++;
}
else if (tetWedge.isA(mesh, celli))
{
nTetWedge++;
}
else
{
nUnknown++;
polyhedralFaces(mesh.cells()[celli].size())++;
}
}
reduce(nHex,sumOp<label>());
reduce(nPrism,sumOp<label>());
reduce(nWedge,sumOp<label>());
reduce(nPyr,sumOp<label>());
reduce(nTetWedge,sumOp<label>());
reduce(nTet,sumOp<label>());
reduce(nUnknown,sumOp<label>());
Info<< "Overall number of cells of each type:" << nl
<< " hexahedra: " << nHex << nl
<< " prisms: " << nPrism << nl
<< " wedges: " << nWedge << nl
<< " pyramids: " << nPyr << nl
<< " tet wedges: " << nTetWedge << nl
<< " tetrahedra: " << nTet << nl
<< " polyhedra: " << nUnknown
<< endl;
if (nUnknown > 0)
{
Pstream::mapCombineGather(polyhedralFaces, plusEqOp<label>());
Info<< " Breakdown of polyhedra by number of faces:" << nl
<< " faces" << " number of cells" << endl;
const labelList sortedKeys = polyhedralFaces.sortedToc();
forAll(sortedKeys, keyI)
{
const label nFaces = sortedKeys[keyI];
Info<< setf(std::ios::right) << setw(13)
<< nFaces << " " << polyhedralFaces[nFaces] << nl;
}
}
Info<< endl;
}
void Foam::mergeAndWrite
(
const polyMesh& mesh,
const surfaceWriter& writer,
const word& name,
const indirectPrimitivePatch setPatch,
const fileName& outputDir
)
{
if (Pstream::parRun())
{
labelList pointToGlobal;
labelList uniqueMeshPointLabels;
autoPtr<globalIndex> globalPoints;
autoPtr<globalIndex> globalFaces;
faceList mergedFaces;
pointField mergedPoints;
Foam::PatchTools::gatherAndMerge
(
mesh,
setPatch.localFaces(),
setPatch.meshPoints(),
setPatch.meshPointMap(),
pointToGlobal,
uniqueMeshPointLabels,
globalPoints,
globalFaces,
mergedFaces,
mergedPoints
);
// Write
if (Pstream::master())
{
writer.write(outputDir, name, mergedPoints, mergedFaces);
}
}
else
{
writer.write
(
outputDir,
name,
setPatch.localPoints(),
setPatch.localFaces()
);
}
}
Foam::fileName Foam::checkMeshOutputDir(const polyMesh& mesh)
{
return
mesh.time().globalPath()
/functionObjects::writeFile::outputPrefix
/(
(mesh.name() != polyMesh::defaultRegion)
? mesh.name()
: word::null
)
/"checkMesh"
/mesh.pointsInstance();
}
void Foam::mergeAndWrite
(
const surfaceWriter& writer,
const faceSet& set
)
{
const polyMesh& mesh = refCast<const polyMesh>(set.db());
const indirectPrimitivePatch setPatch
(
IndirectList<face>(mesh.faces(), set.sortedToc()),
mesh.points()
);
mergeAndWrite(mesh, writer, set.name(), setPatch, checkMeshOutputDir(mesh));
}
void Foam::mergeAndWrite
(
const surfaceWriter& writer,
const cellSet& set
)
{
const polyMesh& mesh = refCast<const polyMesh>(set.db());
const polyBoundaryMesh& pbm = mesh.boundaryMesh();
// Determine faces on outside of cellSet
PackedBoolList isInSet(mesh.nCells());
forAllConstIter(cellSet, set, iter)
{
isInSet[iter.key()] = true;
}
boolList bndInSet(mesh.nFaces()-mesh.nInternalFaces());
forAll(pbm, patchi)
{
const polyPatch& pp = pbm[patchi];
const labelList& fc = pp.faceCells();
forAll(fc, i)
{
bndInSet[pp.start()+i-mesh.nInternalFaces()] = isInSet[fc[i]];
}
}
syncTools::swapBoundaryFaceList(mesh, bndInSet);
DynamicList<label> outsideFaces(3*set.size());
for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
{
bool ownVal = isInSet[mesh.faceOwner()[facei]];
bool neiVal = isInSet[mesh.faceNeighbour()[facei]];
if (ownVal != neiVal)
{
outsideFaces.append(facei);
}
}
forAll(pbm, patchi)
{
const polyPatch& pp = pbm[patchi];
const labelList& fc = pp.faceCells();
if (pp.coupled())
{
forAll(fc, i)
{
label facei = pp.start()+i;
bool neiVal = bndInSet[facei-mesh.nInternalFaces()];
if (isInSet[fc[i]] && !neiVal)
{
outsideFaces.append(facei);
}
}
}
else
{
forAll(fc, i)
{
if (isInSet[fc[i]])
{
outsideFaces.append(pp.start()+i);
}
}
}
}
const indirectPrimitivePatch setPatch
(
IndirectList<face>(mesh.faces(), outsideFaces),
mesh.points()
);
const fileName outputDir(checkMeshOutputDir(mesh));
outputDir.clean();
mergeAndWrite(mesh, writer, set.name(), setPatch, outputDir);
}
void Foam::mergeAndWrite
(
const setWriter& writer,
const pointSet& set
)
{
const polyMesh& mesh = refCast<const polyMesh>(set.db());
pointField mergedPts;
labelList mergedIDs;
if (Pstream::parRun())
{
// Note: we explicitly do not merge the points
// (mesh.globalData().mergePoints etc) since this might
// hide any synchronisation problem
globalIndex globalNumbering(mesh.nPoints());
mergedPts.setSize(returnReduce(set.size(), sumOp<label>()));
mergedIDs.setSize(mergedPts.size());
labelList setPointIDs(set.sortedToc());
// Get renumbered local data
pointField myPoints(mesh.points(), setPointIDs);
labelList myIDs(setPointIDs.size());
forAll(setPointIDs, i)
{
myIDs[i] = globalNumbering.toGlobal(setPointIDs[i]);
}
if (Pstream::master())
{
// Insert master data first
label pOffset = 0;
SubList<point>(mergedPts, myPoints.size(), pOffset) = myPoints;
SubList<label>(mergedIDs, myIDs.size(), pOffset) = myIDs;
pOffset += myPoints.size();
// Receive slave ones
for (int slave=1; slave<Pstream::nProcs(); slave++)
{
IPstream fromSlave(Pstream::commsTypes::scheduled, slave);
pointField slavePts(fromSlave);
labelList slaveIDs(fromSlave);
SubList<point>(mergedPts, slavePts.size(), pOffset) = slavePts;
SubList<label>(mergedIDs, slaveIDs.size(), pOffset) = slaveIDs;
pOffset += slaveIDs.size();
}
}
else
{
// Construct processor stream with estimate of size. Could
// be improved.
OPstream toMaster
(
Pstream::commsTypes::scheduled,
Pstream::masterNo(),
myPoints.byteSize() + myIDs.byteSize()
);
toMaster << myPoints << myIDs;
}
}
else
{
mergedIDs = set.sortedToc();
mergedPts = pointField(mesh.points(), mergedIDs);
}
// Write with scalar pointID
if (Pstream::master())
{
writer.write
(
checkMeshOutputDir(mesh),
set.name(),
coordSet(false, word::null, mergedPts),
"pointID",
scalarField(scalarList(mergedIDs))
);
}
}
// ************************************************************************* //

View File

@ -1,85 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2015-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Description
Tools for checking the mesh
SourceFiles
checkTools.C
\*---------------------------------------------------------------------------*/
#ifndef checkTools_H
#define checkTools_H
#include "scalar.H"
#include "indirectPrimitivePatch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class polyMesh;
class surfaceWriter;
class setWriter;
class pointSet;
class faceSet;
class cellSet;
class fileName;
class polyMesh;
//- Output directory for sets and surfaces
fileName checkMeshOutputDir(const polyMesh& mesh);
//- Print mesh statistics
void printMeshStats(const polyMesh& mesh, const bool allTopology);
//- Generate merged surface on master and write. Needs input patch
// to be of mesh faces.
void mergeAndWrite
(
const polyMesh& mesh,
const surfaceWriter& setWriter,
const word& name,
const indirectPrimitivePatch setPatch,
const fileName& outputDir
);
//- Write representation of (assembled) faceSet to surface file in
// postProcessing/ directory
void mergeAndWrite(const surfaceWriter&, const faceSet&);
//- Write representation of (assembled) cellSet to surface file in
// postProcessing/ directory
void mergeAndWrite(const surfaceWriter&, const cellSet&);
//- Write representation of (assembled) pointSet to 'set' file in
// postProcessing/ directory
void mergeAndWrite(const setWriter&, const pointSet&);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -1,594 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "checkTopology.H"
#include "polyMeshCheck.H"
#include "Time.H"
#include "regionSplit.H"
#include "cellSet.H"
#include "faceSet.H"
#include "pointSet.H"
#include "IOmanip.H"
#include "emptyPolyPatch.H"
#include "processorPolyPatch.H"
#include "surfaceWriter.H"
#include "checkTools.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Foam::label Foam::checkTopology
(
const polyMesh& mesh,
const bool allTopology,
const autoPtr<surfaceWriter>& surfWriter,
const autoPtr<Foam::setWriter>& setWriter
)
{
label noFailedChecks = 0;
Info<< "Checking topology..." << endl;
// Check if the boundary definition is unique
mesh.boundaryMesh().checkDefinition(true);
// Check that empty patches cover all sides of the mesh
{
label nEmpty = 0;
forAll(mesh.boundaryMesh(), patchi)
{
if (isA<emptyPolyPatch>(mesh.boundaryMesh()[patchi]))
{
nEmpty += mesh.boundaryMesh()[patchi].size();
}
}
reduce(nEmpty, sumOp<label>());
label nTotCells = returnReduce(mesh.cells().size(), sumOp<label>());
// These are actually warnings, not errors.
if (nTotCells && (nEmpty % nTotCells))
{
Info<< " ***Total number of faces on empty patches"
<< " is not divisible by the number of cells in the mesh."
<< " Hence this mesh is not 1D or 2D."
<< endl;
}
}
// Check if the boundary processor patches are correct
mesh.boundaryMesh().checkParallelSync(true);
// Check names of zones are equal
mesh.cellZones().checkDefinition(true);
if (mesh.cellZones().checkParallelSync(true))
{
noFailedChecks++;
}
mesh.faceZones().checkDefinition(true);
if (mesh.faceZones().checkParallelSync(true))
{
noFailedChecks++;
}
mesh.pointZones().checkDefinition(true);
if (mesh.pointZones().checkParallelSync(true))
{
noFailedChecks++;
}
{
cellSet cells(mesh, "illegalCells", mesh.nCells()/100);
forAll(mesh.cells(), celli)
{
const cell& cFaces = mesh.cells()[celli];
if (cFaces.size() <= 3)
{
cells.insert(celli);
}
forAll(cFaces, i)
{
if (cFaces[i] < 0 || cFaces[i] >= mesh.nFaces())
{
cells.insert(celli);
break;
}
}
}
label nCells = returnReduce(cells.size(), sumOp<label>());
if (nCells > 0)
{
Info<< " Illegal cells (less than 4 faces or out of range faces)"
<< " found, number of cells: " << nCells << endl;
noFailedChecks++;
Info<< " <<Writing " << nCells
<< " illegal cells to set " << cells.name() << endl;
cells.instance() = mesh.pointsInstance();
cells.write();
if (surfWriter.valid())
{
mergeAndWrite(surfWriter(), cells);
}
}
else
{
Info<< " Cell to face addressing OK." << endl;
}
}
{
pointSet points(mesh, "unusedPoints", mesh.nPoints()/100);
if (meshCheck::checkPoints(mesh, true, &points))
{
noFailedChecks++;
label nPoints = returnReduce(points.size(), sumOp<label>());
Info<< " <<Writing " << nPoints
<< " unused points to set " << points.name() << endl;
points.instance() = mesh.pointsInstance();
points.write();
if (setWriter.valid())
{
mergeAndWrite(setWriter, points);
}
}
}
{
faceSet faces(mesh, "upperTriangularFace", mesh.nFaces()/100);
if (meshCheck::checkUpperTriangular(mesh, true, &faces))
{
noFailedChecks++;
}
label nFaces = returnReduce(faces.size(), sumOp<label>());
if (nFaces > 0)
{
Info<< " <<Writing " << nFaces
<< " unordered faces to set " << faces.name() << endl;
faces.instance() = mesh.pointsInstance();
faces.write();
if (surfWriter.valid())
{
mergeAndWrite(surfWriter(), faces);
}
}
}
{
faceSet faces(mesh, "outOfRangeFaces", mesh.nFaces()/100);
if (meshCheck::checkFaceVertices(mesh, true, &faces))
{
noFailedChecks++;
label nFaces = returnReduce(faces.size(), sumOp<label>());
Info<< " <<Writing " << nFaces
<< " faces with out-of-range or duplicate vertices to set "
<< faces.name() << endl;
faces.instance() = mesh.pointsInstance();
faces.write();
if (surfWriter.valid())
{
mergeAndWrite(surfWriter(), faces);
}
}
}
if (allTopology)
{
cellSet cells(mesh, "zipUpCells", mesh.nCells()/100);
if (meshCheck::checkCellsZipUp(mesh, true, &cells))
{
noFailedChecks++;
label nCells = returnReduce(cells.size(), sumOp<label>());
Info<< " <<Writing " << nCells
<< " cells with over used edges to set " << cells.name()
<< endl;
cells.instance() = mesh.pointsInstance();
cells.write();
if (surfWriter.valid())
{
mergeAndWrite(surfWriter(), cells);
}
}
}
if (allTopology)
{
faceSet faces(mesh, "edgeFaces", mesh.nFaces()/100);
if (meshCheck::checkFaceFaces(mesh, true, &faces))
{
noFailedChecks++;
}
label nFaces = returnReduce(faces.size(), sumOp<label>());
if (nFaces > 0)
{
Info<< " <<Writing " << nFaces
<< " faces with non-standard edge connectivity to set "
<< faces.name() << endl;
faces.instance() = mesh.pointsInstance();
faces.write();
if (surfWriter.valid())
{
mergeAndWrite(surfWriter(), faces);
}
}
}
if (allTopology)
{
labelList nInternalFaces(mesh.nCells(), 0);
for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
{
nInternalFaces[mesh.faceOwner()[facei]]++;
nInternalFaces[mesh.faceNeighbour()[facei]]++;
}
const polyBoundaryMesh& patches = mesh.boundaryMesh();
forAll(patches, patchi)
{
if (patches[patchi].coupled())
{
const labelUList& owners = patches[patchi].faceCells();
forAll(owners, i)
{
nInternalFaces[owners[i]]++;
}
}
}
cellSet oneCells(mesh, "oneInternalFaceCells", mesh.nCells()/100);
cellSet twoCells(mesh, "twoInternalFacesCells", mesh.nCells()/100);
forAll(nInternalFaces, celli)
{
if (nInternalFaces[celli] <= 1)
{
oneCells.insert(celli);
}
else if (nInternalFaces[celli] == 2)
{
twoCells.insert(celli);
}
}
label nOneCells = returnReduce(oneCells.size(), sumOp<label>());
if (nOneCells > 0)
{
Info<< " <<Writing " << nOneCells
<< " cells with zero or one non-boundary face to set "
<< oneCells.name()
<< endl;
oneCells.instance() = mesh.pointsInstance();
oneCells.write();
if (surfWriter.valid())
{
mergeAndWrite(surfWriter(), oneCells);
}
}
label nTwoCells = returnReduce(twoCells.size(), sumOp<label>());
if (nTwoCells > 0)
{
Info<< " <<Writing " << nTwoCells
<< " cells with two non-boundary faces to set "
<< twoCells.name()
<< endl;
twoCells.instance() = mesh.pointsInstance();
twoCells.write();
if (surfWriter.valid())
{
mergeAndWrite(surfWriter(), twoCells);
}
}
}
{
regionSplit rs(mesh);
if (rs.nRegions() <= 1)
{
Info<< " Number of regions: " << rs.nRegions() << " (OK)."
<< endl;
}
else
{
Info<< " *Number of regions: " << rs.nRegions() << endl;
Info<< " The mesh has multiple regions which are not connected "
"by any face." << endl
<< " <<Writing region information to "
<< mesh.time().name()/"cellToRegion"
<< endl;
labelIOList ctr
(
IOobject
(
"cellToRegion",
mesh.time().name(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
rs
);
ctr.write();
// Points in multiple regions
pointSet points
(
mesh,
"multiRegionPoints",
mesh.nPoints()/1000
);
// Is region disconnected
boolList regionDisconnected(rs.nRegions(), true);
if (allTopology)
{
// -1 : not assigned
// -2 : multiple regions
// >= 0 : single region
labelList pointToRegion(mesh.nPoints(), -1);
for
(
label faceI = mesh.nInternalFaces();
faceI < mesh.nFaces();
faceI++
)
{
label regionI = rs[mesh.faceOwner()[faceI]];
const face& f = mesh.faces()[faceI];
forAll(f, fp)
{
label& pRegion = pointToRegion[f[fp]];
if (pRegion == -1)
{
pRegion = regionI;
}
else if (pRegion == -2)
{
// Already marked
regionDisconnected[regionI] = false;
}
else if (pRegion != regionI)
{
// Multiple regions
regionDisconnected[regionI] = false;
regionDisconnected[pRegion] = false;
pRegion = -2;
points.insert(f[fp]);
}
}
}
Pstream::listCombineGather(regionDisconnected, andEqOp<bool>());
Pstream::listCombineScatter(regionDisconnected);
}
// write cellSet for each region
PtrList<cellSet> cellRegions(rs.nRegions());
for (label i = 0; i < rs.nRegions(); i++)
{
cellRegions.set
(
i,
new cellSet
(
mesh,
"region" + Foam::name(i),
mesh.nCells()/100
)
);
}
forAll(rs, i)
{
cellRegions[rs[i]].insert(i);
}
for (label i = 0; i < rs.nRegions(); i++)
{
Info<< " <<Writing region " << i;
if (allTopology)
{
if (regionDisconnected[i])
{
Info<< " (fully disconnected)";
}
else
{
Info<< " (point connected)";
}
}
Info<< " with "
<< returnReduce(cellRegions[i].size(), sumOp<scalar>())
<< " cells to cellSet " << cellRegions[i].name() << endl;
cellRegions[i].write();
}
label nPoints = returnReduce(points.size(), sumOp<label>());
if (nPoints)
{
Info<< " <<Writing " << nPoints
<< " points that are in multiple regions to set "
<< points.name() << endl;
points.write();
if (setWriter.valid())
{
mergeAndWrite(setWriter, points);
}
}
}
}
{
if (!Pstream::parRun())
{
Info<< "\nChecking patch topology for multiply connected"
<< " surfaces..." << endl;
}
else
{
Info<< "\nChecking basic patch addressing..." << endl;
}
const polyBoundaryMesh& patches = mesh.boundaryMesh();
// Non-manifold points
pointSet points
(
mesh,
"nonManifoldPoints",
mesh.nPoints()/1000
);
Pout.setf(ios_base::left);
Info<< " "
<< setw(20) << "Patch"
<< setw(9) << "Faces"
<< setw(9) << "Points";
if (!Pstream::parRun())
{
Info<< setw(34) << "Surface topology";
}
Info<< endl;
forAll(patches, patchi)
{
const polyPatch& pp = patches[patchi];
if (!isA<processorPolyPatch>(pp))
{
Info<< " "
<< setw(20) << pp.name()
<< setw(9) << returnReduce(pp.size(), sumOp<label>())
<< setw(9) << returnReduce(pp.nPoints(), sumOp<label>());
if (!Pstream::parRun())
{
primitivePatch::surfaceTopo pTyp = pp.surfaceType();
if (pp.empty())
{
Info<< setw(34) << "ok (empty)";
}
else if (pTyp == primitivePatch::MANIFOLD)
{
if (pp.checkPointManifold(true, &points))
{
Info<< setw(34)
<< "multiply connected (shared point)";
}
else
{
Info<< setw(34) << "ok (closed singly connected)";
}
// Add points on non-manifold edges to make set complete
pp.checkTopology(false, &points);
}
else
{
pp.checkTopology(false, &points);
if (pTyp == primitivePatch::OPEN)
{
Info<< setw(34)
<< "ok (non-closed singly connected)";
}
else
{
Info<< setw(34)
<< "multiply connected (shared edge)";
}
}
}
Info<< endl;
}
}
if (points.size())
{
Info<< " <<Writing " << returnReduce(points.size(), sumOp<label>())
<< " conflicting points to set "
<< points.name() << endl;
points.instance() = mesh.pointsInstance();
points.write();
}
// Info.setf(ios_base::right);
}
// Force creation of all addressing if requested.
// Errors will be reported as required
if (allTopology)
{
mesh.cells();
mesh.faces();
mesh.edges();
mesh.points();
mesh.faceOwner();
mesh.faceNeighbour();
mesh.cellCells();
mesh.edgeCells();
mesh.pointCells();
mesh.edgeFaces();
mesh.pointFaces();
mesh.cellEdges();
mesh.faceEdges();
mesh.pointEdges();
}
return noFailedChecks;
}
// ************************************************************************* //

View File

@ -1,59 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Description
Tools for checking the mesh
SourceFiles
checkTopology.C
\*---------------------------------------------------------------------------*/
#ifndef checkTopology_H
#define checkTopology_H
#include "label.H"
#include "autoPtr.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class polyMesh;
class surfaceWriter;
class setWriter;
label checkTopology
(
const polyMesh& mesh,
const bool allTopology,
const autoPtr<surfaceWriter>& surfWriter,
const autoPtr<Foam::setWriter>& setWriter
);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //