ENH: improve surfaceSplitByPatch controls (#1600)

- uses MeshedSurface instead of triSurface to prevent automatic
  triangulation.

- supports '-patches' and '-excludePatches' controls as per foamToVTK.
  For example,

    surfaceSplitByPatch -patches '( ".*rider.*" )'  motorBike.obj

ENH: use MeshedSurface for surfaceSubset

applications/utilities/surface/surfaceSplitByPatch/surfaceSplitByPatch.C

@@ -31,12 +31,33 @@ Group
     grpSurfaceUtilities
 
 Description
-    Writes regions of triSurface to separate files.
+    Writes surface regions to separate files.
+
+Usage
+    \b surfaceSplitByPatch [OPTION]
+
+    Options:
+      - \par -patches NAME | LIST
+        Specify single patch or multiple patches (name or regex) to extract
+        For example,
+        \verbatim
+          -patches top
+          -patches '( front \".*back\" )'
+        \endverbatim
+
+      - \par -excludePatches NAME | LIST
+        Specify single or multiple patches (name or regex) not to extract.
+        For example,
+        \verbatim
+          -excludePatches '( inlet_1 inlet_2 "proc.*")'
+        \endverbatim
 
 \*---------------------------------------------------------------------------*/
 
 #include "argList.H"
-#include "triSurface.H"
+#include "MeshedSurfaces.H"
+#include "stringListOps.H"
+#include "geometricSurfacePatch.H"
 
 using namespace Foam;
 
@@ -48,65 +69,103 @@ int main(int argc, char *argv[])
     (
         "Write surface mesh regions to separate files"
     );
-
     argList::noParallel();
+    argList::addOption
+    (
+        "patches",
+        "wordRes",
+        "Specify single patch or multiple patches to write\n"
+        "Eg, 'top' or '( front \".*back\" )'"
+    );
+    argList::addOption
+    (
+        "excludePatches",
+        "wordRes",
+        "Specify single patch or multiple patches to exclude from writing."
+        " Eg, 'outlet' or '( inlet \".*Wall\" )'"
+    );
+
     argList::addArgument("input", "The input surface file");
+
     argList args(argc, argv);
 
     const fileName surfName = args[1];
 
-    Info<< "Reading surf from " << surfName << " ..." << nl << endl;
+    const fileName surfBase(surfName.lessExt());
 
-    fileName surfBase = surfName.lessExt();
+    const word extension(surfName.ext());
-
-    word extension = surfName.ext();
-
-    triSurface surf(surfName);
-
-    Info<< "Writing regions to separate files ..."
-        << nl << endl;
 
 
-    const geometricSurfacePatchList& patches = surf.patches();
+    Info<< nl
+        << "Read surface from " << surfName << " ..." << nl << endl;
 
-    forAll(patches, patchi)
+    meshedSurface surf(surfName);
+
+    const surfZoneList& zones = surf.surfZones();
+
+    Info<< "    " << surf.size() << " faces with "
+        << zones.size() << " zones" << nl << nl;
+
+
+    wordRes includePatches, excludePatches;
+
+    if (args.readListIfPresent<wordRe>("patches", includePatches))
     {
-        const geometricSurfacePatch& pp = patches[patchi];
+        Info<< "Including patches " << flatOutput(includePatches)
+            << nl << endl;
+    }
+    if (args.readListIfPresent<wordRe>("excludePatches", excludePatches))
+    {
+        Info<< "Excluding patches " << flatOutput(excludePatches)
+            << nl << endl;
+    }
 
-        word patchName(pp.name());
+    // Identity if both whitelist and blacklist are empty
+    const labelList zoneIndices
+    (
+        stringListOps::findMatching
+        (
+            zones,
+            includePatches,
+            excludePatches,
+            nameOp<surfZone>()
+        )
+    );
+
+
+    Info<< "Writing regions to "
+        << zoneIndices.size() << " separate files ..." << nl << endl;
+
+
+    // Faces to subset
+    bitSet includeMap(surf.size());
+
+    for (const label zonei : zoneIndices)
+    {
+        const surfZone& zn = zones[zonei];
+
+        includeMap.reset();
+        includeMap.set(zn.range());
+
+        word patchName(zn.name());
 
         if (patchName.empty())
         {
-            patchName = geometricSurfacePatch::defaultName(patchi);
+            // In case people expect the same names as with triSurface
+            patchName = geometricSurfacePatch::defaultName(zonei);
         }
 
         fileName outFile(surfBase + '_' + patchName + '.' + extension);
 
-        Info<< "   Writing patch " << patchName << " to file " << outFile
+        Info<< "   Zone " << zonei << " (" << zn.size() << " faces) "
-            << endl;
+            << patchName
+            << " to file " << outFile << nl;
 
-
+        // Subset and write
-        // Collect faces of region
+        surf.subsetMesh(includeMap).write(outFile);
-        boolList includeMap(surf.size(), false);
-
-        forAll(surf, facei)
-        {
-            const labelledTri& f = surf[facei];
-
-            if (f.region() == patchi)
-            {
-                includeMap[facei] = true;
-            }
     }
 
-        // Subset triSurface
+    Info<< "\nEnd\n" << endl;
-
-        triSurface subSurf(surf.subsetMesh(includeMap));
-
-        subSurf.write(outFile);
-    }
-
-    Info<< "End\n" << endl;
 
     return 0;
 }

applications/utilities/surface/surfaceSubset/surfaceSubset.C

@@ -6,7 +6,7 @@
      \\/     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.
@@ -36,8 +36,8 @@ Description
 
 \*---------------------------------------------------------------------------*/
 
-#include "triSurface.H"
 #include "triSurfaceSearch.H"
+#include "MeshedSurfaces.H"
 #include "argList.H"
 #include "Fstream.H"
 #include "IOdictionary.H"
@@ -56,7 +56,7 @@ int main(int argc, char *argv[])
 {
     argList::addNote
     (
-        "A surface analysis tool that subsets the triSurface to choose a"
+        "A surface analysis tool that subsets the surface to choose a"
         " region of interest."
     );
 
@@ -71,7 +71,8 @@ int main(int argc, char *argv[])
     dictionary meshSubsetDict(dictFile);
 
     Info<< "Reading surface " << args[2] << " ..." << endl;
-    triSurface surf1(args[2]);
+
+    meshedSurface surf1(args[2]);
 
     const fileName outFileName(args[3]);
 
@@ -101,7 +102,12 @@ int main(int argc, char *argv[])
         meshSubsetDict.lookup("zone")
     );
 
-    if (markedZone.size() && markedZone.size() != 2)
+
+    boundBox zoneBb;
+
+    if (markedZone.size())
+    {
+        if (markedZone.size() != 2)
         {
             FatalErrorInFunction
                 << "zone specification should be two points, min and max of "
@@ -110,125 +116,101 @@ int main(int argc, char *argv[])
                 << exit(FatalError);
         }
 
+        zoneBb.min() = markedZone[0];
+        zoneBb.max() = markedZone[1];
+
+        if (!zoneBb.valid())
+        {
+            WarningInFunction
+                << "Defined zone is invalid: " << zoneBb << nl;
+        }
+    }
+
+
     const bool addFaceNeighbours =
         meshSubsetDict.get<bool>("addFaceNeighbours");
 
     const bool invertSelection =
-        meshSubsetDict.lookupOrDefault("invertSelection", false);
+        meshSubsetDict.getOrDefault("invertSelection", false);
 
     // Mark the cells for the subset
 
     // Faces to subset
-    boolList facesToSubset(surf1.size(), false);
+    bitSet facesToSubset(surf1.size(), false);
 
 
     //
-    // pick up faces connected to "localPoints"
+    // Faces connected to "localPoints"
     //
 
     if (markedPoints.size())
     {
         Info<< "Found " << markedPoints.size() << " marked point(s)." << endl;
 
-        // pick up cells sharing the point
+        for (const label pointi : markedPoints)
-
-        forAll(markedPoints, pointi)
         {
-            if
+            if (pointi < 0 || pointi >= surf1.nPoints())
-            (
-                markedPoints[pointi] < 0
-             || markedPoints[pointi] >= surf1.nPoints()
-            )
             {
                 FatalErrorInFunction
-                    << "localPoint label " << markedPoints[pointi]
+                    << "localPoint label " << pointi << "out of range."
-                    << "out of range."
+                    << " Surface has " << surf1.nPoints() << " localPoints."
-                    << " The mesh has got "
-                    << surf1.nPoints() << " localPoints."
                     << exit(FatalError);
             }
 
-            const labelList& curFaces =
+            const labelList& curFaces = surf1.pointFaces()[pointi];
-                surf1.pointFaces()[markedPoints[pointi]];
 
-            forAll(curFaces, i)
+            facesToSubset.set(curFaces);
-            {
-                facesToSubset[curFaces[i]] =  true;
         }
     }
-    }
-
 
 
     //
-    // pick up faces connected to "edges"
+    // Faces connected to "edges"
     //
 
     if (markedEdges.size())
     {
         Info<< "Found " << markedEdges.size() << " marked edge(s)." << endl;
 
-        // pick up cells sharing the edge
+        for (const label edgei : markedEdges)
-
-        forAll(markedEdges, edgeI)
         {
-            if
+            if (edgei < 0 || edgei >= surf1.nEdges())
-            (
-                markedEdges[edgeI] < 0
-             || markedEdges[edgeI] >= surf1.nEdges()
-            )
             {
                 FatalErrorInFunction
-                    << "edge label " << markedEdges[edgeI]
+                    << "edge label " << edgei << "out of range."
-                    << "out of range."
+                    << " Surface has " << surf1.nEdges() << " edges."
-                    << " The mesh has got "
-                    << surf1.nEdges() << " edges."
                     << exit(FatalError);
             }
 
-            const labelList& curFaces = surf1.edgeFaces()[markedEdges[edgeI]];
+            const labelList& curFaces = surf1.edgeFaces()[edgei];
 
-            forAll(curFaces, i)
+            facesToSubset.set(curFaces);
-            {
-                facesToSubset[curFaces[i]] =  true;
-            }
         }
     }
 
 
     //
-    // pick up faces with centre inside "zone"
+    // Faces with centre inside "zone"
     //
 
-    if (markedZone.size() == 2)
+    if (zoneBb.valid())
     {
-        const point& min = markedZone[0];
+        Info<< "Using zone " << zoneBb << endl;
-        const point& max = markedZone[1];
-
-        Info<< "Using zone min:" << min << " max:" << max << endl;
 
         forAll(surf1, facei)
         {
             const point centre = surf1[facei].centre(surf1.points());
 
-            if
+            if (zoneBb.contains(centre))
-            (
-                (centre.x() >= min.x())
-             && (centre.y() >= min.y())
-             && (centre.z() >= min.z())
-             && (centre.x() <= max.x())
-             && (centre.y() <= max.y())
-             && (centre.z() <= max.z())
-            )
             {
-                facesToSubset[facei] = true;
+                facesToSubset.set(facei);
             }
         }
     }
 
 
     //
-    // pick up faces on certain side of surface
+    // Faces on certain side of surface
     //
 
     if (meshSubsetDict.found("surface"))
@@ -237,18 +219,16 @@ int main(int argc, char *argv[])
 
         const fileName surfName(surfDict.get<fileName>("name"));
 
-        const bool outside(surfDict.get<bool>("outside"));
+        const volumeType::type volType =
+        (
+            surfDict.getOrDefault("outside", false)
+          ? volumeType::OUTSIDE
+          : volumeType::INSIDE
+        );
 
-        if (outside)
+        Info<< "Selecting faces with centre located "
-        {
+            << volumeType::names[volType] << " of surface "
-            Info<< "Selecting all triangles with centre outside surface "
             << surfName << endl;
-        }
-        else
-        {
-            Info<< "Selecting all triangles with centre inside surface "
-                << surfName << endl;
-        }
 
         // Read surface to select on
         triSurface selectSurf(surfName);
@@ -264,22 +244,15 @@ int main(int argc, char *argv[])
             searchSelectSurf.tree();
 
         // Check if face (centre) is in outside or inside.
-        forAll(facesToSubset, facei)
+        forAll(surf1, facei)
         {
             if (!facesToSubset[facei])
             {
                 const point fc(surf1[facei].centre(surf1.points()));
 
-                volumeType t = selectTree.getVolumeType(fc);
+                if (volType == selectTree.getVolumeType(fc))
-
-                if
-                (
-                    outside
-                  ? (t == volumeType::OUTSIDE)
-                  : (t == volumeType::INSIDE)
-                )
                 {
-                    facesToSubset[facei] = true;
+                    facesToSubset.set(facei);
                 }
             }
         }
@@ -304,7 +277,7 @@ int main(int argc, char *argv[])
 
             if (pl.distance(fc) < distance && mag(pl.normal() & nf) > cosAngle)
             {
-                facesToSubset[facei] = true;
+                facesToSubset.set(facei);
             }
         }
     }
@@ -323,38 +296,29 @@ int main(int argc, char *argv[])
         Info<< "Found " << markedFaces.size() << " marked face(s)." << endl;
 
         // Check and mark faces to pick up
-        forAll(markedFaces, facei)
+        for (const label facei : markedFaces)
         {
-            if
+            if (facei < 0 || facei >= surf1.size())
-            (
-                markedFaces[facei] < 0
-             || markedFaces[facei] >= surf1.size()
-            )
             {
                 FatalErrorInFunction
-                    << "Face label " << markedFaces[facei] << "out of range."
+                    << "Face label " << facei << "out of range."
-                    << " The mesh has got "
+                    << " Surface has " << surf1.size() << " faces."
-                    << surf1.size() << " faces."
                     << exit(FatalError);
             }
 
             // Mark the face
-            facesToSubset[markedFaces[facei]] = true;
+            facesToSubset.set(facei);
 
-            // mark its neighbours if requested
+            // Mark its neighbours if requested
             if (addFaceNeighbours)
             {
-                const labelList& curFaces =
+                const labelList& curFaces = surf1.faceFaces()[facei];
-                    surf1.faceFaces()[markedFaces[facei]];
 
-                forAll(curFaces, i)
+                for (const label neiFacei : curFaces)
                 {
-                    label facei = curFaces[i];
+                    if (facesToSubset.set(neiFacei))
-
-                    if (!facesToSubset[facei])
                     {
-                        facesToSubset[facei] =  true;
+                        ++nFaceNeighbours;
-                        nFaceNeighbours++;
                     }
                 }
             }
@@ -372,15 +336,12 @@ int main(int argc, char *argv[])
     {
         Info<< "Inverting selection." << endl;
 
-        forAll(facesToSubset, i)
+        facesToSubset.flip();
-        {
-            facesToSubset[i] = facesToSubset[i] ? false : true;
-        }
     }
 
 
     // Create subsetted surface
-    triSurface surf2(surf1.subsetMesh(facesToSubset));
+    meshedSurface surf2(surf1.subsetMesh(facesToSubset));
 
     Info<< "Subset:" << endl;
     surf2.writeStats(Info);