triSurfaceTools.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 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 Class
25  Foam::triSurfaceTools
26 
27 Description
28  A collection of tools for triSurface.
29 
30 SourceFiles
31  triSurfaceTools.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef triSurfaceTools_H
36 #define triSurfaceTools_H
37 
38 #include "boolList.H"
39 #include "pointField.H"
40 #include "DynamicList.H"
41 #include "HashSet.H"
42 #include "FixedList.H"
43 #include "vector2D.H"
44 #include "triPointRef.H"
45 #include "surfaceLocation.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 // Forward declaration of classes
53 class triSurface;
54 class edge;
55 class labelledTri;
56 class polyBoundaryMesh;
57 class plane;
58 
59 /*---------------------------------------------------------------------------*\
60  Class triSurfaceTools Declaration
61 \*---------------------------------------------------------------------------*/
62 
63 class triSurfaceTools
64 {
65  // Private Member Functions
66 
67  // Refinement
68 
69  enum refineType
70  {
71  NONE,
72  RED,
73  GREEN
74  };
75  static void calcRefineStatus
76  (
77  const triSurface& surf,
78  const label facei,
79  List<refineType>& refine
80  );
81  static void greenRefine
82  (
83  const triSurface& surf,
84  const label facei,
85  const label edgeI,
86  const label newPointi,
87  DynamicList<labelledTri>& newFaces
88  );
89  static triSurface doRefine
90  (
91  const triSurface& surf,
92  const List<refineType>& refineStatus
93  );
94 
95 
96  // Coarsening
97 
98  static scalar faceCosAngle
99  (
100  const point& pStart,
101  const point& pEnd,
102  const point& pLeft,
103  const point& pRight
104  );
105 
106  static void protectNeighbours
107  (
108  const triSurface& surf,
109  const label vertI,
110  labelList& faceStatus
111  );
112 
113  //- Faces to collapse because of edge collapse
114  static labelHashSet getCollapsedFaces
115  (
116  const triSurface& surf,
117  label edgeI
118  );
119 
120  // Return value of faceUsed for faces using vertI (local numbering).
121  // Used internally.
122  static label vertexUsesFace
123  (
124  const triSurface& surf,
125  const labelHashSet& faceUsed,
126  const label vertI
127  );
128 
129  // Get new connections between faces (because of edge collapse)
130  // in form of tables:
131  // - given edge get other edge
132  // - given edge get other face
133  // A face using point v1 on edge will get connected to a face using
134  // point v2 if they share a common vertex
135  // (but not a common edge since then the triangles collapse to
136  // nothing)
137  static void getMergedEdges
138  (
139  const triSurface& surf,
140  const label edgeI,
141  const labelHashSet& collapsedFaces,
142  HashTable<label, label, Hash<label>>& edgeToEdge,
143  HashTable<label, label, Hash<label>>& edgeToFace
144  );
145 
146  //- Calculates (cos of) angle across edgeI of facei,
147  // taking into account updated addressing (resulting from edge
148  // collapse)
149  static scalar edgeCosAngle
150  (
151  const triSurface& surf,
152  const label v1,
153  const point& pt,
154  const labelHashSet& collapsedFaces,
155  const HashTable<label, label, Hash<label>>& edgeToEdge,
156  const HashTable<label, label, Hash<label>>& edgeToFace,
157  const label facei,
158  const label edgeI
159  );
160 
161  //- Calculate minimum (cos of) edge angle using addressing from
163  // edge to v1 at pt. Returns 1 if v1 is on edge without neighbours
164  // (and hence no edge angle can be defined)
165  static scalar collapseMinCosAngle
166  (
167  const triSurface& surf,
168  const label v1,
169  const point& pt,
170  const labelHashSet& collapsedFaces,
171  const HashTable<label, label, Hash<label>>& edgeToEdge,
172  const HashTable<label, label, Hash<label>>& edgeToFace
173  );
174 
175  //- Like collapseMinCosAngle but return true for value < minCos
176  bool collapseCreatesFold
177  (
178  const triSurface& surf,
179  const label v1,
180  const point& pt,
181  const labelHashSet& collapsedFaces,
182  const HashTable<label, label, Hash<label>>& edgeToEdge,
183  const HashTable<label, label, Hash<label>>& edgeToFace,
184  const scalar minCos
185  );
186 
189  //static bool collapseCreatesDuplicates
190  //(
191  // const triSurface& surf,
192  // const label edgeI,
193  // const HashTable<bool, label, Hash<label>>& collapsedFaces
194  //);
195 
196  // Tracking
197 
198  //- Finds the triangle edge/point cut by the plane between
199  // a point inside/on edge of a triangle and a point outside.
200  // Returns
201  // - location on edge/point and hit()
202  // - or miss() if no intersection found
203  static surfaceLocation cutEdge
204  (
205  const triSurface& s,
206  const label triI,
207  const label excludeEdgeI,
208  const label excludePointi,
209  const point& triPoint,
210  const plane& cutPlane,
211  const point& toPoint
212  );
213 
214  //- Checks if current is on the same triangle as the endpoint
215  // and shifts it there. If so updates current and sets a hit.
216  static void snapToEnd
217  (
218  const triSurface& s,
219  const surfaceLocation& endInfo,
220  surfaceLocation& current
221  );
222 
223  //- Visits faces eFaces around start. Does not visit triangle
224  // start.triangle() nor edge excludeEdgeI.
225  // Returns edge, triangle (if more than one choice) which gets
226  // us nearer endpoint.
227  // Returns
228  // - hit() if triangle contains endpoint
229  // - triangle()=-1 if no triangle found
230  // - nearest triangle/edge otherwise
231  static surfaceLocation visitFaces
232  (
233  const triSurface& s,
234  const labelList& eFaces,
235  const surfaceLocation& start,
236  const label excludeEdgeI,
237  const label excludePointi,
238  const surfaceLocation& end,
239  const plane& cutPlane
240  );
241 
242 
243 public:
244 
245  // OBJ writing
246 
247  //- Write pointField to OBJ format file
248  static void writeOBJ
249  (
250  const fileName& fName,
251  const pointField& pts
252  );
253 
254  //- Write vertex subset to OBJ format file
255  static void writeOBJ
256  (
257  const triSurface& surf,
258  const fileName& fName,
259  const boolList& markedVerts
260  );
261 
262 
263  // Additional addressing
264 
265  //- Get all triangles using edge endpoint
266  static void getVertexTriangles
267  (
268  const triSurface& surf,
269  const label edgeI,
270  labelList& edgeTris
271  );
272 
273  //- Get all vertices (local numbering) connected to vertices of edge
275  (
276  const triSurface& surf,
277  const edge& e
278  );
279 
281  //static void orderVertices
282  //(
283  // const labelledTri& f,
284  // const label v1,
285  // const label v2,
286  // label& vA,
287  // label& vB
288  //);
289 
290  //- Get face connected to edge not facei
291  static label otherFace
292  (
293  const triSurface& surf,
294  const label facei,
295  const label edgeI
296  );
297 
298  //- Get the two edges on facei counterclockwise after edgeI
299  static void otherEdges
300  (
301  const triSurface& surf,
302  const label facei,
303  const label edgeI,
304  label& e1,
305  label& e2
306  );
307 
308  //- Get the two vertices (local numbering) on facei counterclockwise
309  // vertI
310  static void otherVertices
311  (
312  const triSurface& surf,
313  const label facei,
314  const label vertI,
315  label& vert1I,
316  label& vert2I
317  );
318 
319  //- Get edge opposite vertex (local numbering)
320  static label oppositeEdge
321  (
322  const triSurface& surf,
323  const label facei,
324  const label vertI
325  );
326 
327  //- Get vertex (local numbering) opposite edge
328  static label oppositeVertex
329  (
330  const triSurface& surf,
331  const label facei,
332  const label edgeI
333  );
334 
335  //- Returns edge label connecting v1, v2 (local numbering)
336  static label getEdge
337  (
338  const triSurface& surf,
339  const label vert1I,
340  const label vert2I
341  );
342 
343  //- Return index of triangle (or -1) using all three edges
344  static label getTriangle
345  (
346  const triSurface& surf,
347  const label e0I,
348  const label e1I,
349  const label e2I
350  );
351 
352  // Coarsening
353 
354  //- Create new triSurface by collapsing edges to edge mids.
356  (
357  const triSurface& surf,
358  const labelList& collapsableEdges
359  );
360 
361 
362  //- Face collapse status.
363  // anyEdge: any edge can be collapsed
364  // noEdge: no edge can be collapsed
365  // collapsed: already collapsed
366  // >0: edge label that can be collapsed
367  static const label ANYEDGE;
368  static const label NOEDGE;
369  static const label COLLAPSED;
370 
371  //- Create new triSurface by collapsing edges to specified
372  // positions. faceStatus allows
373  // explicit control over which faces need to be protected (see above).
374  // faceStatus gets updated to protect collapsing already collapsed
375  // faces.
377  (
378  const triSurface& surf,
379  const labelList& collapsableEdges,
380  const pointField& edgeMids,
381  labelList& faceStatus
382  );
383 
384 
385  // Refinement
386 
387  //- Refine edges by splitting to opposite vertex
388  static triSurface greenRefine
389  (
390  const triSurface& surf,
391  const labelList& refineEdges
392  );
393 
394  //- Refine face by splitting all edges. Neighbouring face is
395  // greenRefine'd.
397  (
398  const triSurface& surf,
399  const labelList& refineFaces
400  );
401 
402 
403  // Geometric
404 
405  //- Returns element in edgeIndices with minimum length
406  static label minEdge
407  (
408  const triSurface& surf,
409  const labelList& edgeIndices
410  );
411 
412  //- Returns element in edgeIndices with minimum length
413  static label maxEdge
414  (
415  const triSurface& surf,
416  const labelList& edgeIndices
417  );
418 
419  //- Merge points within distance
420  static triSurface mergePoints
421  (
422  const triSurface& surf,
423  const scalar mergeTol
424  );
425 
426  //- Triangle (unit) normal. If nearest point to triangle on edge use
427  // edge normal (calculated on the fly); if on vertex use vertex normal.
428  // Uses planarTol.
429  static vector surfaceNormal
430  (
431  const triSurface& surf,
432  const label nearestFacei,
433  const point& nearestPt
434  );
435 
436  //- On which side of surface
437  enum sideType
438  {
439  UNKNOWN, // cannot be determined (e.g. non-manifold)
440  INSIDE, // inside
441  OUTSIDE // outside
442  };
443 
444  //- If nearest point is on edgeI, determine on which side of surface
445  // sample is.
446  static sideType edgeSide
447  (
448  const triSurface& surf,
449  const point& sample,
450  const point& nearestPoint,
451  const label edgeI
452  );
453 
454  //- Given nearest point (to sample) on surface determines which side
455  // sample is. Uses either face normal, edge normal or point normal
456  // (non-trivial). Uses triangle::classify.
457  static sideType surfaceSide
458  (
459  const triSurface& surf,
460  const point& sample,
461  const label nearestFacei
462  );
463 
464  // Triangulation of faces
465 
466  //- Simple triangulation of (selected patches of) boundaryMesh. Needs
467  // polyMesh (or polyBoundaryMesh) since only at this level are the
468  // triangles on neighbouring patches connected.
469  static triSurface triangulate
470  (
471  const polyBoundaryMesh& mBesh,
473  const bool verbose = false
474  );
475 
476  //- Face-centre triangulation of (selected patches of) boundaryMesh.
477  // Needs
478  // polyMesh (or polyBoundaryMesh) since only at this level are the
479  // triangles on neighbouring patches connected.
481  (
482  const polyBoundaryMesh& mBesh,
484  const bool verbose = false
485  );
486 
487 
488  // Triangulation and interpolation
489 
490  //- Calculate linear interpolation weights for point (guaranteed to be
491  // inside triangle)
492  static void calcInterpolationWeights
493  (
494  const triPointRef&,
495  const point&,
496  FixedList<scalar, 3>& weights
497  );
498 
499  // Calculate weighting factors from samplePts to triangle it is in.
500  // Uses linear search to find triangle.
501  // Vertices are:
502  // (a b c) : vertices of the triangle abc the point is in
503  // or if the point is outside all triangles:
504  // (a b -1) : the edge ab the point is nearest to.
505  // (a -1 -1) : the vertex a the point is nearest to
506  static void calcInterpolationWeights
507  (
508  const triSurface& s,
509  const pointField& samplePts,
510  List<FixedList<label, 3>>& verts,
511  List<FixedList<scalar, 3>>& weights
512  );
513 
514  //- Do unconstrained Delaunay of points. Returns triSurface with 3D
515  // points with z=0. All triangles in region 0.
516  static triSurface delaunay2D(const List<vector2D>&);
517 
518 
519  // Tracking
520 
521  //- Test point on plane of triangle to see if on edge or point or inside
523  (
524  const triSurface&,
525  const label triI,
526  const point& trianglePoint
527  );
528 
529  //- Track on surface to get closer to point.
530  // Possible situations:
531  // - 1. reached endpoint
532  // - 2. reached edge (normal situation)
533  // - 3. reached end of surface (edge on single face)
534  // Input:
535  // - starting position+triangle/edge/point (so has to be on surface!)
536  // - (optional) previous position+triangle/edge/point to prevent
537  // going back. Set index (of triangle/edge/point) to -1 if not
538  // used.
539  // - end position+triangle/edge/point (so has to be on surface!)
540  // - plane to follow. Has to go through end point!
541  // Returns:
542  // - true if end point reached (situation 1)
543  // - new position+triangle/edge/point
544  // Caller has to check for situation 3 by checking that triangle()
545  // is not set.
547  (
548  const triSurface&,
549  const surfaceLocation& start,
550  const surfaceLocation& end,
551  const plane& cutPlane
552  );
553 
554  //- Track from edge to edge across surface. Uses trackToEdge.
555  // Not really useful by itself, more example of how to use trackToEdge.
556  // endInfo should be location on surface.
557  // hitInfo should be initialised to starting location (on surface as
558  // well). Upon return is set to end location.
559  static void track
560  (
561  const triSurface&,
562  const surfaceLocation& endInfo,
563  const plane& cutPlane,
564  surfaceLocation& hitInfo
565  );
566 };
567 
568 
569 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
570 
571 } // End namespace Foam
572 
573 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
574 
575 #endif
576 
577 // ************************************************************************* //
static labelList getVertexVertices(const triSurface &surf, const edge &e)
Get all vertices (local numbering) connected to vertices of edge.
A triangle primitive used to calculate face normals and swept volumes.
Definition: triangle.H:57
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A class for handling file names.
Definition: fileName.H:69
const double e
Elementary charge.
Definition: doubleFloat.H:78
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
static void track(const triSurface &, const surfaceLocation &endInfo, const plane &cutPlane, surfaceLocation &hitInfo)
Track from edge to edge across surface. Uses trackToEdge.
static triSurface redGreenRefine(const triSurface &surf, const labelList &refineFaces)
Refine face by splitting all edges. Neighbouring face is.
static label maxEdge(const triSurface &surf, const labelList &edgeIndices)
Returns element in edgeIndices with minimum length.
Geometric class that creates a 2D plane and can return the intersection point between a line and the ...
Definition: plane.H:60
static vector surfaceNormal(const triSurface &surf, const label nearestFacei, const point &nearestPt)
Triangle (unit) normal. If nearest point to triangle on edge use.
static const label NOEDGE
static triSurface delaunay2D(const List< vector2D > &)
Do unconstrained Delaunay of points. Returns triSurface with 3D.
static triSurface collapseEdges(const triSurface &surf, const labelList &collapsableEdges)
Create new triSurface by collapsing edges to edge mids.
static label otherFace(const triSurface &surf, const label facei, const label edgeI)
Get face connected to edge not facei.
triSurface triangulateFaceCentre(const polyBoundaryMesh &mBesh, const labelHashSet &includePatches, const bool verbose=false)
Face-centre triangulation of (selected patches of) boundaryMesh.
static label oppositeVertex(const triSurface &surf, const label facei, const label edgeI)
Get vertex (local numbering) opposite edge.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
static triSurface mergePoints(const triSurface &surf, const scalar mergeTol)
Merge points within distance.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
static sideType edgeSide(const triSurface &surf, const point &sample, const point &nearestPoint, const label edgeI)
If nearest point is on edgeI, determine on which side of surface.
static label minEdge(const triSurface &surf, const labelList &edgeIndices)
Returns element in edgeIndices with minimum length.
static void otherVertices(const triSurface &surf, const label facei, const label vertI, label &vert1I, label &vert2I)
Get the two vertices (local numbering) on facei counterclockwise.
PointFrompoint toPoint(const Foam::point &p)
static const label COLLAPSED
An STL-conforming hash table.
Definition: HashTable.H:61
static surfaceLocation trackToEdge(const triSurface &, const surfaceLocation &start, const surfaceLocation &end, const plane &cutPlane)
Track on surface to get closer to point.
Foam::polyBoundaryMesh.
A collection of tools for triSurface.
sideType
On which side of surface.
labelHashSet includePatches
static void calcInterpolationWeights(const triPointRef &, const point &, FixedList< scalar, 3 > &weights)
Calculate linear interpolation weights for point (guaranteed to be.
static label getTriangle(const triSurface &surf, const label e0I, const label e1I, const label e2I)
Return index of triangle (or -1) using all three edges.
static sideType surfaceSide(const triSurface &surf, const point &sample, const label nearestFacei)
Given nearest point (to sample) on surface determines which side.
label newPointi
Definition: readKivaGrid.H:501
static void otherEdges(const triSurface &surf, const label facei, const label edgeI, label &e1, label &e2)
Get the two edges on facei counterclockwise after edgeI.
Hash function class for primitives. All non-primitives used to hash entries on hash tables likely nee...
Definition: Hash.H:54
static triSurface triangulate(const polyBoundaryMesh &mBesh, const labelHashSet &includePatches, const bool verbose=false)
Simple triangulation of (selected patches of) boundaryMesh. Needs.
Contains information about location on a triSurface.
static void writeOBJ(const fileName &fName, const pointField &pts)
Write pointField to OBJ format file.
Triangulated surface description with patch information.
Definition: triSurface.H:65
static const label ANYEDGE
Face collapse status.
static label getEdge(const triSurface &surf, const label vert1I, const label vert2I)
Returns edge label connecting v1, v2 (local numbering)
static void getVertexTriangles(const triSurface &surf, const label edgeI, labelList &edgeTris)
Get all triangles using edge endpoint.
static surfaceLocation classify(const triSurface &, const label triI, const point &trianglePoint)
Test point on plane of triangle to see if on edge or point or inside.
static label oppositeEdge(const triSurface &surf, const label facei, const label vertI)
Get edge opposite vertex (local numbering)
Namespace for OpenFOAM.