src/meshTools/triSurface/booleanOps/intersectedSurface/edgeSurface.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
   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 \*---------------------------------------------------------------------------*/
  25
  26 #include "edgeSurface.H"
  27 #include "triSurface.H"
  28 #include "surfaceIntersection.H"
  29 #include "meshTools.H"
  30 #include "OFstream.H"
  31
  32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  33
  34 defineTypeNameAndDebug(Foam::edgeSurface, 0);
  35
  36
  37 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  38
  39 // Write whole pointField and edges to stream
  40 void Foam::edgeSurface::writeOBJ
  41 (
  42     const pointField& points,
  43     const edgeList& edges,
  44     Ostream& os
  45 )
  46 {
  47     forAll(points, pointI)
  48     {
  49         const point& pt = points[pointI];
  50
  51         os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << endl;
  52     }
  53     forAll(edges, edgeI)
  54     {
  55         const edge& e = edges[edgeI];
  56
  57         os << "l " << e.start()+1 << ' ' << e.end()+1 << endl;
  58     }
  59 }
  60
  61
  62 // Write whole pointField and selected edges to stream
  63 void Foam::edgeSurface::writeOBJ
  64 (
  65     const pointField& points,
  66     const edgeList& edges,
  67     const labelList& edgeLabels,
  68     Ostream& os
  69 )
  70 {
  71     forAll(points, pointI)
  72     {
  73         const point& pt = points[pointI];
  74
  75         os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << endl;
  76     }
  77     forAll(edgeLabels, i)
  78     {
  79         const edge& e = edges[edgeLabels[i]];
  80
  81         os << "l " << e.start()+1 << ' ' << e.end()+1 << endl;
  82     }
  83 }
  84
  85
  86 // Pointedges in edgeSurface indices only.
  87 void Foam::edgeSurface::calcPointEdges()
  88 {
  89     pointEdges_.setSize(points_.size());
  90
  91     labelList pointNEdges(points_.size(), 0);
  92
  93     forAll(edges_, edgeI)
  94     {
  95         const edge& e = edges_[edgeI];
  96
  97         pointNEdges[e[0]]++;
  98         pointNEdges[e[1]]++;
  99     }
 100
 101     forAll(pointEdges_, pointI)
 102     {
 103         pointEdges_[pointI].setSize(pointNEdges[pointI]);
 104     }
 105
 106     pointNEdges = 0;
 107
 108     forAll(edges_, edgeI)
 109     {
 110         const edge& e = edges_[edgeI];
 111
 112         labelList& pEdges0 = pointEdges_[e[0]];
 113         pEdges0[pointNEdges[e[0]]++] = edgeI;
 114
 115         labelList& pEdges1 = pointEdges_[e[1]];
 116         pEdges1[pointNEdges[e[1]]++] = edgeI;
 117     }
 118 }
 119
 120
 121 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 122
 123 // Construct from surface and intersection description
 124 Foam::edgeSurface::edgeSurface
 125 (
 126     const triSurface& surf,
 127     const bool isFirstSurface,
 128     const surfaceIntersection& inter
 129 )
 130 :
 131     points_(surf.nPoints() + inter.cutPoints().size()),
 132     nSurfacePoints_(surf.nPoints()),
 133     edges_(),
 134     nSurfaceEdges_(surf.nEdges()),
 135     parentEdges_(0),
 136     faceEdges_(surf.size()),
 137     pointEdges_(points_.size())
 138 {
 139     // Copy points (surface ones first)
 140     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 141
 142     label pointI = 0;
 143
 144     const pointField& surfPoints = surf.localPoints();
 145
 146     forAll(surfPoints, i)
 147     {
 148         points_[pointI++] = surfPoints[i];
 149     }
 150
 151     const pointField& cutPoints = inter.cutPoints();
 152
 153     forAll(cutPoints, i)
 154     {
 155         points_[pointI++] = cutPoints[i];
 156     }
 157
 158
 159     // Copy edges (surface ones first)
 160     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 161
 162     DynamicList<edge> allEdges(surf.nEdges() + inter.cutEdges().size());
 163     DynamicList<label> allParentEdges(surf.nEdges());
 164     List<DynamicList<label> > allFaceEdges(surf.size());
 165
 166
 167     // Copy surface edges (can be split!)
 168
 169     const edgeList& surfEdges = surf.edges();
 170
 171     forAll(surfEdges, edgeI)
 172     {
 173         const edge& e = surfEdges[edgeI];
 174
 175         // Get additional vertices for this edge.
 176         const labelList& extraVerts = inter.edgeCuts(isFirstSurface)[edgeI];
 177
 178         // Store current top of allEdges.
 179         label freeNewEdgeI = allEdges.size();
 180
 181         if (extraVerts.empty())
 182         {
 183             // No cuts across this edge. Note that vertices do not need to be
 184             // renumbered.
 185             allEdges.append(e);
 186         }
 187         else
 188         {
 189             // Edge is cut. From e.start() to extraVerts[0],
 190             // from extraVerts[i] to i+1 and finally to e.end().
 191             allEdges.append
 192             (
 193                 edge
 194                 (
 195                     e.start(),
 196                     extraVerts[0] + nSurfacePoints_
 197                 )
 198             );
 199
 200             for (label extraI = 1; extraI < extraVerts.size(); extraI++)
 201             {
 202                 allEdges.append
 203                 (
 204                     edge
 205                     (
 206                         extraVerts[extraI-1] + nSurfacePoints_,
 207                         extraVerts[extraI] + nSurfacePoints_
 208                     )
 209                 );
 210             }
 211             allEdges.append
 212             (
 213                 edge
 214                 (
 215                     extraVerts.last() + nSurfacePoints_,
 216                     e.end()
 217                 )
 218             );
 219         }
 220
 221         // Update allFaceEdges, parentEdges_ for the newly added edges.
 222
 223         // Add each edge label to all face neighbours of edgeI
 224         const labelList& myFaces = surf.edgeFaces()[edgeI];
 225
 226         for (label eI = freeNewEdgeI; eI < allEdges.size(); eI++)
 227         {
 228             allParentEdges.append(edgeI);
 229
 230             forAll(myFaces, myFaceI)
 231             {
 232                 allFaceEdges[myFaces[myFaceI]].append(eI);
 233             }
 234         }
 235     }
 236
 237     // Done all (possibly split) surface edges by now.
 238     nSurfaceEdges_ = allEdges.size();
 239
 240
 241     // Copy intersection edges
 242     // (note no parentEdges)
 243     const edgeList& cutEdges = inter.cutEdges();
 244
 245     forAll(cutEdges, i)
 246     {
 247         const edge& e = cutEdges[i];
 248
 249         allEdges.append(edge(e[0] + nSurfacePoints_, e[1] + nSurfacePoints_));
 250     }
 251
 252
 253
 254
 255     // Add intersection edges to faceEdges
 256     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 257
 258     forAllConstIter(labelPairLookup, inter.facePairToEdge(), iter)
 259     {
 260         // Edge label in intersection
 261         const label edgeI = iter();
 262
 263         // Get the face from the correct surface
 264         const FixedList<label, 2>& twoFaces = iter.key();
 265
 266         label faceI;
 267
 268         if (isFirstSurface)
 269         {
 270             faceI = twoFaces[0];
 271         }
 272         else
 273         {
 274             faceI = twoFaces[1];
 275         }
 276
 277         // Store on face-edge addressing. (note: offset edge)
 278         allFaceEdges[faceI].append(edgeI + nSurfaceEdges_);
 279     }
 280
 281     // Transfer.
 282     edges_.transfer(allEdges);
 283     parentEdges_.transfer(allParentEdges);
 284
 285     forAll(allFaceEdges, faceI)
 286     {
 287         faceEdges_[faceI].transfer(allFaceEdges[faceI]);
 288     }
 289
 290
 291     // Additional addressing
 292     // ~~~~~~~~~~~~~~~~~~~~~
 293
 294     calcPointEdges();
 295
 296
 297     if (debug & 4)
 298     {
 299         Pout<< "edgeSurface : Dumping faceEdges to files" << endl;
 300
 301         forAll(faceEdges_, faceI)
 302         {
 303             const labelList& fEdges = faceEdges_[faceI];
 304
 305             if (fEdges.size() != 3)
 306             {
 307                 fileName faceFName("face_" + name(faceI) + ".obj");
 308                 Pout<< "edgeSurface : Dumping faceEdges for face " << faceI
 309                     << " to " << faceFName << endl;
 310
 311                 OFstream fStream(faceFName);
 312                 writeOBJ(points_, edges_, fEdges, fStream);
 313             }
 314         }
 315
 316         Pout<< "edgeSurface : Dumping edges to edges.obj" << endl;
 317         OFstream eStream("edges.obj");
 318         writeOBJ(points_, edges_, eStream);
 319
 320         Pout<< "edgeSurface : Dumping intersectionEdges to"
 321             << " intersectionEdges.obj" << endl;
 322         OFstream intEdgesStream("intersectionEdges.obj");
 323
 324         labelList edgeLabels(edges_.size() - nSurfaceEdges_);
 325
 326         label i = 0;
 327         for (label edgeI = nSurfaceEdges_; edgeI < edges_.size(); edgeI++)
 328         {
 329             edgeLabels[i++] = edgeI;
 330         }
 331
 332         writeOBJ(points_, edges_, edgeLabels, intEdgesStream);
 333     }
 334 }
 335
 336
 337 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 338
 339 void Foam::edgeSurface::addIntersectionEdges
 340 (
 341     const label faceI,
 342     const edgeList& additionalEdges
 343 )
 344 {
 345     if (debug & 2)
 346     {
 347         Pout<< "Old face consisted of edges:" << endl;
 348
 349         const labelList& fEdges = faceEdges_[faceI];
 350         forAll(fEdges, i)
 351         {
 352             const edge& e = edges_[fEdges[i]];
 353
 354             Pout<< "    " << fEdges[i] << ' ' << e
 355                 << points_[e.start()] << ' ' << points_[e.end()] << endl;
 356         }
 357     }
 358
 359     // Make space for additional intersection edges (copies old ones)
 360     const label oldNEdges = edges_.size();
 361
 362     edges_.setSize(oldNEdges + additionalEdges.size());
 363
 364     // Append new intersection edges
 365     label newEdgeI = oldNEdges;
 366
 367     forAll(additionalEdges, i)
 368     {
 369         edges_[newEdgeI] = additionalEdges[i];  // Vertices already in eSurf
 370                                                 // indices.
 371         newEdgeI++;
 372     }
 373
 374     // Append to faceEdges.
 375     labelList& fEdges = faceEdges_[faceI];
 376
 377     label nFEdges = fEdges.size();
 378
 379     fEdges.setSize(nFEdges + additionalEdges.size());
 380
 381     forAll(additionalEdges, i)
 382     {
 383         fEdges[nFEdges++] = oldNEdges + i;
 384     }
 385
 386
 387     // Update pointEdge addressing
 388     calcPointEdges();
 389
 390
 391     if (debug & 2)
 392     {
 393         const labelList& fEdges = faceEdges_[faceI];
 394
 395         Pout<< "New face consists of edges:" << endl;
 396         forAll(fEdges, i)
 397         {
 398             const edge& e = edges_[fEdges[i]];
 399
 400             Pout<< "    " << fEdges[i] << ' ' << e
 401                 << points_[e.start()] << ' ' << points_[e.end()] << endl;
 402         }
 403     }
 404 }
 405
 406
 407 // ************************************************************************* //