applications/utilities/surface/surfaceAdd/surfaceAdd.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   6      \\/     M anipulation  |
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of OpenFOAM.
  10
  11     OpenFOAM is free software: you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by
  13     the Free Software Foundation, either version 3 of the License, or
  14     (at your option) any later version.
  15
  16     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
  17     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  18     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  19     for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Description
  25     Add two surfaces. Does geometric merge on points. Does not check for
  26     overlapping/intersecting triangles.
  27
  28     Keeps patches separate by renumbering.
  29
  30 \*---------------------------------------------------------------------------*/
  31
  32 #include "argList.H"
  33 #include "fileName.H"
  34 #include "triSurface.H"
  35 #include "OFstream.H"
  36 #include "IFstream.H"
  37 #include "triFace.H"
  38 #include "triFaceList.H"
  39
  40 using namespace Foam;
  41
  42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  43
  44 // Main program:
  45
  46 int main(int argc, char *argv[])
  47 {
  48     argList::addNote
  49     (
  50         "add two surfaces via a geometric merge on points."
  51     );
  52
  53     argList::noParallel();
  54     argList::validArgs.append("surfaceFile");
  55     argList::validArgs.append("surfaceFile");
  56     argList::validArgs.append("output surfaceFile");
  57
  58     argList::addOption
  59     (
  60         "points",
  61         "file",
  62         "provide additional points"
  63     );
  64     argList::addBoolOption
  65     (
  66         "mergeRegions",
  67         "combine regions from both surfaces"
  68     );
  69
  70     argList args(argc, argv);
  71
  72     const fileName inFileName1 = args[1];
  73     const fileName inFileName2 = args[2];
  74     const fileName outFileName = args[3];
  75
  76     const bool addPoint     = args.optionFound("points");
  77     const bool mergeRegions = args.optionFound("mergeRegions");
  78
  79     if (addPoint)
  80     {
  81         Info<< "Reading a surface and adding points from a file"
  82             << "; merging the points and writing the surface to another file"
  83             << nl << endl;
  84
  85         Info<< "Surface  : " << inFileName1<< nl
  86             << "Points   : " << args["points"] << nl
  87             << "Writing  : " << outFileName << nl << endl;
  88     }
  89     else
  90     {
  91         Info<< "Reading two surfaces"
  92             << "; merging points and writing the surface to another file"
  93             << nl << endl;
  94
  95         if (mergeRegions)
  96         {
  97             Info<< "Regions from the two files will get merged" << nl
  98                 << "Do not use this option if you want to keep the regions"
  99                 << " separate" << nl << endl;
 100         }
 101         else
 102         {
 103             Info<< "Regions from the two files will not get merged" << nl
 104                 << "Regions from " << inFileName2 << " will get offset so"
 105                 << " as not to overlap with the regions in " << inFileName1
 106                 << nl << endl;
 107         }
 108
 109
 110         Info<< "Surface1 : " << inFileName1<< nl
 111             << "Surface2 : " << inFileName2<< nl
 112             << "Writing  : " << outFileName << nl << endl;
 113     }
 114
 115     const triSurface surface1(inFileName1);
 116
 117     Info<< "Surface1:" << endl;
 118     surface1.writeStats(Info);
 119     Info<< endl;
 120
 121     const pointField& points1 = surface1.points();
 122
 123     // Final surface
 124     triSurface combinedSurf;
 125
 126     if (addPoint)
 127     {
 128         IFstream pointsFile(args["points"]);
 129         pointField extraPoints(pointsFile);
 130
 131         Info<< "Additional Points:" << extraPoints.size() << endl;
 132
 133         vectorField pointsAll(points1);
 134         label pointI = pointsAll.size();
 135         pointsAll.setSize(pointsAll.size() + extraPoints.size());
 136
 137         forAll(extraPoints, i)
 138         {
 139             pointsAll[pointI++] = extraPoints[i];
 140         }
 141
 142         combinedSurf = triSurface(surface1, surface1.patches(), pointsAll);
 143     }
 144     else
 145     {
 146         const triSurface surface2(inFileName2);
 147
 148         Info<< "Surface2:" << endl;
 149         surface2.writeStats(Info);
 150         Info<< endl;
 151
 152
 153         // Make new storage
 154         List<labelledTri> facesAll(surface1.size() + surface2.size());
 155
 156         const pointField& points2 = surface2.points();
 157
 158         vectorField pointsAll(points1.size() + points2.size());
 159
 160
 161         label pointi = 0;
 162         // Copy points1 into pointsAll
 163         forAll(points1, point1i)
 164         {
 165             pointsAll[pointi++] = points1[point1i];
 166         }
 167         // Add surface2 points
 168         forAll(points2, point2i)
 169         {
 170             pointsAll[pointi++] = points2[point2i];
 171         }
 172
 173
 174         label trianglei = 0;
 175
 176         // Copy triangles1 into trianglesAll
 177         forAll(surface1, faceI)
 178         {
 179             facesAll[trianglei++] = surface1[faceI];
 180         }
 181         label nRegions1 = surface1.patches().size();
 182
 183
 184         if (!mergeRegions)
 185         {
 186             Info<< "Surface " << inFileName1 << " has " << nRegions1
 187                 << " regions"
 188                 << nl
 189                 << "All region numbers in " << inFileName2 << " will be offset"
 190                 << " by this amount" << nl << endl;
 191         }
 192
 193         // Add (renumbered) surface2 triangles
 194         forAll(surface2, faceI)
 195         {
 196             const labelledTri& tri = surface2[faceI];
 197
 198             labelledTri& destTri = facesAll[trianglei++];
 199             destTri[0] = tri[0] + points1.size();
 200             destTri[1] = tri[1] + points1.size();
 201             destTri[2] = tri[2] + points1.size();
 202             if (mergeRegions)
 203             {
 204                 destTri.region() = tri.region();
 205             }
 206             else
 207             {
 208                 destTri.region() = tri.region() + nRegions1;
 209             }
 210         }
 211
 212         label nRegions2 = surface2.patches().size();
 213
 214         geometricSurfacePatchList newPatches;
 215
 216         if (mergeRegions)
 217         {
 218             // Overwrite
 219             newPatches.setSize(max(nRegions1, nRegions2));
 220
 221             forAll(surface1.patches(), patchI)
 222             {
 223                 newPatches[patchI] = surface1.patches()[patchI];
 224             }
 225             forAll(surface2.patches(), patchI)
 226             {
 227                 newPatches[patchI] = surface2.patches()[patchI];
 228             }
 229         }
 230         else
 231         {
 232             Info<< "Regions from " << inFileName2 << " have been renumbered:"
 233                 << nl
 234                 << "    old\tnew" << nl;
 235
 236             for (label regionI = 0; regionI < nRegions2; regionI++)
 237             {
 238                 Info<< "    " << regionI << '\t' << regionI+nRegions1
 239                     << nl;
 240             }
 241             Info<< nl;
 242
 243             newPatches.setSize(nRegions1 + nRegions2);
 244
 245             label newPatchI = 0;
 246
 247             forAll(surface1.patches(), patchI)
 248             {
 249                 newPatches[newPatchI++] = surface1.patches()[patchI];
 250             }
 251
 252             forAll(surface2.patches(), patchI)
 253             {
 254                 newPatches[newPatchI++] = surface2.patches()[patchI];
 255             }
 256         }
 257
 258
 259         Info<< "New patches:" << nl;
 260         forAll(newPatches, patchI)
 261         {
 262             Info<< "    " << patchI << '\t' << newPatches[patchI].name() << nl;
 263         }
 264         Info<< endl;
 265
 266
 267         // Construct new surface mesh
 268         combinedSurf = triSurface(facesAll, newPatches, pointsAll);
 269     }
 270
 271     // Merge all common points and do some checks
 272     combinedSurf.cleanup(true);
 273
 274     Info<< "Merged surface:" << endl;
 275
 276     combinedSurf.writeStats(Info);
 277
 278     Info<< endl;
 279
 280     Info<< "Writing : " << outFileName << endl;
 281
 282     // No need to 'group' while writing since all in correct order anyway.
 283     combinedSurf.write(outFileName);
 284
 285     Info<< "End\n" << endl;
 286
 287     return 0;
 288 }
 289
 290
 291 // ************************************************************************* //