From 7fc943c178a36fcebdff05479d420239c32271c2 Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@esi-group.com>
Date: Tue, 12 Oct 2021 08:47:18 +0200
Subject: [PATCH] ENH: improvements for makeFaMesh, checkFaMesh
- added -dry-run, -write-vtk options.
Additional mesh information after creation.
- add parallel reductions and more information for checkFaMesh
ENH: minor cleanup of faPatch internals
- align pointLabels and pointEdges creation more closely with coding
patterns used in PrimitivePatch
- use fileHandler when loading "S0" field.
---
.../finiteArea/checkFaMesh/Make/options | 8 +-
.../finiteArea/checkFaMesh/checkFaMesh.C | 60 ++++----
.../finiteArea/checkFaMesh/faMeshWriteVTK.H | 47 ++++++
.../finiteArea/checkFaMesh/printMeshSummary.H | 134 ++++++++++++++++++
.../finiteArea/makeFaMesh/Make/options | 2 +
.../finiteArea/makeFaMesh/decomposeFaFields.H | 49 +++++--
.../makeFaMesh/faMeshWriteEdgesOBJ.H | 11 +-
.../finiteArea/makeFaMesh/faMeshWriteVTK.H | 47 ++++++
.../finiteArea/makeFaMesh/makeFaMesh.C | 53 ++++---
.../finiteArea/makeFaMesh/printMeshSummary.H | 99 +++++++++++++
.../faMesh/faBoundaryMesh/faBoundaryMesh.C | 9 +-
src/finiteArea/faMesh/faMesh.C | 4 +-
.../constraint/cyclic/cyclicFaPatch.C | 37 +++--
.../constraint/processor/processorFaPatch.C | 7 +-
.../constraint/processor/processorFaPatch.H | 19 +--
.../faMesh/faPatches/faPatch/faPatch.C | 132 +++++++++--------
.../faMesh/faPatches/faPatch/faPatch.H | 15 +-
17 files changed, 561 insertions(+), 172 deletions(-)
create mode 100644 applications/utilities/finiteArea/checkFaMesh/faMeshWriteVTK.H
create mode 100644 applications/utilities/finiteArea/checkFaMesh/printMeshSummary.H
create mode 100644 applications/utilities/finiteArea/makeFaMesh/faMeshWriteVTK.H
create mode 100644 applications/utilities/finiteArea/makeFaMesh/printMeshSummary.H
diff --git a/applications/utilities/finiteArea/checkFaMesh/Make/options b/applications/utilities/finiteArea/checkFaMesh/Make/options
index b69f32d758c..2a0d2e96cef 100644
--- a/applications/utilities/finiteArea/checkFaMesh/Make/options
+++ b/applications/utilities/finiteArea/checkFaMesh/Make/options
@@ -1,9 +1,11 @@
EXE_INC = \
-I$(LIB_SRC)/finiteArea/lnInclude \
- -I$(LIB_SRC)/finiteVolume/lnInclude \
- -I$(LIB_SRC)/meshTools/lnInclude \
+ -I$(LIB_SRC)/fileFormats/lnInclude \
+ -I$(LIB_SRC)/surfMesh/lnInclude \
+ -I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lfiniteArea \
- -lfiniteVolume \
+ -lfileFormats \
+ -lsurfMesh \
-lmeshTools
diff --git a/applications/utilities/finiteArea/checkFaMesh/checkFaMesh.C b/applications/utilities/finiteArea/checkFaMesh/checkFaMesh.C
index d83aaeb7f2e..e15b830083c 100644
--- a/applications/utilities/finiteArea/checkFaMesh/checkFaMesh.C
+++ b/applications/utilities/finiteArea/checkFaMesh/checkFaMesh.C
@@ -6,6 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2016-2017 Wikki Ltd
+ Copyright (C) 2021 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -29,14 +30,22 @@ Application
Description
Check a finiteArea mesh
-Author
+Original Authors
Zeljko Tukovic, FAMENA
Hrvoje Jasak, Wikki Ltd.
\*---------------------------------------------------------------------------*/
-#include "fvCFD.H"
-#include "faCFD.H"
+#include "Time.H"
+#include "argList.H"
+#include "faMesh.H"
+#include "polyMesh.H"
+#include "areaFaMesh.H"
+#include "edgeFaMesh.H"
+#include "areaFields.H"
+#include "edgeFields.H"
+#include "processorFaPatch.H"
+#include "foamVtkUIndPatchWriter.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -51,37 +60,32 @@ int main(int argc, char *argv[])
"Check a finiteArea mesh"
);
- #include "addRegionOption.H"
+ argList::addBoolOption
+ (
+ "write-vtk",
+ "Write mesh as a vtp (vtk) file for display or debugging"
+ );
+ #include "addRegionOption.H"
#include "setRootCase.H"
#include "createTime.H"
- #include "createNamedMesh.H"
- #include "createFaMesh.H"
+ #include "createNamedPolyMesh.H"
+
+
+ // Create
+ faMesh aMesh(mesh);
Info<< "Time = " << runTime.timeName() << nl << endl;
- // General mesh statistics
- Info<< "Number of points: " << aMesh.nPoints() << nl
- << "Number of internal edges: " << aMesh.nInternalEdges() << nl
- << "Number of edges: " << aMesh.nEdges() << nl
- << "Number of faces: " << aMesh.nFaces() << nl
- << endl;
-
- // Check geometry
- Info<< "Face area: min = " << min(aMesh.S().field())
- << " max = " << max(aMesh.S().field()) << nl
- << "Internal edge length: min = "
- << min(aMesh.magLe().internalField()) << nl
- << " max = " << max(aMesh.magLe().internalField()) << nl
- << "Edge length: min = "
- << min(aMesh.magLe()).value() << nl
- << " max = " << max(aMesh.magLe()).value() << nl
- << "Face area normals: min = " << min(aMesh.faceAreaNormals().field())
- << " max = " << max(aMesh.faceAreaNormals().field()) << nl
- << endl;
-
-
- Info << "\nEnd" << endl;
+ #include "printMeshSummary.H"
+
+ if (args.found("write-vtk"))
+ {
+ #include "faMeshWriteVTK.H"
+ }
+
+ Info<< "\nEnd\n" << endl;
+
return 0;
}
diff --git a/applications/utilities/finiteArea/checkFaMesh/faMeshWriteVTK.H b/applications/utilities/finiteArea/checkFaMesh/faMeshWriteVTK.H
new file mode 100644
index 00000000000..1df8b4d3125
--- /dev/null
+++ b/applications/utilities/finiteArea/checkFaMesh/faMeshWriteVTK.H
@@ -0,0 +1,47 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+ VTK output of faMesh with some geometric or debug fields
+
+\*---------------------------------------------------------------------------*/
+
+{
+ // finiteArea
+ vtk::uindirectPatchWriter writer
+ (
+ aMesh.patch(),
+ // vtk::formatType::INLINE_ASCII,
+ fileName
+ (
+ aMesh.mesh().time().globalPath() / "finiteArea"
+ )
+ );
+
+ writer.writeGeometry();
+
+ // CellData
+ writer.beginCellData(3);
+ writer.writeProcIDs();
+ writer.write("area", aMesh.S().field());
+ writer.write("normal", aMesh.faceAreaNormals());
+
+ // PointData
+ writer.beginPointData(1);
+ writer.write("normal", aMesh.pointAreaNormals());
+
+ Info<< nl
+ << "Wrote faMesh in vtk format: " << writer.output().name() << nl;
+}
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/finiteArea/checkFaMesh/printMeshSummary.H b/applications/utilities/finiteArea/checkFaMesh/printMeshSummary.H
new file mode 100644
index 00000000000..256d2b1371f
--- /dev/null
+++ b/applications/utilities/finiteArea/checkFaMesh/printMeshSummary.H
@@ -0,0 +1,134 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+ Summary of faMesh information
+
+\*---------------------------------------------------------------------------*/
+
+{
+ const faBoundaryMesh& patches = aMesh.boundary();
+ const label nNonProcessor = patches.nNonProcessor();
+ const label nPatches = patches.size();
+
+ Info<< "----------------" << nl
+ << "Mesh Information" << nl
+ << "----------------" << nl
+ << " " << "boundingBox: " << boundBox(aMesh.points()) << nl;
+
+ Info<< " Number of points: "
+ << returnReduce(aMesh.nPoints(), sumOp<label>()) << nl
+ << " Number of faces: "
+ << returnReduce(aMesh.nFaces(), sumOp<label>()) << nl;
+
+ Info<< " Number of edges: "
+ << returnReduce(aMesh.nEdges(), sumOp<label>()) << nl
+ << " Number of internal edges: "
+ << returnReduce(aMesh.nInternalEdges(), sumOp<label>()) << nl;
+
+
+ label nProcEdges = 0;
+
+ if (Pstream::parRun())
+ {
+ for (const faPatch& fap : patches)
+ {
+ const auto* cpp = isA<processorFaPatch>(fap);
+
+ if (cpp)
+ {
+ nProcEdges += fap.nEdges();
+ }
+ }
+ }
+
+ const label nBoundEdges = aMesh.nBoundaryEdges() - nProcEdges;
+
+ Info<< " Number of boundary edges: "
+ << returnReduce(nBoundEdges - nProcEdges, sumOp<label>()) << nl;
+
+ if (Pstream::parRun())
+ {
+ Info<< " Number of processor edges: "
+ << returnReduce(nProcEdges, sumOp<label>()) << nl;
+ }
+
+
+ Info<< "----------------" << nl
+ << "Patches" << nl
+ << "----------------" << nl;
+
+ for (label patchi = 0; patchi < nNonProcessor; ++patchi)
+ {
+ const faPatch& p = patches[patchi];
+
+ Info<< " " << "patch " << p.index()
+ << " (size: " << returnReduce(p.size(), sumOp<label>())
+ << ") name: " << p.name()
+ << nl;
+ }
+
+
+ // Geometry information
+ Info<< nl;
+ {
+ scalarMinMax limit(gMinMax(aMesh.S().field()));
+ Info<< "Face area:" << nl
+ << " min = " << limit.min() << " max = " << limit.max() << nl;
+ }
+
+ {
+ scalarMinMax limit(minMax(aMesh.magLe().primitiveField()));
+
+ // Include processor boundaries into 'internal' edges
+ if (Pstream::parRun())
+ {
+ for (label patchi = nNonProcessor; patchi < nPatches; ++patchi)
+ {
+ limit.add(minMax(aMesh.magLe().boundaryField()[patchi]));
+ }
+
+ reduce(limit, minMaxOp<scalar>());
+ }
+
+ Info<< "Edge length (internal):" << nl
+ << " min = " << limit.min() << " max = " << limit.max() << nl;
+
+
+ // Include (non-processor) boundaries
+ for (label patchi = 0; patchi < nNonProcessor; ++patchi)
+ {
+ limit.add(minMax(aMesh.magLe().boundaryField()[patchi]));
+ }
+
+ if (Pstream::parRun())
+ {
+ reduce(limit, minMaxOp<scalar>());
+ }
+
+ Info<< "Edge length:" << nl
+ << " min = " << limit.min() << " max = " << limit.max() << nl;
+ }
+
+ // Not particularly meaningful
+ #if 0
+ {
+ MinMax<vector> limit(gMinMax(aMesh.faceAreaNormals().field()));
+
+ Info<< "Face area normals:" << nl
+ << " min = " << limit.min() << " max = " << limit.max() << nl;
+ }
+ #endif
+}
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/finiteArea/makeFaMesh/Make/options b/applications/utilities/finiteArea/makeFaMesh/Make/options
index 7a28b2c4b67..3f949d53443 100644
--- a/applications/utilities/finiteArea/makeFaMesh/Make/options
+++ b/applications/utilities/finiteArea/makeFaMesh/Make/options
@@ -1,6 +1,7 @@
EXE_INC = \
-I$(LIB_SRC)/finiteArea/lnInclude \
-I$(LIB_SRC)/fileFormats/lnInclude \
+ -I$(LIB_SRC)/surfMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/parallel/decompose/faDecompose/lnInclude \
-I$(LIB_SRC)/parallel/reconstruct/faReconstruct/lnInclude
@@ -8,6 +9,7 @@ EXE_INC = \
EXE_LIBS = \
-lfiniteArea \
-lfileFormats \
+ -lsurfMesh \
-lmeshTools \
-lfaDecompose \
-lfaReconstruct
diff --git a/applications/utilities/finiteArea/makeFaMesh/decomposeFaFields.H b/applications/utilities/finiteArea/makeFaMesh/decomposeFaFields.H
index 5c38c9b4f13..3b91b260ca3 100644
--- a/applications/utilities/finiteArea/makeFaMesh/decomposeFaFields.H
+++ b/applications/utilities/finiteArea/makeFaMesh/decomposeFaFields.H
@@ -17,9 +17,11 @@ Description
if (Pstream::parRun())
{
- faMeshReconstructor reconstructor(areaMesh);
+ faMeshReconstructor reconstructor(aMesh);
reconstructor.writeAddressing();
+ Info<< "Wrote proc-addressing" << nl << endl;
+
// Handle area fields
// ------------------
@@ -72,22 +74,53 @@ if (Pstream::parRun())
if (nAreaFields)
{
- Info<< "Decomposing " << nAreaFields << " area fields" << nl << endl;
+ Info<< "Decomposing " << nAreaFields << " area fields" << nl;
faFieldDecomposer fieldDecomposer
(
fullMesh,
- areaMesh,
+ aMesh,
reconstructor.edgeProcAddressing(),
reconstructor.faceProcAddressing(),
reconstructor.boundaryProcAddressing()
);
- fieldDecomposer.decomposeFields(areaScalarFields);
- fieldDecomposer.decomposeFields(areaVectorFields);
- fieldDecomposer.decomposeFields(areaSphTensorFields);
- fieldDecomposer.decomposeFields(areaSymmTensorFields);
- fieldDecomposer.decomposeFields(areaTensorFields);
+ if (areaScalarFields.size())
+ {
+ Info<< " scalars: "
+ << flatOutput(PtrListOps::names(areaScalarFields)) << nl;
+ fieldDecomposer.decomposeFields(areaScalarFields);
+ }
+
+ if (areaVectorFields.size())
+ {
+ Info<< " vectors: "
+ << flatOutput(PtrListOps::names(areaVectorFields)) << nl;
+ fieldDecomposer.decomposeFields(areaVectorFields);
+ }
+
+ if (areaSphTensorFields.size())
+ {
+ Info<< " sphTensors: "
+ << flatOutput(PtrListOps::names(areaSphTensorFields)) << nl;
+ fieldDecomposer.decomposeFields(areaSphTensorFields);
+ }
+
+ if (areaSymmTensorFields.size())
+ {
+ Info<< " symmTensors: "
+ << flatOutput(PtrListOps::names(areaSymmTensorFields)) << nl;
+ fieldDecomposer.decomposeFields(areaSymmTensorFields);
+ }
+
+ if (areaTensorFields.size())
+ {
+ Info<< " tensors: "
+ << flatOutput(PtrListOps::names(areaTensorFields)) << nl;
+ fieldDecomposer.decomposeFields(areaTensorFields);
+ }
+
+ Info<< endl;
}
}
diff --git a/applications/utilities/finiteArea/makeFaMesh/faMeshWriteEdgesOBJ.H b/applications/utilities/finiteArea/makeFaMesh/faMeshWriteEdgesOBJ.H
index a50d18df3c3..a5f4192d392 100644
--- a/applications/utilities/finiteArea/makeFaMesh/faMeshWriteEdgesOBJ.H
+++ b/applications/utilities/finiteArea/makeFaMesh/faMeshWriteEdgesOBJ.H
@@ -16,15 +16,18 @@ Description
\*---------------------------------------------------------------------------*/
{
- Info<< "Writing edges in obj format" << endl;
+ Info<< nl
+ << "Writing edges in obj format" << endl;
- word outputName("faMesh-edges.obj");
+ word outputName("finiteArea-edges.obj");
if (Pstream::parRun())
{
outputName = word
(
- "faMesh-edges-" + Foam::name(Pstream::myProcNo()) + ".obj"
+ "finiteArea-edges-proc"
+ + Foam::name(Pstream::myProcNo())
+ + ".obj"
);
}
@@ -36,7 +39,7 @@ Description
false
);
- os.write(areaMesh.patch().edges(), areaMesh.patch().localPoints());
+ os.write(aMesh.patch().edges(), aMesh.patch().localPoints());
}
diff --git a/applications/utilities/finiteArea/makeFaMesh/faMeshWriteVTK.H b/applications/utilities/finiteArea/makeFaMesh/faMeshWriteVTK.H
new file mode 100644
index 00000000000..1df8b4d3125
--- /dev/null
+++ b/applications/utilities/finiteArea/makeFaMesh/faMeshWriteVTK.H
@@ -0,0 +1,47 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+ VTK output of faMesh with some geometric or debug fields
+
+\*---------------------------------------------------------------------------*/
+
+{
+ // finiteArea
+ vtk::uindirectPatchWriter writer
+ (
+ aMesh.patch(),
+ // vtk::formatType::INLINE_ASCII,
+ fileName
+ (
+ aMesh.mesh().time().globalPath() / "finiteArea"
+ )
+ );
+
+ writer.writeGeometry();
+
+ // CellData
+ writer.beginCellData(3);
+ writer.writeProcIDs();
+ writer.write("area", aMesh.S().field());
+ writer.write("normal", aMesh.faceAreaNormals());
+
+ // PointData
+ writer.beginPointData(1);
+ writer.write("normal", aMesh.pointAreaNormals());
+
+ Info<< nl
+ << "Wrote faMesh in vtk format: " << writer.output().name() << nl;
+}
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/finiteArea/makeFaMesh/makeFaMesh.C b/applications/utilities/finiteArea/makeFaMesh/makeFaMesh.C
index e80ed303d9f..f93685ea48d 100644
--- a/applications/utilities/finiteArea/makeFaMesh/makeFaMesh.C
+++ b/applications/utilities/finiteArea/makeFaMesh/makeFaMesh.C
@@ -32,7 +32,7 @@ Description
When called in parallel, it will also try to act like decomposePar,
create procAddressing and decompose serial finite-area fields.
-Author
+Original Authors
Zeljko Tukovic, FAMENA
Hrvoje Jasak, Wikki Ltd.
@@ -48,6 +48,8 @@ Author
#include "areaFields.H"
#include "faFieldDecomposer.H"
#include "faMeshReconstructor.H"
+#include "PtrListOps.H"
+#include "foamVtkUIndPatchWriter.H"
#include "OBJstream.H"
using namespace Foam;
@@ -69,11 +71,21 @@ int main(int argc, char *argv[])
);
argList::addOption("dict", "file", "Alternative faMeshDefinition");
+ argList::addBoolOption
+ (
+ "dry-run",
+ "Create but do not write"
+ );
+ argList::addBoolOption
+ (
+ "write-vtk",
+ "Write mesh as a vtp (vtk) file for display or debugging"
+ );
argList::addBoolOption
(
"write-edges-obj",
- "Write mesh edges as obj files and exit",
- false // could make an advanced option
+ "Write mesh edges as obj files (one per processor)",
+ true // advanced option (debugging only)
);
#include "addRegionOption.H"
@@ -81,6 +93,8 @@ int main(int argc, char *argv[])
#include "createTime.H"
#include "createNamedPolyMesh.H"
+ const bool dryrun = args.found("dry-run");
+
// Reading faMeshDefinition dictionary
#include "findMeshDefinitionDict.H"
@@ -91,33 +105,40 @@ int main(int argc, char *argv[])
meshDefDict.add("emptyPatch", patchName, true);
}
+
// Create
- faMesh areaMesh(mesh, meshDefDict);
+ faMesh aMesh(mesh, meshDefDict);
- bool quickExit = false;
+ // Mesh information
+ #include "printMeshSummary.H"
if (args.found("write-edges-obj"))
{
- quickExit = true;
#include "faMeshWriteEdgesOBJ.H"
}
- if (quickExit)
+ if (args.found("write-vtk"))
{
- Info<< "\nEnd\n" << endl;
- return 0;
+ #include "faMeshWriteVTK.H"
}
- // Set the precision of the points data to 10
- IOstream::defaultPrecision(10);
+ if (dryrun)
+ {
+ Info<< "\ndry-run: not writing mesh or decomposing fields\n" << nl;
+ }
+ else
+ {
+ // Set the precision of the points data to 10
+ IOstream::defaultPrecision(10);
- Info<< nl << "Write finite area mesh." << nl;
- areaMesh.write();
- Info<< endl;
+ Info<< nl << "Write finite area mesh." << nl;
+ aMesh.write();
- #include "decomposeFaFields.H"
+ Info<< endl;
+ #include "decomposeFaFields.H"
+ }
- Info << "\nEnd\n" << endl;
+ Info<< "\nEnd\n" << endl;
return 0;
}
diff --git a/applications/utilities/finiteArea/makeFaMesh/printMeshSummary.H b/applications/utilities/finiteArea/makeFaMesh/printMeshSummary.H
new file mode 100644
index 00000000000..5dd75c8804d
--- /dev/null
+++ b/applications/utilities/finiteArea/makeFaMesh/printMeshSummary.H
@@ -0,0 +1,99 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+ Summary of faMesh information
+
+\*---------------------------------------------------------------------------*/
+
+{
+ const faBoundaryMesh& patches = aMesh.boundary();
+ const label nNonProcessor = patches.nNonProcessor();
+ const label nPatches = patches.size();
+
+ Info<< "----------------" << nl
+ << "Mesh Information" << nl
+ << "----------------" << nl
+ << " " << "boundingBox: " << boundBox(aMesh.points()) << nl
+ << " " << "nFaces: " << returnReduce(aMesh.nFaces(), sumOp<label>())
+ << nl;
+
+
+ Info<< "----------------" << nl
+ << "Patches" << nl
+ << "----------------" << nl;
+
+ for (label patchi = 0; patchi < nNonProcessor; ++patchi)
+ {
+ const faPatch& p = patches[patchi];
+
+ Info<< " " << "patch " << p.index()
+ << " (size: " << returnReduce(p.size(), sumOp<label>())
+ << ") name: " << p.name()
+ << nl;
+ }
+
+ // Geometry information
+ Info<< nl;
+ {
+ scalarMinMax limit(gMinMax(aMesh.S().field()));
+ Info<< "Face area:" << nl
+ << " min = " << limit.min() << " max = " << limit.max() << nl;
+ }
+
+ {
+ scalarMinMax limit(minMax(aMesh.magLe().primitiveField()));
+
+ // Include processor boundaries into 'internal' edges
+ if (Pstream::parRun())
+ {
+ for (label patchi = nNonProcessor; patchi < nPatches; ++patchi)
+ {
+ limit.add(minMax(aMesh.magLe().boundaryField()[patchi]));
+ }
+
+ reduce(limit, minMaxOp<scalar>());
+ }
+
+ Info<< "Edge length (internal):" << nl
+ << " min = " << limit.min() << " max = " << limit.max() << nl;
+
+
+ // Include (non-processor) boundaries
+ for (label patchi = 0; patchi < nNonProcessor; ++patchi)
+ {
+ limit.add(minMax(aMesh.magLe().boundaryField()[patchi]));
+ }
+
+ if (Pstream::parRun())
+ {
+ reduce(limit, minMaxOp<scalar>());
+ }
+
+ Info<< "Edge length:" << nl
+ << " min = " << limit.min()
+ << " max = " << limit.max() << nl;
+ }
+
+ // Not particularly meaningful
+ #if 0
+ {
+ MinMax<vector> limit(gMinMax(aMesh.faceAreaNormals().field()));
+
+ Info<< "Face area normals:" << nl
+ << " min = " << limit.min() << " max = " << limit.max() << nl;
+ }
+ #endif
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteArea/faMesh/faBoundaryMesh/faBoundaryMesh.C b/src/finiteArea/faMesh/faBoundaryMesh/faBoundaryMesh.C
index c40ff59d4e6..cd4eee3537a 100644
--- a/src/finiteArea/faMesh/faBoundaryMesh/faBoundaryMesh.C
+++ b/src/finiteArea/faMesh/faBoundaryMesh/faBoundaryMesh.C
@@ -75,8 +75,15 @@ Foam::faBoundaryMesh::faBoundaryMesh
regIOobject(io),
mesh_(mesh)
{
- if (readOpt() == IOobject::MUST_READ)
+ if
+ (
+ readOpt() == IOobject::MUST_READ
+ || readOpt() == IOobject::MUST_READ_IF_MODIFIED
+ )
{
+ // Warn for MUST_READ_IF_MODIFIED
+ warnNoRereading<faBoundaryMesh>();
+
faPatchList& patches = *this;
// Read faPatch list
diff --git a/src/finiteArea/faMesh/faMesh.C b/src/finiteArea/faMesh/faMesh.C
index fea1ea36e86..a98c49060d4 100644
--- a/src/finiteArea/faMesh/faMesh.C
+++ b/src/finiteArea/faMesh/faMesh.C
@@ -343,7 +343,7 @@ Foam::faMesh::faMesh(const polyMesh& pMesh)
// Calculate the geometry for the patches (transformation tensors etc.)
boundary_.calcGeometry();
- if (isFile(pMesh.time().timePath()/mesh().dbDir()/"S0"))
+ if (fileHandler().isFile(pMesh.time().timePath()/"S0"))
{
S0Ptr_ = new DimensionedField<scalar, areaMesh>
(
@@ -508,7 +508,7 @@ Foam::faMesh::faMesh
// Calculate the geometry for the patches (transformation tensors etc.)
boundary_.calcGeometry();
- if (isFile(mesh().time().timePath()/"S0"))
+ if (fileHandler().isFile(pMesh.time().timePath()/"S0"))
{
S0Ptr_ = new DimensionedField<scalar, areaMesh>
(
diff --git a/src/finiteArea/faMesh/faPatches/constraint/cyclic/cyclicFaPatch.C b/src/finiteArea/faMesh/faPatches/constraint/cyclic/cyclicFaPatch.C
index 0bee133c300..32f3e4ca5ce 100644
--- a/src/finiteArea/faMesh/faPatches/constraint/cyclic/cyclicFaPatch.C
+++ b/src/finiteArea/faMesh/faPatches/constraint/cyclic/cyclicFaPatch.C
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2016-2017 Wikki Ltd
- Copyright (C) 2019 OpenCFD Ltd.
+ Copyright (C) 2019-2021 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -47,34 +47,33 @@ const Foam::scalar Foam::cyclicFaPatch::matchTol_ = 1e-3;
void Foam::cyclicFaPatch::calcTransforms()
{
+ const label sizeby2 = this->size()/2;
if (size() > 0)
{
- // const label sizeby2 = this->size()/2;
- pointField half0Ctrs(size()/2);
- pointField half1Ctrs(size()/2);
- for (label i=0; i<size()/2; ++i)
+ pointField half0Ctrs(sizeby2);
+ pointField half1Ctrs(sizeby2);
+ for (label edgei=0; edgei < sizeby2; ++edgei)
{
- half0Ctrs[i] = this->edgeCentres()[i];
- half1Ctrs[i] = this->edgeCentres()[i+size()/2];
+ half0Ctrs[edgei] = this->edgeCentres()[edgei];
+ half1Ctrs[edgei] = this->edgeCentres()[edgei + sizeby2];
}
- vectorField half0Normals(size()/2);
- vectorField half1Normals(size()/2);
+ vectorField half0Normals(sizeby2);
+ vectorField half1Normals(sizeby2);
const vectorField eN(edgeNormals()*magEdgeLengths());
scalar maxMatchError = 0;
label errorEdge = -1;
- for (label edgei = 0; edgei < size()/2; ++edgei)
+ for (label edgei = 0; edgei < sizeby2; ++edgei)
{
half0Normals[edgei] = eN[edgei];
- label nbrEdgei = edgei + size()/2;
- half1Normals[edgei] = eN[nbrEdgei];
+ half1Normals[edgei] = eN[edgei + sizeby2];
scalar magLe = mag(half0Normals[edgei]);
scalar nbrMagLe = mag(half1Normals[edgei]);
- scalar avLe = (magLe + nbrMagLe)/2.0;
+ scalar avLe = 0.5*(magLe + nbrMagLe);
if (magLe < ROOTVSMALL && nbrMagLe < ROOTVSMALL)
{
@@ -101,10 +100,10 @@ void Foam::cyclicFaPatch::calcTransforms()
// Check for error in edge matching
if (maxMatchError > matchTol_)
{
- label nbrEdgei = errorEdge + size()/2;
+ label nbrEdgei = errorEdge + sizeby2;
scalar magLe = mag(half0Normals[errorEdge]);
scalar nbrMagLe = mag(half1Normals[errorEdge]);
- scalar avLe = (magLe + nbrMagLe)/2.0;
+ scalar avLe = 0.5*(magLe + nbrMagLe);
FatalErrorInFunction
<< "edge " << errorEdge
@@ -162,7 +161,7 @@ void Foam::cyclicFaPatch::makeWeights(scalarField& w) const
for (label edgei = 0; edgei < sizeby2; ++edgei)
{
- scalar avL = (magL[edgei] + magL[edgei + sizeby2])/2.0;
+ scalar avL = 0.5*(magL[edgei] + magL[edgei + sizeby2]);
if
(
@@ -191,7 +190,7 @@ void Foam::cyclicFaPatch::makeWeights(scalarField& w) const
// Check for error in matching
if (maxMatchError > matchTol_)
{
- scalar avL = (magL[errorEdge] + magL[errorEdge + sizeby2])/2.0;
+ scalar avL = 0.5*(magL[errorEdge] + magL[errorEdge + sizeby2]);
FatalErrorInFunction
<< "edge " << errorEdge << " and " << errorEdge + sizeby2
@@ -313,7 +312,7 @@ Foam::tmp<Foam::labelField> Foam::cyclicFaPatch::transfer
const label sizeby2 = this->size()/2;
- for (label edgei=0; edgei<sizeby2; ++edgei)
+ for (label edgei=0; edgei < sizeby2; ++edgei)
{
pnf[edgei] = interfaceData[edgei + sizeby2];
pnf[edgei + sizeby2] = interfaceData[edgei];
@@ -345,7 +344,7 @@ Foam::tmp<Foam::labelField> Foam::cyclicFaPatch::internalFieldTransfer
const label sizeby2 = this->size()/2;
- for (label edgei=0; edgei<sizeby2; ++edgei)
+ for (label edgei=0; edgei < sizeby2; ++edgei)
{
pnf[edgei] = iF[edgeCells[edgei + sizeby2]];
pnf[edgei + sizeby2] = iF[edgeCells[edgei]];
diff --git a/src/finiteArea/faMesh/faPatches/constraint/processor/processorFaPatch.C b/src/finiteArea/faMesh/faPatches/constraint/processor/processorFaPatch.C
index d257e4a199c..4ac1c79f227 100644
--- a/src/finiteArea/faMesh/faPatches/constraint/processor/processorFaPatch.C
+++ b/src/finiteArea/faMesh/faPatches/constraint/processor/processorFaPatch.C
@@ -80,12 +80,7 @@ void Foam::processorFaPatch::makeNonGlobalPatchPoints() const
|| !boundaryMesh().mesh()().globalData().nGlobalPoints()
)
{
- nonGlobalPatchPointsPtr_ = new labelList(nPoints());
- labelList& ngpp = *nonGlobalPatchPointsPtr_;
- forAll(ngpp, i)
- {
- ngpp[i] = i;
- }
+ nonGlobalPatchPointsPtr_ = new labelList(identity(nPoints()));
}
else
{
diff --git a/src/finiteArea/faMesh/faPatches/constraint/processor/processorFaPatch.H b/src/finiteArea/faMesh/faPatches/constraint/processor/processorFaPatch.H
index 5924eaae017..b21be3af932 100644
--- a/src/finiteArea/faMesh/faPatches/constraint/processor/processorFaPatch.H
+++ b/src/finiteArea/faMesh/faPatches/constraint/processor/processorFaPatch.H
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2016-2017 Wikki Ltd
- Copyright (C) 2019 OpenCFD Ltd.
+ Copyright (C) 2019-2021 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -177,19 +177,19 @@ public:
// Member Functions
//- Return interface size
- virtual label interfaceSize() const
+ virtual label interfaceSize() const noexcept
{
return size();
}
//- Return processor number
- int myProcNo() const
+ int myProcNo() const noexcept
{
return myProcNo_;
}
//- Return neighbour processor number
- int neighbProcNo() const
+ int neighbProcNo() const noexcept
{
return neighbProcNo_;
}
@@ -197,18 +197,11 @@ public:
//- Return true if running parallel
virtual bool coupled() const
{
- if (Pstream::parRun())
- {
- return true;
- }
- else
- {
- return false;
- }
+ return Pstream::parRun();
}
//- Is this the master side?
- virtual bool master() const
+ virtual bool master() const noexcept
{
return (myProcNo_ < neighbProcNo_);
}
diff --git a/src/finiteArea/faMesh/faPatches/faPatch/faPatch.C b/src/finiteArea/faMesh/faPatches/faPatch/faPatch.C
index 81ddba7ee6d..f6e8b98aae9 100644
--- a/src/finiteArea/faMesh/faPatches/faPatch/faPatch.C
+++ b/src/finiteArea/faMesh/faPatches/faPatch/faPatch.C
@@ -210,96 +210,102 @@ const Foam::labelList& Foam::faPatch::pointLabels() const
}
-void Foam::faPatch::calcPointLabels() const
+const Foam::labelListList& Foam::faPatch::pointEdges() const
{
- SLList<label> labels;
+ if (!pointEdgesPtr_)
+ {
+ calcPointEdges();
+ }
+
+ return *pointEdgesPtr_;
+}
+
+void Foam::faPatch::calcPointLabels() const
+{
const edgeList::subList edges = patchSlice(boundaryMesh().mesh().edges());
- forAll(edges, edgeI)
- {
- bool existStart = false;
- bool existEnd = false;
+ // Walk boundary edges.
+ // The edge orientation corresponds to the face orientation
+ // (outwards normal).
- forAllIters(labels, iter)
- {
- if (*iter == edges[edgeI].start())
- {
- existStart = true;
- }
-
- if (*iter == edges[edgeI].end())
- {
- existEnd = true;
- }
- }
+ // Note: could combine this with calcPointEdges for more efficiency
+
+ // Map<label> markedPoints(4*edges.size());
+ labelHashSet markedPoints(4*edges.size());
+ DynamicList<label> dynEdgePoints(2*edges.size());
- if (!existStart)
+ for (const edge& e : edges)
+ {
+ // if (markedPoints.insert(e.first(), markedPoints.size()))
+ if (markedPoints.insert(e.first()))
{
- labels.append(edges[edgeI].start());
+ dynEdgePoints.append(e.first());
}
-
- if (!existEnd)
+ // if (markedPoints.insert(e.second(), markedPoints.size()))
+ if (markedPoints.insert(e.second()))
{
- labels.append(edges[edgeI].end());
+ dynEdgePoints.append(e.second());
}
}
- pointLabelsPtr_ = new labelList(labels);
+ // Transfer to plain list (reuse storage)
+ pointLabelsPtr_ = new labelList(std::move(dynEdgePoints));
+ /// const labelList& edgePoints = *pointLabelsPtr_;
+ ///
+ /// // Cannot use invertManyToMany - we have non-local edge numbering
+ ///
+ /// // Intermediate storage for pointEdges.
+ /// // Points on the boundary will normally connect 1 or 2 edges only.
+ /// List<DynamicList<label,2>> dynPointEdges(edgePoints.size());
+ ///
+ /// forAll(edges, edgei)
+ /// {
+ /// const edge& e = edges[edgei];
+ ///
+ /// dynPointEdges[markedPoints[e.first()]].append(edgei);
+ /// dynPointEdges[markedPoints[e.second()]].append(edgei);
+ /// }
+ ///
+ /// // Flatten to regular list
+ /// pointEdgesPtr_ = new labelListList(edgePoints.size());
+ /// auto& pEdges = *pointEdgesPtr_;
+ ///
+ /// forAll(pEdges, pointi)
+ /// {
+ /// pEdges[pointi] = std::move(dynPointEdges[pointi]);
+ /// }
}
void Foam::faPatch::calcPointEdges() const
{
- const labelList& points = pointLabels();
+ const edgeList::subList edges = patchSlice(boundaryMesh().mesh().edges());
- const edgeList::subList e = patchSlice(boundaryMesh().mesh().edges());
+ const labelList& edgePoints = pointLabels();
- // set up storage for pointEdges
- List<SLList<label>> pointEdgs(points.size());
+ // Cannot use invertManyToMany - we have non-local edge numbering
- forAll(e, edgeI)
- {
- const edge& curPoints = e[edgeI];
+ // Intermediate storage for pointEdges.
+ // Points on the boundary will normally connect 1 or 2 edges only.
+ List<DynamicList<label,2>> dynPointEdges(edgePoints.size());
- forAll(curPoints, pointI)
- {
- const label localPointIndex = points.find(curPoints[pointI]);
-
- pointEdgs[localPointIndex].append(edgeI);
- }
- }
-
- // sort out the list
- pointEdgesPtr_ = new labelListList(pointEdgs.size());
- labelListList& pEdges = *pointEdgesPtr_;
-
- forAll(pointEdgs, pointI)
+ forAll(edges, edgei)
{
- pEdges[pointI].setSize(pointEdgs[pointI].size());
+ const edge& e = edges[edgei];
- label i = 0;
- for
- (
- SLList<label>::iterator curEdgesIter = pointEdgs[pointI].begin();
- curEdgesIter != pointEdgs[pointI].end();
- ++curEdgesIter, ++i
- )
- {
- pEdges[pointI][i] = curEdgesIter();
- }
+ dynPointEdges[edgePoints.find(e.first())].append(edgei);
+ dynPointEdges[edgePoints.find(e.second())].append(edgei);
}
-}
+ // Flatten to regular list
+ pointEdgesPtr_ = new labelListList(edgePoints.size());
+ auto& pEdges = *pointEdgesPtr_;
-const Foam::labelListList& Foam::faPatch::pointEdges() const
-{
- if (!pointEdgesPtr_)
+ forAll(pEdges, pointi)
{
- calcPointEdges();
+ pEdges[pointi] = std::move(dynPointEdges[pointi]);
}
-
- return *pointEdgesPtr_;
}
@@ -383,7 +389,7 @@ const Foam::scalarField& Foam::faPatch::magEdgeLengths() const
Foam::tmp<Foam::vectorField> Foam::faPatch::edgeNormals() const
{
- tmp<vectorField> tedgeNorm(edgeLengths());
+ auto tedgeNorm = tmp<vectorField>::New(edgeLengths());
for (vector& n : tedgeNorm.ref())
{
diff --git a/src/finiteArea/faMesh/faPatches/faPatch/faPatch.H b/src/finiteArea/faMesh/faPatches/faPatch/faPatch.H
index e2529421185..71df3aa5004 100644
--- a/src/finiteArea/faMesh/faPatches/faPatch/faPatch.H
+++ b/src/finiteArea/faMesh/faPatches/faPatch/faPatch.H
@@ -80,17 +80,14 @@ class faPatch
//- Reference to boundary mesh
const faBoundaryMesh& boundaryMesh_;
+ //- Demand-driven: edge-face addressing
+ mutable labelList::subList* edgeFacesPtr_;
- // Demand-driven private data
+ //- Demand-driven: local points labels
+ mutable labelList* pointLabelsPtr_;
- //- Edge-face addressing
- mutable labelList::subList* edgeFacesPtr_;
-
- //- Local points labels
- mutable labelList* pointLabelsPtr_;
-
- //- Point-edge addressing
- mutable labelListList* pointEdgesPtr_;
+ //- Demand-driven: point-edge addressing
+ mutable labelListList* pointEdgesPtr_;
// Private Member Functions
--
GitLab