ENH: generate local orientation vectors for blockMesh VTK output (#1903)

- can help when designing/debugging blockMesh layouts

- propagate low-level cellModel methods face() and edge() to cellShape

STYLE: relocate blockMesh OBJ output to application only

- remove blockTopology files in cleanCase function

- improve code consistency in top-level blockMesh, PDRblockMesh
  generation.

applications/utilities/mesh/generation/PDRblockMesh/PDRblockMesh.C

@@ -43,8 +43,8 @@ Usage
       - \par -dict \<filename\>
         Alternative dictionary for the mesh description.
 
-      - \par -noClean
+      - \par -no-clean
-        Do not remove any existing polyMesh/ directory or files
+        Do not remove polyMesh/ directory or files
 
       - \par -time
         Write resulting mesh to a time directory (instead of constant)
@@ -92,9 +92,11 @@ int main(int argc, char *argv[])
 
     argList::addBoolOption
     (
-        "noClean",
+        "no-clean",
-        "Do not remove any existing polyMesh/ directory or files"
+        "Do not remove polyMesh/ directory or files"
     );
+    argList::addOptionCompat("no-clean", {"noClean", -2006});
+
     argList::addOption("dict", "file", "Alternative PDRblockMeshDict");
     argList::addOption
     (
@@ -107,7 +109,11 @@ int main(int argc, char *argv[])
     #include "createTime.H"
 
     // Remove old files, unless disabled
-    const bool removeOldFiles = !args.found("noClean");
+    const bool removeOldFiles = !args.found("no-clean");
+
+    const word regionName(polyMesh::defaultRegion);
+    const word regionPath;
+
 
     // Instance for resulting mesh
     bool useTime = false;
@@ -134,8 +140,8 @@ int main(int argc, char *argv[])
     }
 
 
+    // Locate appropriate PDRblockMeshDict
     const word dictName("PDRblockMeshDict");
-
     #include "setSystemRunTimeDictionaryIO.H"
 
     IOdictionary meshDict(dictIO);
@@ -159,35 +165,18 @@ int main(int argc, char *argv[])
 
     if (removeOldFiles)
     {
-        const fileName polyMeshPath
+        #include "cleanMeshDirectory.H"
-        (
-            runTime.path()/meshInstance/polyMesh::meshSubDir
-        );
-
-        if (exists(polyMeshPath))
-        {
-            Info<< "Deleting polyMesh directory "
-                << runTime.relativePath(polyMeshPath) << endl;
-            rmDir(polyMeshPath);
-        }
     }
 
 
     Info<< nl << "Creating polyMesh from PDRblockMesh" << endl;
 
-    auto meshPtr = blkMesh.mesh
+    autoPtr<polyMesh> meshPtr =
+        blkMesh.mesh
         (
-        IOobject
+            IOobject(regionName, meshInstance, runTime)
-        (
-            polyMesh::defaultRegion,
-            meshInstance,
-            runTime,
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        )
         );
 
-
     const polyMesh& mesh = *meshPtr;
 
     // Set the precision of the points data to 10
@@ -204,30 +193,7 @@ int main(int argc, char *argv[])
             << exit(FatalError);
     }
 
-    // Mesh summary
+    #include "printMeshSummary.H"
-    {
-        Info<< "----------------" << nl
-            << "Mesh Information" << nl
-            << "----------------" << nl
-            << "  " << "boundingBox: " << boundBox(mesh.points()) << nl
-            << "  " << "nPoints: " << mesh.nPoints() << nl
-            << "  " << "nCells: " << mesh.nCells() << nl
-            << "  " << "nFaces: " << mesh.nFaces() << nl
-            << "  " << "nInternalFaces: " << mesh.nInternalFaces() << nl;
-
-        Info<< "----------------" << nl
-            << "Patches" << nl
-            << "----------------" << nl;
-
-        for (const polyPatch& p : mesh.boundaryMesh())
-        {
-            Info<< "  " << "patch " << p.index()
-                << " (start: " << p.start()
-                << " size: " << p.size()
-                << ") name: " << p.name()
-                << nl;
-        }
-    }
 
     Info<< "\nEnd\n" << endl;
 

applications/utilities/mesh/generation/PDRblockMesh/cleanMeshDirectory.H

@@ -0,0 +1,45 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2020 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Removal of polyMesh directory
+
+\*---------------------------------------------------------------------------*/
+
+{
+    // Shadows enclosing parameter (dictName)
+    const word blockMeshDictName("blockMeshDict");
+
+    const fileName polyMeshPath
+    (
+        runTime.path()/meshInstance/regionPath/polyMesh::meshSubDir
+    );
+
+    if (exists(polyMeshPath))
+    {
+        if (exists(polyMeshPath/blockMeshDictName))
+        {
+            Info<< "Not deleting polyMesh directory "
+                << runTime.relativePath(polyMeshPath) << nl
+                << "    because it contains " << blockMeshDictName << endl;
+        }
+        else
+        {
+            Info<< "Deleting polyMesh directory "
+                << runTime.relativePath(polyMeshPath) << endl;
+            rmDir(polyMeshPath);
+        }
+    }
+}
+
+
+// ************************************************************************* //

applications/utilities/mesh/generation/PDRblockMesh/printMeshSummary.H

@@ -0,0 +1,85 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2020 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Summary of mesh information (eg, after blockMesh)
+
+\*---------------------------------------------------------------------------*/
+
+{
+    Info<< "----------------" << nl
+        << "Mesh Information" << nl
+        << "----------------" << nl
+        << "  " << "boundingBox: " << boundBox(mesh.points()) << nl
+        << "  " << "nPoints: " << mesh.nPoints() << nl
+        << "  " << "nCells: " << mesh.nCells() << nl
+        << "  " << "nFaces: " << mesh.nFaces() << nl
+        << "  " << "nInternalFaces: " << mesh.nInternalFaces() << nl;
+
+    const auto printZone =
+        [](const Foam::zone& zn)
+        {
+            Info<< "  " << "zone " << zn.index()
+                << " (size: " << zn.size()
+                << ") name: " << zn.name() << nl;
+        };
+
+    if (mesh.cellZones().size())
+    {
+        Info<< "----------------" << nl
+            << "Cell Zones" << nl
+            << "----------------" << nl;
+
+        for (const cellZone& zn : mesh.cellZones())
+        {
+            printZone(zn);
+        }
+    }
+    if (mesh.faceZones().size())
+    {
+        Info<< "----------------" << nl
+            << "Face Zones" << nl
+            << "----------------" << nl;
+
+        for (const faceZone& zn : mesh.faceZones())
+        {
+            printZone(zn);
+        }
+    }
+    if (mesh.pointZones().size())
+    {
+        Info<< "----------------" << nl
+            << "Point Zones" << nl
+            << "----------------" << nl;
+
+        for (const pointZone& zn : mesh.pointZones())
+        {
+            printZone(zn);
+        }
+    }
+
+    Info<< "----------------" << nl
+        << "Patches" << nl
+        << "----------------" << nl;
+
+    for (const polyPatch& p : mesh.boundaryMesh())
+    {
+        Info<< "  " << "patch " << p.index()
+            << " (start: " << p.start()
+            << " size: " << p.size()
+            << ") name: " << p.name()
+            << nl;
+    }
+}
+
+
+// ************************************************************************* //

applications/utilities/mesh/generation/blockMesh/addCellZones.H

@@ -1 +0,0 @@
-#warning File removed - left for old dependency check only

applications/utilities/mesh/generation/blockMesh/blockMesh.C

@@ -61,8 +61,8 @@ Usage
       - \par -sets
         Write cellZones as cellSets too (for processing purposes)
 
-      - \par -noClean
+      - \par -no-clean
-        Do not remove any existing polyMesh/ directory or files
+        Do not remove polyMesh/ directory or files
 
       - \par -time
         Write resulting mesh to a time directory (instead of constant)
@@ -141,9 +141,11 @@ int main(int argc, char *argv[])
     );
     argList::addBoolOption
     (
-        "noClean",
+        "no-clean",
-        "Do not remove any existing polyMesh/ directory or files"
+        "Do not remove polyMesh/ directory or files"
     );
+    argList::addOptionCompat("no-clean", {"noClean", -2006});
+
     argList::addOption("dict", "file", "Alternative blockMeshDict");
     argList::addBoolOption
     (
@@ -162,7 +164,7 @@ int main(int argc, char *argv[])
     #include "createTime.H"
 
     // Remove old files, unless disabled
-    const bool removeOldFiles = !args.found("noClean");
+    const bool removeOldFiles = !args.found("no-clean");
 
     // Write cellSets
     const bool writeCellSets = args.found("sets");
@@ -232,66 +234,20 @@ int main(int argc, char *argv[])
     if (args.found("write-obj"))
     {
         quickExit = true;
-
         Info<< nl;
-
+        #include "blockMeshOBJ.H"
-        // Write mesh as edges
-        {
-            OFstream os(runTime.path()/"blockTopology.obj");
-
-            Info<< "Writing block structure in obj format: "
-                << os.name().name() << endl;
-
-            blocks.writeTopology(os);
-        }
-
-        // Write centres of blocks
-        {
-            OFstream os(runTime.path()/"blockCentres.obj");
-
-            Info<< "Writing block centres in obj format: "
-                << os.name().name() << endl;
-
-            for (const point& cc : blocks.topology().cellCentres())
-            {
-                os << "v " << cc.x() << ' ' << cc.y() << ' ' << cc.z() << nl;
-            }
-        }
     }
 
     if (args.found("write-vtk"))
     {
         quickExit = true;
-
-        // non-legacy and ASCII (mesh is small, want readable output)
-        const vtk::outputOptions writeOpts = vtk::formatType::INLINE_ASCII;
-
         Info<< nl;
-
+        #include "blockMeshVTK.H"
-        const polyMesh& topoMesh = blocks.topology();
-        const vtk::vtuCells topoCells(topoMesh, writeOpts);
-
-        vtk::internalMeshWriter writer
-        (
-            topoMesh,
-            topoCells,
-            writeOpts,
-            runTime.path()/"blockTopology"
-        );
-
-        Info<< "Writing block topology in vtk format: "
-            << args.relativePath(writer.output()).c_str() << endl;
-
-        writer.writeGeometry();
-        writer.beginCellData();
-        writer.writeCellIDs();
     }
 
-
     if (quickExit)
     {
         Info<< "\nEnd\n" << endl;
-
         return 0;
     }
 
@@ -309,26 +265,7 @@ int main(int argc, char *argv[])
 
     if (removeOldFiles)
     {
-        const fileName polyMeshPath
+        #include "cleanMeshDirectory.H"
-        (
-            runTime.path()/meshInstance/regionPath/polyMesh::meshSubDir
-        );
-
-        if (exists(polyMeshPath))
-        {
-            if (exists(polyMeshPath/dictName))
-            {
-                Info<< "Not deleting polyMesh directory "
-                    << runTime.relativePath(polyMeshPath) << nl
-                    << "    because it contains " << dictName << endl;
-            }
-            else
-            {
-                Info<< "Deleting polyMesh directory "
-                    << runTime.relativePath(polyMeshPath) << endl;
-                rmDir(polyMeshPath);
-            }
-        }
     }
 
 
@@ -377,30 +314,7 @@ int main(int argc, char *argv[])
         }
     }
 
-    // Write summary
+    #include "printMeshSummary.H"
-    {
-        Info<< "----------------" << nl
-            << "Mesh Information" << nl
-            << "----------------" << nl
-            << "  " << "boundingBox: " << boundBox(mesh.points()) << nl
-            << "  " << "nPoints: " << mesh.nPoints() << nl
-            << "  " << "nCells: " << mesh.nCells() << nl
-            << "  " << "nFaces: " << mesh.nFaces() << nl
-            << "  " << "nInternalFaces: " << mesh.nInternalFaces() << nl;
-
-        Info<< "----------------" << nl
-            << "Patches" << nl
-            << "----------------" << nl;
-
-        for (const polyPatch& p : mesh.boundaryMesh())
-        {
-            Info<< "  " << "patch " << p.index()
-                << " (start: " << p.start()
-                << " size: " << p.size()
-                << ") name: " << p.name()
-                << nl;
-        }
-    }
 
     Info<< "\nEnd\n" << endl;
 

applications/utilities/mesh/generation/blockMesh/blockMeshOBJ.H

@@ -0,0 +1,54 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2020 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    OBJ output of blockMesh topology blocks
+
+\*---------------------------------------------------------------------------*/
+
+{
+    const polyMesh& topoMesh = blocks.topology();
+
+    // Write mesh as edges
+    {
+        OFstream os(runTime.path()/"blockTopology.obj");
+
+        Info<< "Writing block structure in obj format: "
+            << os.name().name() << endl;
+
+        for (const point& p : topoMesh.points())
+        {
+            os << "v " << p.x() << ' ' << p.y() << ' ' << p.z() << nl;
+        }
+
+        for (const edge& e : topoMesh.edges())
+        {
+            os << "l " << e.start() + 1 << ' ' << e.end() + 1 << nl;
+        }
+    }
+
+    // Write centres of blocks
+    {
+        OFstream os(runTime.path()/"blockCentres.obj");
+
+        Info<< "Writing block centres in obj format: "
+            << os.name().name() << endl;
+
+        for (const point& p : topoMesh.cellCentres())
+        {
+            os << "v " << p.x() << ' ' << p.y() << ' ' << p.z() << nl;
+        }
+    }
+}
+
+
+// ************************************************************************* //

applications/utilities/mesh/generation/blockMesh/blockMeshVTK.H

@@ -0,0 +1,100 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2020 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    VTK output of blockMesh topology blocks
+
+\*---------------------------------------------------------------------------*/
+
+{
+    // Non-legacy and ASCII (mesh is small, want readable output)
+    const vtk::outputOptions writeOpts = vtk::formatType::INLINE_ASCII;
+
+    const polyMesh& topoMesh = blocks.topology();
+    const vtk::vtuCells topoCells(topoMesh, writeOpts);
+
+    vtk::internalMeshWriter writer
+    (
+        topoMesh,
+        topoCells,
+        writeOpts,
+        runTime.path()/"blockTopology"
+    );
+
+    Info<< "Writing block topology in vtk format: "
+        << args.relativePath(writer.output()).c_str() << endl;
+
+    writer.writeGeometry();
+    writer.beginCellData();
+    writer.writeCellIDs();
+
+    // No cell decomposition, so there is a 1-to-1 correspondence between
+    // input blocks and VTK output cells.
+
+    vectorField localNormal(blocks.size());
+
+    // Generate local normals as fields for visualisation
+    for (direction cmpt = 0; cmpt < vector::nComponents; ++cmpt)
+    {
+        const label faceMin = label(2*cmpt);
+        const label faceMax = faceMin+1;
+
+        localNormal.resize(blocks.size());
+
+        forAll(blocks, blocki)
+        {
+            const cellShape& shape = blocks[blocki].blockShape();
+            const pointField& verts = blocks[blocki].vertices();
+
+            if (shape.model() == cellModel::ref(cellModel::HEX))
+            {
+                localNormal[blocki] =
+                (
+                    // 30% cell-width as arbitrary value for vector length
+                    0.3*mag
+                    (
+                        shape.face(faceMax).centre(verts)
+                      - shape.face(faceMin).centre(verts)
+                    )
+                  * normalised
+                    (
+                        // Weigh by area to avoid influence of zero-faces
+                        shape.face(faceMax).areaNormal(verts)
+                      - shape.face(faceMin).areaNormal(verts)
+                    )
+                );
+            }
+            else
+            {
+                // Could handle other shapes (how?)
+                localNormal[blocki] = Zero;
+            }
+        }
+
+        // Additional safety (should not happen)
+        localNormal.resize(topoMesh.nCells(), Zero);
+
+        writer.writeCellData
+        (
+            word("local-direction" + name(cmpt)),
+            localNormal
+        );
+    }
+
+    // if (topoMesh.nCells() != blocks.size())
+    // {
+    //     Info<< "Warning: indicated normals may be incorrect" << nl;
+    // }
+}
+
+
+// ************************************************************************* //

applications/utilities/mesh/generation/blockMesh/cleanMeshDirectory.H

@@ -0,0 +1,45 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2020 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Removal of polyMesh directory
+
+\*---------------------------------------------------------------------------*/
+
+{
+    // Shadows enclosing parameter (dictName)
+    const word blockMeshDictName("blockMeshDict");
+
+    const fileName polyMeshPath
+    (
+        runTime.path()/meshInstance/regionPath/polyMesh::meshSubDir
+    );
+
+    if (exists(polyMeshPath))
+    {
+        if (exists(polyMeshPath/blockMeshDictName))
+        {
+            Info<< "Not deleting polyMesh directory "
+                << runTime.relativePath(polyMeshPath) << nl
+                << "    because it contains " << blockMeshDictName << endl;
+        }
+        else
+        {
+            Info<< "Deleting polyMesh directory "
+                << runTime.relativePath(polyMeshPath) << endl;
+            rmDir(polyMeshPath);
+        }
+    }
+}
+
+
+// ************************************************************************* //

applications/utilities/mesh/generation/blockMesh/printMeshSummary.H

@@ -0,0 +1,85 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2020 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Summary of mesh information (eg, after blockMesh)
+
+\*---------------------------------------------------------------------------*/
+
+{
+    Info<< "----------------" << nl
+        << "Mesh Information" << nl
+        << "----------------" << nl
+        << "  " << "boundingBox: " << boundBox(mesh.points()) << nl
+        << "  " << "nPoints: " << mesh.nPoints() << nl
+        << "  " << "nCells: " << mesh.nCells() << nl
+        << "  " << "nFaces: " << mesh.nFaces() << nl
+        << "  " << "nInternalFaces: " << mesh.nInternalFaces() << nl;
+
+    const auto printZone =
+        [](const Foam::zone& zn)
+        {
+            Info<< "  " << "zone " << zn.index()
+                << " (size: " << zn.size()
+                << ") name: " << zn.name() << nl;
+        };
+
+    if (mesh.cellZones().size())
+    {
+        Info<< "----------------" << nl
+            << "Cell Zones" << nl
+            << "----------------" << nl;
+
+        for (const cellZone& zn : mesh.cellZones())
+        {
+            printZone(zn);
+        }
+    }
+    if (mesh.faceZones().size())
+    {
+        Info<< "----------------" << nl
+            << "Face Zones" << nl
+            << "----------------" << nl;
+
+        for (const faceZone& zn : mesh.faceZones())
+        {
+            printZone(zn);
+        }
+    }
+    if (mesh.pointZones().size())
+    {
+        Info<< "----------------" << nl
+            << "Point Zones" << nl
+            << "----------------" << nl;
+
+        for (const pointZone& zn : mesh.pointZones())
+        {
+            printZone(zn);
+        }
+    }
+
+    Info<< "----------------" << nl
+        << "Patches" << nl
+        << "----------------" << nl;
+
+    for (const polyPatch& p : mesh.boundaryMesh())
+    {
+        Info<< "  " << "patch " << p.index()
+            << " (start: " << p.start()
+            << " size: " << p.size()
+            << ") name: " << p.name()
+            << nl;
+    }
+}
+
+
+// ************************************************************************* //

applications/utilities/surface/surfaceClean/surfaceClean.C

@@ -69,9 +69,11 @@ int main(int argc, char *argv[])
 
     argList::addBoolOption
     (
-        "noClean",
+        "no-clean",
         "Suppress surface checking/cleanup on the input surface"
     );
+    argList::addOptionCompat("no-clean", {"noClean", -2006});
+
     argList::addOption
     (
         "scale",
@@ -101,7 +103,7 @@ int main(int argc, char *argv[])
     );
     surf.writeStats(Info);
 
-    if (!args.found("noClean"))
+    if (!args.found("no-clean"))
     {
         Info<< "Removing duplicate and illegal triangles ..." << nl << endl;
         surf.cleanup(true);

bin/tools/CleanFunctions

@@ -6,11 +6,10 @@
 #    \\/     M anipulation  |
 #------------------------------------------------------------------------------
 #     Copyright (C) 2011-2016 OpenFOAM Foundation
-#     Copyright (C) 2015-2019 OpenCFD Ltd.
+#     Copyright (C) 2015-2020 OpenCFD Ltd.
 #------------------------------------------------------------------------------
 # License
-#     This file is part of OpenFOAM, licensed under GNU General Public License
+#     This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
-#     <http://www.gnu.org/licenses/>.
 #
 # Script
 #     CleanFunctions
@@ -116,6 +115,9 @@ cleanCase()
     rm -rf sets
     rm -rf system/machines
 
+    # Possible blockMesh output
+    rm -f blockTopology.vtu blockTopology.obj blockCentres.obj
+
     # From mpirunDebug
     rm -f gdbCommands  mpirun.schema
 

src/OpenFOAM/meshes/meshShapes/cellModel/cellModel.H

@@ -214,7 +214,7 @@ public:
         ) const;
 
 
-        //- Vector centroid
+        //- Centroid of the cell
         vector centre
         (
             const labelList& pointLabels,

src/OpenFOAM/meshes/meshShapes/cellShape/cellShape.H

@@ -126,12 +126,21 @@ public:
 
     // Member Functions
 
-        //- Return the points corresponding to this cellShape
-        inline pointField points(const UList<point>& meshPoints) const;
-
         //- Model reference
         inline const cellModel& model() const;
 
+        //- Number of points
+        inline label nPoints() const;
+
+        //- Number of edges
+        inline label nEdges() const;
+
+        //- Number of faces
+        inline label nFaces() const;
+
+        //- The points corresponding to this shape
+        inline pointField points(const UList<point>& meshPoints) const;
+
         //- Mesh face labels of this cell (in order of model)
         inline labelList meshFaces
         (
@@ -146,29 +155,25 @@ public:
             const labelList& cEdges
         ) const;
 
+        //- The face for the specified model face
+        inline Foam::face face(const label modelFacei) const;
+
         //- Faces of this cell
         inline faceList faces() const;
 
         //- Collapsed faces of this cell
         inline faceList collapsedFaces() const;
 
-        //- Number of faces
+        //- The edge for the specified model edge
-        inline label nFaces() const;
+        inline Foam::edge edge(const label modelEdgei) const;
 
-        //- Edges of this cellShape
+        //- Edges of this shape
         inline edgeList edges() const;
 
-        //- Number of edges
-        inline label nEdges() const;
-
-        //- Number of points
-        inline label nPoints() const;
-
         //- Centroid of the cell
         inline point centre(const UList<point>& points) const;
 
-        //- Return info proxy.
+        //- Return info proxy, to print information to a stream
-        //  Used to print token information to a stream
         Foam::InfoProxy<cellShape> info() const
         {
             return *this;

src/OpenFOAM/meshes/meshShapes/cellShape/cellShapeI.H

@@ -29,6 +29,7 @@ License
 #include "Istream.H"
 #include "cell.H"
 #include "cellModel.H"
+#include "UIndirectList.H"
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
@@ -122,29 +123,36 @@ inline Foam::autoPtr<Foam::cellShape> Foam::cellShape::clone() const
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+inline const Foam::cellModel& Foam::cellShape::model() const
+{
+    return *m;
+}
+
+
+inline Foam::label Foam::cellShape::nPoints() const
+{
+    return size();
+}
+
+
+inline Foam::label Foam::cellShape::nEdges() const
+{
+    return m->nEdges();
+}
+
+
+inline Foam::label Foam::cellShape::nFaces() const
+{
+    return m->nFaces();
+}
+
+
 inline Foam::pointField Foam::cellShape::points
 (
     const UList<point>& meshPoints
 ) const
 {
-    // There are as many points as there labels for them
+    return pointField(UIndirectList<point>(meshPoints, *this));
-    pointField p(size());
-
-    // For each point in list, set it to the point in 'pnts' addressed
-    // by 'labs'
-    forAll(p, i)
-    {
-        p[i] = meshPoints[operator[](i)];
-    }
-
-    // Return list
-    return p;
-}
-
-
-inline const Foam::cellModel& Foam::cellShape::model() const
-{
-    return *m;
 }
 
 
@@ -163,16 +171,13 @@ inline Foam::labelList Foam::cellShape::meshFaces
 
     forAll(localFaces, i)
     {
-        const face& localF = localFaces[i];
+        const Foam::face& localF = localFaces[i];
 
-        forAll(cFaces, j)
+        for (const label meshFacei : cFaces)
         {
-            label meshFacei = cFaces[j];
-
             if (allFaces[meshFacei] == localF)
             {
                 modelToMesh[i] = meshFacei;
-
                 break;
             }
         }
@@ -197,16 +202,13 @@ inline Foam::labelList Foam::cellShape::meshEdges
 
     forAll(localEdges, i)
     {
-        const edge& e = localEdges[i];
+        const Foam::edge& e = localEdges[i];
 
-        forAll(cEdges, j)
+        for (const label edgei : cEdges)
         {
-            label edgeI = cEdges[j];
+            if (allEdges[edgei] == e)
-
-            if (allEdges[edgeI] == e)
             {
-                modelToMesh[i] = edgeI;
+                modelToMesh[i] = edgei;
-
                 break;
             }
         }
@@ -216,6 +218,12 @@ inline Foam::labelList Foam::cellShape::meshEdges
 }
 
 
+inline Foam::face Foam::cellShape::face(const label modelFacei) const
+{
+    return m->face(modelFacei, *this);
+}
+
+
 inline Foam::faceList Foam::cellShape::faces() const
 {
     return m->faces(*this);
@@ -224,31 +232,25 @@ inline Foam::faceList Foam::cellShape::faces() const
 
 inline Foam::faceList Foam::cellShape::collapsedFaces() const
 {
-    faceList oldFaces(faces());
+    const faceList oldFaces(faces());
 
     faceList newFaces(oldFaces.size());
     label newFacei = 0;
 
-    forAll(oldFaces, oldFacei)
+    for (const Foam::face& f : oldFaces)
     {
-        const face& f = oldFaces[oldFacei];
+        Foam::face& newF = newFaces[newFacei];
-
+        newF.resize(f.size());
-        face& newF = newFaces[newFacei];
-
-        newF.setSize(f.size());
 
         label newFp = 0;
-        label prevVertI = -1;
+        label prevVerti = -1;
 
-        forAll(f, fp)
+        for (const label verti : f)
         {
-            label vertI = f[fp];
+            if (verti != prevVerti)
-
-            if (vertI != prevVertI)
             {
-                newF[newFp++] = vertI;
+                newF[newFp++] = verti;
-
+                prevVerti = verti;
-                prevVertI = vertI;
             }
         }
 
@@ -260,20 +262,19 @@ inline Foam::faceList Foam::cellShape::collapsedFaces() const
         if (newFp > 2)
         {
             // Size face and go to next one
-            newF.setSize(newFp);
+            newF.resize(newFp);
-
+            ++newFacei;
-            newFacei++;
         }
     }
-    newFaces.setSize(newFacei);
+    newFaces.resize(newFacei);
 
     return newFaces;
 }
 
 
-inline Foam::label Foam::cellShape::nFaces() const
+inline Foam::edge Foam::cellShape::edge(const label modelEdgei) const
 {
-    return m->nFaces();
+    return m->edge(modelEdgei, *this);
 }
 
 
@@ -283,18 +284,6 @@ inline Foam::edgeList Foam::cellShape::edges() const
 }
 
 
-inline Foam::label Foam::cellShape::nEdges() const
-{
-    return m->nEdges();
-}
-
-
-inline Foam::label Foam::cellShape::nPoints() const
-{
-    return size();
-}
-
-
 inline Foam::point Foam::cellShape::centre(const UList<point>& points) const
 {
     return m->centre(*this, points);

src/mesh/blockMesh/blockMesh/blockMesh.C

@@ -240,17 +240,4 @@ Foam::label Foam::blockMesh::numZonedBlocks() const
 }
 
 
-void Foam::blockMesh::writeTopology(Ostream& os) const
-{
-    for (const point& p : topology().points())
-    {
-        os << "v " << p.x() << ' ' << p.y() << ' ' << p.z() << nl;
-    }
-
-    for (const edge& e : topology().edges())
-    {
-        os << "l " << e.start() + 1 << ' ' << e.end() + 1 << nl;
-    }
-}
-
 // ************************************************************************* //

src/mesh/blockMesh/blockMesh/blockMesh.H

@@ -313,12 +313,6 @@ public:
 
         //- Create polyMesh, with cell zones
         autoPtr<polyMesh> mesh(const IOobject& io) const;
-
-
-    // Write
-
-        //- Writes edges of blockMesh in OBJ format.
-        void writeTopology(Ostream& os) const;
 };