src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersectionFuncs.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2009 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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 Description
  26
  27 \*---------------------------------------------------------------------------*/
  28
  29 #include "surfaceIntersection.H"
  30 #include "triSurfaceSearch.H"
  31 #include "labelPairLookup.H"
  32 #include "OFstream.H"
  33 #include "HashSet.H"
  34 #include "triSurface.H"
  35 #include "pointIndexHit.H"
  36 #include "octreeDataTriSurface.H"
  37 #include "octree.H"
  38 #include "meshTools.H"
  39
  40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  41
  42
  43 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  44
  45 void Foam::surfaceIntersection::writeOBJ(const point& pt, Ostream& os)
  46 {
  47     os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << endl;
  48 }
  49
  50
  51 void Foam::surfaceIntersection::writeOBJ
  52 (
  53     const List<point>& pts,
  54     const List<edge>& edges,
  55     Ostream& os
  56 )
  57 {
  58     forAll(pts, i)
  59     {
  60         writeOBJ(pts[i], os);
  61     }
  62     forAll(edges, i)
  63     {
  64         const edge& e = edges[i];
  65
  66         os << "l " << e.start()+1 << ' ' << e.end()+1 << endl;
  67     }
  68 }
  69
  70
  71 // Get minimum length of all edges connected to point
  72 Foam::scalar Foam::surfaceIntersection::minEdgeLen
  73 (
  74     const triSurface& surf,
  75     const label pointI
  76 )
  77 {
  78     const labelList& pEdges = surf.pointEdges()[pointI];
  79
  80     scalar minLen = GREAT;
  81
  82     forAll(pEdges, pEdgeI)
  83     {
  84         const edge& e = surf.edges()[pEdges[pEdgeI]];
  85
  86         minLen = min(minLen, e.mag(surf.localPoints()));
  87     }
  88
  89     return minLen;
  90 }
  91
  92
  93 // Get edge between fp and fp+1 on faceI.
  94 Foam::label Foam::surfaceIntersection::getEdge
  95 (
  96     const triSurface& surf,
  97     const label faceI,
  98     const label fp
  99 )
 100 {
 101     const labelledTri& f = surf.localFaces()[faceI];
 102
 103     edge faceEdge(f[fp], f[(fp+1) % 3]);
 104
 105     const labelList& eLabels = surf.faceEdges()[faceI];
 106
 107     forAll(eLabels, eI)
 108     {
 109         const label edgeI = eLabels[eI];
 110
 111         if (surf.edges()[edgeI] == faceEdge)
 112         {
 113             return edgeI;
 114         }
 115     }
 116
 117     FatalErrorIn
 118     (
 119         "surfaceIntersection::getEdge(const triSurface&"
 120         ", const label, const label"
 121     )   << "Problem:: Cannot find edge with vertices " << faceEdge
 122         << " in face " << faceI
 123         << abort(FatalError);
 124
 125     return -1;
 126 }
 127
 128
 129 // Given a map remove all consecutive duplicate elements.
 130 void Foam::surfaceIntersection::removeDuplicates
 131 (
 132     const labelList& map,
 133     labelList& elems
 134 )
 135 {
 136     bool hasDuplicate = false;
 137
 138     label prevVertI = -1;
 139
 140     forAll(elems, elemI)
 141     {
 142         label newVertI = map[elems[elemI]];
 143
 144         if (newVertI == prevVertI)
 145         {
 146             hasDuplicate = true;
 147
 148             break;
 149         }
 150         prevVertI = newVertI;
 151     }
 152
 153     if (hasDuplicate)
 154     {
 155         // Create copy
 156         labelList oldElems(elems);
 157
 158         label elemI = 0;
 159
 160         // Insert first
 161         elems[elemI++] = map[oldElems[0]];
 162
 163         for(label vertI = 1; vertI < oldElems.size(); vertI++)
 164         {
 165             // Insert others only if they differ from one before
 166             label newVertI = map[oldElems[vertI]];
 167
 168             if (newVertI != elems[elems.size()-1])
 169             {
 170                 elems[elemI++] = newVertI;
 171             }
 172         }
 173         elems.setSize(elemI);
 174     }
 175 }
 176
 177
 178 // Remap.
 179 void Foam::surfaceIntersection::inlineRemap
 180 (
 181     const labelList& map,
 182     labelList& elems
 183 )
 184 {
 185     forAll(elems, elemI)
 186     {
 187         elems[elemI] = map[elems[elemI]];
 188     }
 189 }
 190
 191
 192 // Remove all duplicate and degenerate elements. Return unique elements and
 193 // map from old to new.
 194 Foam::edgeList Foam::surfaceIntersection::filterEdges
 195 (
 196     const edgeList& edges,
 197     labelList& map
 198 )
 199 {
 200     HashSet<edge, Hash<edge> > uniqueEdges(10*edges.size());
 201
 202     edgeList newEdges(edges.size());
 203
 204     map.setSize(edges.size());
 205     map = -1;
 206
 207     label newEdgeI = 0;
 208
 209     forAll(edges, edgeI)
 210     {
 211         const edge& e = edges[edgeI];
 212
 213         if
 214         (
 215             (e.start() != e.end())
 216          && (uniqueEdges.find(e) == uniqueEdges.end())
 217         )
 218         {
 219             // Edge is -non degenerate and -not yet seen.
 220             uniqueEdges.insert(e);
 221
 222             map[edgeI] = newEdgeI;
 223
 224             newEdges[newEdgeI++] = e;
 225         }
 226     }
 227
 228     newEdges.setSize(newEdgeI);
 229
 230     return newEdges;
 231 }
 232
 233
 234 // Remove all duplicate elements.
 235 Foam::labelList Foam::surfaceIntersection::filterLabels
 236 (
 237     const labelList& elems,
 238     labelList& map
 239 )
 240 {
 241     labelHashSet uniqueElems(10*elems.size());
 242
 243     labelList newElems(elems.size());
 244
 245     map.setSize(elems.size());
 246     map = -1;
 247
 248     label newElemI = 0;
 249
 250     forAll(elems, elemI)
 251     {
 252         label elem = elems[elemI];
 253
 254         if (uniqueElems.find(elem) == uniqueElems.end())
 255         {
 256             // First time elem is seen
 257             uniqueElems.insert(elem);
 258
 259             map[elemI] = newElemI;
 260
 261             newElems[newElemI++] = elem;
 262         }
 263     }
 264
 265     newElems.setSize(newElemI);
 266
 267     return newElems;
 268 }
 269
 270
 271 void Foam::surfaceIntersection::writeIntersectedEdges
 272 (
 273     const triSurface& surf,
 274     const labelListList& edgeCutVerts,
 275     Ostream& os
 276 ) const
 277 {
 278     // Dump all points (surface followed by cutPoints)
 279     const pointField& pts = surf.localPoints();
 280
 281     forAll(pts, pointI)
 282     {
 283         writeOBJ(pts[pointI], os);
 284     }
 285     forAll(cutPoints(), cutPointI)
 286     {
 287         writeOBJ(cutPoints()[cutPointI], os);
 288     }
 289
 290     forAll(edgeCutVerts, edgeI)
 291     {
 292         const labelList& extraVerts = edgeCutVerts[edgeI];
 293
 294         if (extraVerts.size())
 295         {
 296             const edge& e = surf.edges()[edgeI];
 297
 298             // Start of original edge to first extra point
 299             os  << "l " << e.start()+1 << ' '
 300                 << extraVerts[0] + surf.nPoints() + 1 << endl;
 301
 302             for(label i = 1; i < extraVerts.size(); i++)
 303             {
 304                 os  << "l " << extraVerts[i-1] + surf.nPoints() + 1  << ' '
 305                     << extraVerts[i] + surf.nPoints() + 1 << endl;
 306             }
 307
 308             os  << "l " << extraVerts[extraVerts.size()-1] + surf.nPoints() + 1
 309                 << ' ' << e.end()+1 << endl;
 310         }
 311     }
 312 }
 313
 314
 315 // Return 0 (p close to start), 1(close to end) or -1.
 316 Foam::label Foam::surfaceIntersection::classify
 317 (
 318     const scalar startTol,
 319     const scalar endTol,
 320     const point& p,
 321     const edge& e,
 322     const pointField& points
 323 )
 324 {
 325     if (mag(p - points[e.start()]) < startTol)
 326     {
 327         return 0;
 328     }
 329     else if (mag(p - points[e.end()]) < endTol)
 330     {
 331         return 1;
 332     }
 333     else
 334     {
 335         return -1;
 336     }
 337 }
 338
 339
 340 // ************************************************************************* //