From 134f7abd57ae1c2f6327ba27b2b3d26daaa8162d Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Mon, 8 May 2017 10:50:59 +0100
Subject: [PATCH] ENH: checkMesh: output vol fields of mesh quality. Fixes
#466.
---
.../mesh/manipulation/checkMesh/Make/files | 1 +
.../mesh/manipulation/checkMesh/checkMesh.C | 55 ++-
.../mesh/manipulation/checkMesh/writeFields.C | 360 ++++++++++++++++++
.../mesh/manipulation/checkMesh/writeFields.H | 10 +
4 files changed, 423 insertions(+), 3 deletions(-)
create mode 100644 applications/utilities/mesh/manipulation/checkMesh/writeFields.C
create mode 100644 applications/utilities/mesh/manipulation/checkMesh/writeFields.H
diff --git a/applications/utilities/mesh/manipulation/checkMesh/Make/files b/applications/utilities/mesh/manipulation/checkMesh/Make/files
index 21febe27891..08cbcaf0885 100644
--- a/applications/utilities/mesh/manipulation/checkMesh/Make/files
+++ b/applications/utilities/mesh/manipulation/checkMesh/Make/files
@@ -1,3 +1,4 @@
+writeFields.C
checkTools.C
checkTopology.C
checkGeometry.C
diff --git a/applications/utilities/mesh/manipulation/checkMesh/checkMesh.C b/applications/utilities/mesh/manipulation/checkMesh/checkMesh.C
index 246d559c547..89aa271bfec 100644
--- a/applications/utilities/mesh/manipulation/checkMesh/checkMesh.C
+++ b/applications/utilities/mesh/manipulation/checkMesh/checkMesh.C
@@ -3,7 +3,7 @@
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
- \\/ M anipulation | Copyright (C) 2015 OpenCFD Ltd.
+ \\/ M anipulation | Copyright (C) 2015-2017 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -57,12 +57,18 @@ Usage
Reconstruct all cellSets and faceSets geometry and write to postProcessing/
directory according to surfaceFormat (e.g. vtk or ensight)
+ \param -writeAllFields \n
+ Writes all mesh quality measures as fields.
+
+ \param -writeFields '(\<fieldName\>)' \n
+ Writes selected mesh quality measures as fields.
+
\*---------------------------------------------------------------------------*/
#include "argList.H"
#include "timeSelector.H"
#include "Time.H"
-#include "polyMesh.H"
+#include "fvMesh.H"
#include "globalMeshData.H"
#include "surfaceWriter.H"
#include "vtkSetWriter.H"
@@ -71,6 +77,7 @@ Usage
#include "checkTopology.H"
#include "checkGeometry.H"
#include "checkMeshQuality.H"
+#include "writeFields.H"
using namespace Foam;
@@ -96,6 +103,17 @@ int main(int argc, char *argv[])
"include extra topology checks"
);
argList::addBoolOption
+ (
+ "writeAllFields",
+ "write volFields with mesh quality parameters"
+ );
+ argList::addOption
+ (
+ "writeFields",
+ "wordList"
+ "write volFields with selected mesh quality parameters"
+ );
+ argList::addBoolOption
(
"meshQuality",
"read user-defined mesh quality criterions from system/meshQualityDict"
@@ -110,7 +128,7 @@ int main(int argc, char *argv[])
#include "setRootCase.H"
#include "createTime.H"
instantList timeDirs = timeSelector::select0(runTime, args);
- #include "createNamedPolyMesh.H"
+ #include "createNamedMesh.H"
const bool noTopology = args.optionFound("noTopology");
const bool allGeometry = args.optionFound("allGeometry");
@@ -119,6 +137,24 @@ int main(int argc, char *argv[])
word surfaceFormat;
const bool writeSets = args.optionReadIfPresent("writeSets", surfaceFormat);
+ HashSet<word> selectedFields;
+ bool writeFields = args.optionReadIfPresent
+ (
+ "writeFields",
+ selectedFields
+ );
+ if (!writeFields && args.optionFound("writeAllFields"))
+ {
+ selectedFields.insert("nonOrthoAngle");
+ selectedFields.insert("faceWeight");
+ selectedFields.insert("skewness");
+ selectedFields.insert("cellDeterminant");
+ selectedFields.insert("aspectRatio");
+ selectedFields.insert("cellShapes");
+ selectedFields.insert("cellVolume");
+ selectedFields.insert("cellVolumeRatio");
+ }
+
if (noTopology)
{
@@ -143,6 +179,11 @@ int main(int argc, char *argv[])
<< " representation"
<< " of all faceSets and cellSets." << nl << endl;
}
+ if (selectedFields.size())
+ {
+ Info<< "Writing mesh quality as fields " << selectedFields << nl
+ << endl;
+ }
autoPtr<IOdictionary> qualDict;
@@ -234,6 +275,10 @@ int main(int argc, char *argv[])
Info<< "\nFailed " << nFailedChecks << " mesh checks.\n"
<< endl;
}
+
+
+ // Write selected fields
+ Foam::writeFields(mesh, selectedFields);
}
else if (state == polyMesh::POINTS_MOVED)
{
@@ -262,6 +307,10 @@ int main(int argc, char *argv[])
{
Info<< "\nMesh OK.\n" << endl;
}
+
+
+ // Write selected fields
+ Foam::writeFields(mesh, selectedFields);
}
}
diff --git a/applications/utilities/mesh/manipulation/checkMesh/writeFields.C b/applications/utilities/mesh/manipulation/checkMesh/writeFields.C
new file mode 100644
index 00000000000..1cdc3dbe61a
--- /dev/null
+++ b/applications/utilities/mesh/manipulation/checkMesh/writeFields.C
@@ -0,0 +1,360 @@
+#include "writeFields.H"
+#include "volFields.H"
+#include "polyMeshTools.H"
+#include "zeroGradientFvPatchFields.H"
+#include "syncTools.H"
+
+using namespace Foam;
+
+void maxFaceToCell
+(
+ const scalarField& faceData,
+ volScalarField& cellData
+)
+{
+ const cellList& cells = cellData.mesh().cells();
+
+ scalarField& cellFld = cellData.ref();
+
+ cellFld = -GREAT;
+ forAll(cells, cellI)
+ {
+ const cell& cFaces = cells[cellI];
+ forAll(cFaces, i)
+ {
+ cellFld[cellI] = max(cellFld[cellI], faceData[cFaces[i]]);
+ }
+ }
+
+ forAll(cellData.boundaryField(), patchI)
+ {
+ fvPatchScalarField& fvp = cellData.boundaryFieldRef()[patchI];
+
+ fvp = fvp.patch().patchSlice(faceData);
+ }
+ cellData.correctBoundaryConditions();
+}
+
+
+void minFaceToCell
+(
+ const scalarField& faceData,
+ volScalarField& cellData
+)
+{
+ const cellList& cells = cellData.mesh().cells();
+
+ scalarField& cellFld = cellData.ref();
+
+ cellFld = GREAT;
+ forAll(cells, cellI)
+ {
+ const cell& cFaces = cells[cellI];
+ forAll(cFaces, i)
+ {
+ cellFld[cellI] = min(cellFld[cellI], faceData[cFaces[i]]);
+ }
+ }
+
+ forAll(cellData.boundaryField(), patchI)
+ {
+ fvPatchScalarField& fvp = cellData.boundaryFieldRef()[patchI];
+
+ fvp = fvp.patch().patchSlice(faceData);
+ }
+ cellData.correctBoundaryConditions();
+}
+
+
+void Foam::writeFields
+(
+ const fvMesh& mesh,
+ const HashSet<word>& selectedFields
+)
+{
+ if (selectedFields.empty())
+ {
+ return;
+ }
+
+ Info<< "Writing fields with mesh quality parameters" << endl;
+
+ if (selectedFields.found("nonOrthoAngle"))
+ {
+ //- Face based orthogonality
+ const scalarField faceOrthogonality
+ (
+ polyMeshTools::faceOrthogonality
+ (
+ mesh,
+ mesh.faceAreas(),
+ mesh.cellCentres()
+ )
+ );
+
+ //- Face based angle
+ const scalarField nonOrthoAngle
+ (
+ radToDeg
+ (
+ Foam::acos(min(1.0, faceOrthogonality))
+ )
+ );
+
+ //- Cell field - max of either face
+ volScalarField cellNonOrthoAngle
+ (
+ IOobject
+ (
+ "nonOrthoAngle",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ mesh,
+ dimensionedScalar("angle", dimless, 0),
+ calculatedFvPatchScalarField::typeName
+ );
+ //- Take max
+ maxFaceToCell(nonOrthoAngle, cellNonOrthoAngle);
+ Info<< " Writing non-orthogonality (angle) to "
+ << cellNonOrthoAngle.name() << endl;
+ cellNonOrthoAngle.write();
+ }
+
+ if (selectedFields.found("faceWeight"))
+ {
+ const scalarField faceWeights
+ (
+ polyMeshTools::faceWeights
+ (
+ mesh,
+ mesh.faceCentres(),
+ mesh.faceAreas(),
+ mesh.cellCentres()
+ )
+ );
+
+ volScalarField cellWeights
+ (
+ IOobject
+ (
+ "faceWeight",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ mesh,
+ dimensionedScalar("weight", dimless, 0),
+ calculatedFvPatchScalarField::typeName
+ );
+ //- Take min
+ minFaceToCell(faceWeights, cellWeights);
+ Info<< " Writing face interpolation weights (0..0.5) to "
+ << cellWeights.name() << endl;
+ cellWeights.write();
+ }
+
+
+ // Skewness
+ // ~~~~~~~~
+
+ if (selectedFields.found("skewness"))
+ {
+ //- Face based skewness
+ const scalarField faceSkewness
+ (
+ polyMeshTools::faceSkewness
+ (
+ mesh,
+ mesh.points(),
+ mesh.faceCentres(),
+ mesh.faceAreas(),
+ mesh.cellCentres()
+ )
+ );
+
+ //- Cell field - max of either face
+ volScalarField cellSkewness
+ (
+ IOobject
+ (
+ "skewness",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ mesh,
+ dimensionedScalar("skewness", dimless, 0),
+ calculatedFvPatchScalarField::typeName
+ );
+ //- Take max
+ maxFaceToCell(faceSkewness, cellSkewness);
+ Info<< " Writing face skewness to " << cellSkewness.name() << endl;
+ cellSkewness.write();
+ }
+
+
+ // cellDeterminant
+ // ~~~~~~~~~~~~~~~
+
+ if (selectedFields.found("cellDeterminant"))
+ {
+ volScalarField cellDeterminant
+ (
+ IOobject
+ (
+ "cellDeterminant",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE,
+ false
+ ),
+ mesh,
+ dimensionedScalar("cellDeterminant", dimless, 0),
+ zeroGradientFvPatchScalarField::typeName
+ );
+ cellDeterminant.primitiveFieldRef() =
+ primitiveMeshTools::cellDeterminant
+ (
+ mesh,
+ mesh.geometricD(),
+ mesh.faceAreas(),
+ syncTools::getInternalOrCoupledFaces(mesh)
+ );
+ cellDeterminant.correctBoundaryConditions();
+ Info<< " Writing cell determinant to "
+ << cellDeterminant.name() << endl;
+ cellDeterminant.write();
+ }
+
+
+ // Aspect ratio
+ // ~~~~~~~~~~~~
+ if (selectedFields.found("aspectRatio"))
+ {
+ volScalarField aspectRatio
+ (
+ IOobject
+ (
+ "aspectRatio",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE,
+ false
+ ),
+ mesh,
+ dimensionedScalar("aspectRatio", dimless, 0),
+ zeroGradientFvPatchScalarField::typeName
+ );
+
+
+ scalarField cellOpenness;
+ polyMeshTools::cellClosedness
+ (
+ mesh,
+ mesh.geometricD(),
+ mesh.faceAreas(),
+ mesh.cellVolumes(),
+ cellOpenness,
+ aspectRatio.ref()
+ );
+
+ aspectRatio.correctBoundaryConditions();
+ Info<< " Writing aspect ratio to " << aspectRatio.name() << endl;
+ aspectRatio.write();
+ }
+
+
+ // cell type
+ // ~~~~~~~~~
+ if (selectedFields.found("cellShapes"))
+ {
+ volScalarField shape
+ (
+ IOobject
+ (
+ "cellShapes",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE,
+ false
+ ),
+ mesh,
+ dimensionedScalar("cellShapes", dimless, 0),
+ zeroGradientFvPatchScalarField::typeName
+ );
+ const cellShapeList& cellShapes = mesh.cellShapes();
+ forAll(cellShapes, cellI)
+ {
+ const cellModel& model = cellShapes[cellI].model();
+ shape[cellI] = model.index();
+ }
+ shape.correctBoundaryConditions();
+ Info<< " Writing cell shape (hex, tet etc.) to " << shape.name()
+ << endl;
+ shape.write();
+ }
+
+ if (selectedFields.found("cellVolume"))
+ {
+ volScalarField V
+ (
+ IOobject
+ (
+ "cellVolume",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE,
+ false
+ ),
+ mesh,
+ dimensionedScalar("cellVolume", dimVolume, 0),
+ calculatedFvPatchScalarField::typeName
+ );
+ V.ref() = mesh.V();
+ Info<< " Writing cell volume to " << V.name() << endl;
+ V.write();
+ }
+
+ if (selectedFields.found("cellVolumeRatio"))
+ {
+ const scalarField faceVolumeRatio
+ (
+ polyMeshTools::volRatio
+ (
+ mesh,
+ mesh.V()
+ )
+ );
+
+ volScalarField cellVolumeRatio
+ (
+ IOobject
+ (
+ "cellVolumeRatio",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ mesh,
+ dimensionedScalar("cellVolumeRatio", dimless, 0),
+ calculatedFvPatchScalarField::typeName
+ );
+ //- Take min
+ minFaceToCell(faceVolumeRatio, cellVolumeRatio);
+ Info<< " Writing cell volume ratio to "
+ << cellVolumeRatio.name() << endl;
+ cellVolumeRatio.write();
+ }
+
+ Info<< endl;
+}
diff --git a/applications/utilities/mesh/manipulation/checkMesh/writeFields.H b/applications/utilities/mesh/manipulation/checkMesh/writeFields.H
new file mode 100644
index 00000000000..09c31b0e62b
--- /dev/null
+++ b/applications/utilities/mesh/manipulation/checkMesh/writeFields.H
@@ -0,0 +1,10 @@
+#include "fvMesh.H"
+
+namespace Foam
+{
+ void writeFields
+ (
+ const fvMesh&,
+ const HashSet<word>& selectedFields
+ );
+}
--
GitLab