collapseBase.C
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 \*---------------------------------------------------------------------------*/
25 
26 #include "collapseBase.H"
27 #include "triSurfaceTools.H"
28 #include "argList.H"
29 #include "OFstream.H"
30 #include "SubList.H"
31 #include "labelPair.H"
32 #include "meshTools.H"
33 #include "OSspecific.H"
34 
35 using namespace Foam;
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
40 //static void writeRegionOBJ
41 //(
42 // const triSurface& surf,
43 // const label regionI,
44 // const labelList& collapseRegion,
45 // const labelList& outsideVerts
46 //)
47 //{
48 // fileName dir("regions");
49 //
50 // mkDir(dir);
51 // fileName regionName(dir / "region_" + name(regionI) + ".obj");
52 //
53 // Pout<< "Dumping region " << regionI << " to file " << regionName << endl;
54 //
55 // boolList include(surf.size(), false);
56 //
57 // forAll(collapseRegion, facei)
58 // {
59 // if (collapseRegion[facei] == regionI)
60 // {
61 // include[facei] = true;
62 // }
63 // }
64 //
65 // labelList pointMap, faceMap;
66 //
67 // triSurface regionSurf(surf.subsetMesh(include, pointMap, faceMap));
68 //
69 // Pout<< "Region " << regionI << " surface:" << nl;
70 // regionSurf.writeStats(Pout);
71 //
72 // regionSurf.write(regionName);
73 //
74 //
75 // // Dump corresponding outside vertices.
76 // fileName pointsName(dir / "regionPoints_" + name(regionI) + ".obj");
77 //
78 // Pout<< "Dumping region " << regionI << " points to file " << pointsName
79 // << endl;
80 //
81 // OFstream str(pointsName);
82 //
83 // forAll(outsideVerts, i)
84 // {
85 // meshTools::writeOBJ(str, surf.localPoints()[outsideVerts[i]]);
86 // }
87 //}
88 
89 
90 // Split triangle into multiple triangles because edge e being split
91 // into multiple edges.
92 static void splitTri
93 (
94  const labelledTri& f,
95  const edge& e,
96  const labelList& splitPoints,
98 )
99 {
100  //label oldNTris = tris.size();
101 
102  label fp = findIndex(f, e[0]);
103  label fp1 = f.fcIndex(fp);
104  label fp2 = f.fcIndex(fp1);
105 
106  if (f[fp1] == e[1])
107  {
108  // Split triangle along fp to fp1
109  tris.append(labelledTri(f[fp2], f[fp], splitPoints[0], f.region()));
110 
111  for (label i = 1; i < splitPoints.size(); i++)
112  {
113  tris.append
114  (
116  (
117  f[fp2],
118  splitPoints[i-1],
119  splitPoints[i],
120  f.region()
121  )
122  );
123  }
124 
125  tris.append
126  (
128  (
129  f[fp2],
130  splitPoints.last(),
131  f[fp1],
132  f.region()
133  )
134  );
135  }
136  else if (f[fp2] == e[1])
137  {
138  // Split triangle along fp2 to fp. (Reverse order of splitPoints)
139 
140  tris.append
141  (
143  (
144  f[fp1],
145  f[fp2],
146  splitPoints.last(),
147  f.region()
148  )
149  );
150 
151  for (label i = splitPoints.size()-1; i > 0; --i)
152  {
153  tris.append
154  (
156  (
157  f[fp1],
158  splitPoints[i],
159  splitPoints[i-1],
160  f.region()
161  )
162  );
163  }
164 
165  tris.append
166  (
168  (
169  f[fp1],
170  splitPoints[0],
171  f[fp],
172  f.region()
173  )
174  );
175  }
176  else
177  {
179  << "Edge " << e << " not part of triangle " << f
180  << " fp:" << fp
181  << " fp1:" << fp1
182  << " fp2:" << fp2
183  << abort(FatalError);
184  }
185 
186  //Pout<< "Split face " << f << " along edge " << e
187  // << " into triangles:" << endl;
188  //
189  //for (label i = oldNTris; i < tris.size(); i++)
190  //{
191  // Pout<< " " << tris[i] << nl;
192  //}
193 }
194 
195 
196 // Insert scalar into sortedVerts/sortedWeights so the weights are in
197 // incrementing order.
198 static bool insertSorted
199 (
200  const label vertI,
201  const scalar weight,
202 
203  labelList& sortedVerts,
204  scalarField& sortedWeights
205 )
206 {
207  if (findIndex(sortedVerts, vertI) != -1)
208  {
210  << " which is already in list of sorted vertices "
211  << sortedVerts << abort(FatalError);
212  }
213 
214  if (weight <= 0 || weight >= 1)
215  {
217  << " with illegal weight " << weight
218  << " into list of sorted vertices "
219  << sortedVerts << abort(FatalError);
220  }
221 
222 
223  label insertI = sortedVerts.size();
224 
225  forAll(sortedVerts, sortedI)
226  {
227  scalar w = sortedWeights[sortedI];
228 
229  if (mag(w - weight) < SMALL)
230  {
232  << "Trying to insert weight " << weight << " which is close to"
233  << " existing weight " << w << " in " << sortedWeights
234  << endl;
235  }
236 
237  if (w > weight)
238  {
239  // Start inserting before sortedI.
240  insertI = sortedI;
241  break;
242  }
243  }
244 
245 
246  label sz = sortedWeights.size();
247 
248  sortedWeights.setSize(sz + 1);
249  sortedVerts.setSize(sz + 1);
250 
251  // Leave everything up to (not including) insertI intact.
252 
253  // Make space by copying from insertI up.
254  for (label i = sz-1; i >= insertI; --i)
255  {
256  sortedWeights[i+1] = sortedWeights[i];
257  sortedVerts[i+1] = sortedVerts[i];
258  }
259  sortedWeights[insertI] = weight;
260  sortedVerts[insertI] = vertI;
261 
262  return true;
263 }
264 
265 
266 // Is triangle candidate for collapse? Small height or small quality
267 bool isSliver
268 (
269  const triSurface& surf,
270  const scalar minLen,
271  const scalar minQuality,
272  const label facei,
273  const label edgeI
274 )
275 {
276  const pointField& localPoints = surf.localPoints();
277 
278  // Check
279  // - opposite vertex projects onto base edge
280  // - normal distance is small
281  // - or triangle quality is small
282 
283  label opposite0 =
285  (
286  surf,
287  facei,
288  edgeI
289  );
290 
291  const edge& e = surf.edges()[edgeI];
292  const labelledTri& f = surf[facei];
293 
294  pointHit pHit =
295  e.line(localPoints).nearestDist
296  (
297  localPoints[opposite0]
298  );
299 
300  if
301  (
302  pHit.hit()
303  && (
304  pHit.distance() < minLen
305  || f.tri(surf.points()).quality() < minQuality
306  )
307  )
308  {
309  // Remove facei and split all other faces using this
310  // edge. This is done by 'replacing' the edgeI with the
311  // opposite0 vertex
312  //Pout<< "Splitting face " << facei << " since distance "
313  // << pHit.distance()
314  // << " from vertex " << opposite0
315  // << " to edge " << edgeI
316  // << " points "
317  // << localPoints[e[0]]
318  // << localPoints[e[1]]
319  // << " is too small or triangle quality "
320  // << f.tri(surf.points()).quality()
321  // << " too small." << endl;
322 
323  return true;
324  }
325  else
326  {
327  return false;
328  }
329 }
330 
331 
332 // Mark all faces that are going to be collapsed.
333 // faceToEdge: per face -1 or the base edge of the face.
334 static void markCollapsedFaces
335 (
336  const triSurface& surf,
337  const scalar minLen,
338  const scalar minQuality,
339  labelList& faceToEdge
340 )
341 {
342  faceToEdge.setSize(surf.size());
343  faceToEdge = -1;
344 
345  const labelListList& edgeFaces = surf.edgeFaces();
346 
347  forAll(edgeFaces, edgeI)
348  {
349  const labelList& eFaces = surf.edgeFaces()[edgeI];
350 
351  forAll(eFaces, i)
352  {
353  label facei = eFaces[i];
354 
355  bool isCandidate = isSliver(surf, minLen, minQuality, facei, edgeI);
356 
357  if (isCandidate)
358  {
359  // Mark face as being collapsed
360  if (faceToEdge[facei] != -1)
361  {
363  << "Cannot collapse face " << facei << " since "
364  << " is marked to be collapsed both to edge "
365  << faceToEdge[facei] << " and " << edgeI
366  << abort(FatalError);
367  }
368 
369  faceToEdge[facei] = edgeI;
370  }
371  }
372  }
373 }
374 
375 
376 // Recurse through collapsed faces marking all of them with regionI (in
377 // collapseRegion)
378 static void markRegion
379 (
380  const triSurface& surf,
381  const labelList& faceToEdge,
382  const label regionI,
383  const label facei,
384  labelList& collapseRegion
385 )
386 {
387  if (faceToEdge[facei] == -1 || collapseRegion[facei] != -1)
388  {
390  << "Problem : crossed into uncollapsed/regionized face"
391  << abort(FatalError);
392  }
393 
394  collapseRegion[facei] = regionI;
395 
396  // Recurse across edges to collapsed neighbours
397 
398  const labelList& fEdges = surf.faceEdges()[facei];
399 
400  forAll(fEdges, fEdgeI)
401  {
402  label edgeI = fEdges[fEdgeI];
403 
404  const labelList& eFaces = surf.edgeFaces()[edgeI];
405 
406  forAll(eFaces, i)
407  {
408  label nbrFacei = eFaces[i];
409 
410  if (faceToEdge[nbrFacei] != -1)
411  {
412  if (collapseRegion[nbrFacei] == -1)
413  {
414  markRegion
415  (
416  surf,
417  faceToEdge,
418  regionI,
419  nbrFacei,
420  collapseRegion
421  );
422  }
423  else if (collapseRegion[nbrFacei] != regionI)
424  {
426  << "Edge:" << edgeI << " between face " << facei
427  << " with region " << regionI
428  << " and face " << nbrFacei
429  << " with region " << collapseRegion[nbrFacei]
430  << endl;
431  }
432  }
433  }
434  }
435 }
436 
437 
438 // Mark every face with region (in collapseRegion) (or -1).
439 // Return number of regions.
440 static label markRegions
441 (
442  const triSurface& surf,
443  const labelList& faceToEdge,
444  labelList& collapseRegion
445 )
446 {
447  label regionI = 0;
448 
449  forAll(faceToEdge, facei)
450  {
451  if (collapseRegion[facei] == -1 && faceToEdge[facei] != -1)
452  {
453  //Pout<< "markRegions : Marking region:" << regionI
454  // << " starting from face " << facei << endl;
455 
456  // Collapsed face. Mark connected region with current region number
457  markRegion(surf, faceToEdge, regionI++, facei, collapseRegion);
458  }
459  }
460  return regionI;
461 }
462 
463 
464 // Type of region.
465 // -1 : edge inbetween uncollapsed faces.
466 // -2 : edge inbetween collapsed faces
467 // >=0 : edge inbetween uncollapsed and collapsed region. Returns region.
468 static label edgeType
469 (
470  const triSurface& surf,
471  const labelList& collapseRegion,
472  const label edgeI
473 )
474 {
475  const labelList& eFaces = surf.edgeFaces()[edgeI];
476 
477  // Detect if edge is inbetween collapseRegion and non-collapse face
478  bool usesUncollapsed = false;
479  label usesRegion = -1;
480 
481  forAll(eFaces, i)
482  {
483  label facei = eFaces[i];
484 
485  label region = collapseRegion[facei];
486 
487  if (region == -1)
488  {
489  usesUncollapsed = true;
490  }
491  else if (usesRegion == -1)
492  {
493  usesRegion = region;
494  }
495  else if (usesRegion != region)
496  {
498  }
499  else
500  {
501  // Equal regions.
502  }
503  }
504 
505  if (usesUncollapsed)
506  {
507  if (usesRegion == -1)
508  {
509  // uncollapsed faces only.
510  return -1;
511  }
512  else
513  {
514  // between uncollapsed and collapsed.
515  return usesRegion;
516  }
517  }
518  else
519  {
520  if (usesRegion == -1)
521  {
523  return -2;
524  }
525  else
526  {
527  return -2;
528  }
529  }
530 }
531 
532 
533 // Get points on outside edge of region (= outside points)
534 static labelListList getOutsideVerts
535 (
536  const triSurface& surf,
537  const labelList& collapseRegion,
538  const label nRegions
539 )
540 {
541  const labelListList& edgeFaces = surf.edgeFaces();
542 
543  // Per region all the outside vertices.
544  labelListList outsideVerts(nRegions);
545 
546  forAll(edgeFaces, edgeI)
547  {
548  // Detect if edge is inbetween collapseRegion and non-collapse face
549  label regionI = edgeType(surf, collapseRegion, edgeI);
550 
551  if (regionI >= 0)
552  {
553  // Edge borders both uncollapsed face and collapsed face on region
554  // usesRegion.
555 
556  const edge& e = surf.edges()[edgeI];
557 
558  labelList& regionVerts = outsideVerts[regionI];
559 
560  // Add both edge points to regionVerts.
561  forAll(e, eI)
562  {
563  label v = e[eI];
564 
565  if (findIndex(regionVerts, v) == -1)
566  {
567  label sz = regionVerts.size();
568  regionVerts.setSize(sz+1);
569  regionVerts[sz] = v;
570  }
571  }
572  }
573  }
574 
575  return outsideVerts;
576 }
577 
578 
579 // n^2 search for furthest removed point pair.
580 static labelPair getSpanPoints
581 (
582  const triSurface& surf,
583  const labelList& outsideVerts
584 )
585 {
586  const pointField& localPoints = surf.localPoints();
587 
588  scalar maxDist = -GREAT;
589  labelPair maxPair;
590 
591  forAll(outsideVerts, i)
592  {
593  label v0 = outsideVerts[i];
594 
595  for (label j = i+1; j < outsideVerts.size(); j++)
596  {
597  label v1 = outsideVerts[j];
598 
599  scalar d = mag(localPoints[v0] - localPoints[v1]);
600 
601  if (d > maxDist)
602  {
603  maxDist = d;
604  maxPair[0] = v0;
605  maxPair[1] = v1;
606  }
607  }
608  }
609 
610  return maxPair;
611 }
612 
613 
614 // Project all non-span points onto the span edge.
615 static void projectNonSpanPoints
616 (
617  const triSurface& surf,
618  const labelList& outsideVerts,
619  const labelPair& spanPair,
620  labelList& sortedVertices,
621  scalarField& sortedWeights
622 )
623 {
624  const point& p0 = surf.localPoints()[spanPair[0]];
625  const point& p1 = surf.localPoints()[spanPair[1]];
626 
627  forAll(outsideVerts, i)
628  {
629  label v = outsideVerts[i];
630 
631  if (v != spanPair[0] && v != spanPair[1])
632  {
633  // Is a non-span point. Project onto spanning edge.
634 
635  pointHit pHit =
636  linePointRef(p0, p1).nearestDist
637  (
638  surf.localPoints()[v]
639  );
640 
641  if (!pHit.hit())
642  {
644  << abort(FatalError);
645  }
646 
647  scalar w = mag(pHit.hitPoint() - p0) / mag(p1 - p0);
648 
649  insertSorted(v, w, sortedVertices, sortedWeights);
650  }
651  }
652 }
653 
654 
655 // Slice part of the orderVertices (and optionally reverse) for this edge.
656 static void getSplitVerts
657 (
658  const triSurface& surf,
659  const label regionI,
660  const labelPair& spanPoints,
661  const labelList& orderedVerts,
662  const scalarField& orderedWeights,
663  const label edgeI,
664 
665  labelList& splitVerts,
666  scalarField& splitWeights
667 )
668 {
669  const edge& e = surf.edges()[edgeI];
670  const label sz = orderedVerts.size();
671 
672  if (e[0] == spanPoints[0])
673  {
674  // Edge in same order as spanPoints&orderedVerts. Keep order.
675 
676  if (e[1] == spanPoints[1])
677  {
678  // Copy all.
679  splitVerts = orderedVerts;
680  splitWeights = orderedWeights;
681  }
682  else
683  {
684  // Copy upto (but not including) e[1]
685  label i1 = findIndex(orderedVerts, e[1]);
686  splitVerts = SubList<label>(orderedVerts, i1, 0);
687  splitWeights = SubList<scalar>(orderedWeights, i1, 0);
688  }
689  }
690  else if (e[0] == spanPoints[1])
691  {
692  // Reverse.
693 
694  if (e[1] == spanPoints[0])
695  {
696  // Copy all.
697  splitVerts = orderedVerts;
698  reverse(splitVerts);
699  splitWeights = orderedWeights;
700  reverse(splitWeights);
701  }
702  else
703  {
704  // Copy downto (but not including) e[1]
705 
706  label i1 = findIndex(orderedVerts, e[1]);
707  splitVerts = SubList<label>(orderedVerts, sz-(i1+1), i1+1);
708  reverse(splitVerts);
709  splitWeights = SubList<scalar>(orderedWeights, sz-(i1+1), i1+1);
710  reverse(splitWeights);
711  }
712  }
713  else if (e[1] == spanPoints[0])
714  {
715  // Reverse.
716 
717  // Copy upto (but not including) e[0]
718 
719  label i0 = findIndex(orderedVerts, e[0]);
720  splitVerts = SubList<label>(orderedVerts, i0, 0);
721  reverse(splitVerts);
722  splitWeights = SubList<scalar>(orderedWeights, i0, 0);
723  reverse(splitWeights);
724  }
725  else if (e[1] == spanPoints[1])
726  {
727  // Copy from (but not including) e[0] to end
728 
729  label i0 = findIndex(orderedVerts, e[0]);
730  splitVerts = SubList<label>(orderedVerts, sz-(i0+1), i0+1);
731  splitWeights = SubList<scalar>(orderedWeights, sz-(i0+1), i0+1);
732  }
733  else
734  {
735  label i0 = findIndex(orderedVerts, e[0]);
736  label i1 = findIndex(orderedVerts, e[1]);
737 
738  if (i0 == -1 || i1 == -1)
739  {
741  << "Did not find edge in projected vertices." << nl
742  << "region:" << regionI << nl
743  << "spanPoints:" << spanPoints
744  << " coords:" << surf.localPoints()[spanPoints[0]]
745  << surf.localPoints()[spanPoints[1]] << nl
746  << "edge:" << edgeI
747  << " verts:" << e
748  << " coords:" << surf.localPoints()[e[0]]
749  << surf.localPoints()[e[1]] << nl
750  << "orderedVerts:" << orderedVerts << nl
751  << abort(FatalError);
752  }
753 
754  if (i0 < i1)
755  {
756  splitVerts = SubList<label>(orderedVerts, i1-i0-1, i0+1);
757  splitWeights = SubList<scalar>(orderedWeights, i1-i0-1, i0+1);
758  }
759  else
760  {
761  splitVerts = SubList<label>(orderedVerts, i0-i1-1, i1+1);
762  reverse(splitVerts);
763  splitWeights = SubList<scalar>(orderedWeights, i0-i1-1, i1+1);
764  reverse(splitWeights);
765  }
766  }
767 }
768 
769 
771 (
772  triSurface& surf,
773  const scalar minLen,
774  const scalar minQuality
775 )
776 {
777  label nTotalSplit = 0;
778 
779  label iter = 0;
780 
781  while (true)
782  {
783  // Detect faces to collapse
784  // ~~~~~~~~~~~~~~~~~~~~~~~~
785 
786  // -1 or edge the face is collapsed onto.
787  labelList faceToEdge(surf.size(), -1);
788 
789  // Calculate faceToEdge (face collapses)
790  markCollapsedFaces(surf, minLen, minQuality, faceToEdge);
791 
792 
793  // Find regions of connected collapsed faces
794  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
795 
796  // per face -1 or region
797  labelList collapseRegion(surf.size(), -1);
798 
799  label nRegions = markRegions(surf, faceToEdge, collapseRegion);
800 
801  //Pout<< "Detected " << nRegions << " regions of faces to be collapsed"
802  // << nl << endl;
803 
804  // Pick up all vertices on outside of region
805  labelListList outsideVerts
806  (
807  getOutsideVerts(surf, collapseRegion, nRegions)
808  );
809 
810  // For all regions determine maximum distance between points
811  List<labelPair> spanPoints(nRegions);
812  labelListList orderedVertices(nRegions);
813  List<scalarField> orderedWeights(nRegions);
814 
815  forAll(spanPoints, regionI)
816  {
817  spanPoints[regionI] = getSpanPoints(surf, outsideVerts[regionI]);
818 
819  //Pout<< "For region " << regionI << " found extrema at points "
820  // << surf.localPoints()[spanPoints[regionI][0]]
821  // << surf.localPoints()[spanPoints[regionI][1]]
822  // << endl;
823 
824  // Project all non-span points onto the span edge.
825  projectNonSpanPoints
826  (
827  surf,
828  outsideVerts[regionI],
829  spanPoints[regionI],
830  orderedVertices[regionI],
831  orderedWeights[regionI]
832  );
833 
834  //Pout<< "For region:" << regionI
835  // << " span:" << spanPoints[regionI]
836  // << " orderedVerts:" << orderedVertices[regionI]
837  // << " orderedWeights:" << orderedWeights[regionI]
838  // << endl;
839 
840  //writeRegionOBJ
841  //(
842  // surf,
843  // regionI,
844  // collapseRegion,
845  // outsideVerts[regionI]
846  //);
847 
848  //Pout<< endl;
849  }
850 
851 
852 
853  // Actually split the edges
854  // ~~~~~~~~~~~~~~~~~~~~~~~~
855 
856 
857  const List<labelledTri>& localFaces = surf.localFaces();
858  const edgeList& edges = surf.edges();
859 
860  label nSplit = 0;
861 
862  // Storage for new triangles.
863  DynamicList<labelledTri> newTris(surf.size());
864 
865  // Whether face has been dealt with (either copied/split or deleted)
866  boolList faceHandled(surf.size(), false);
867 
868 
869  forAll(edges, edgeI)
870  {
871  const edge& e = edges[edgeI];
872 
873  // Detect if edge is inbetween collapseRegion and non-collapse face
874  label regionI = edgeType(surf, collapseRegion, edgeI);
875 
876  if (regionI == -2)
877  {
878  // inbetween collapsed faces. nothing needs to be done.
879  }
880  else if (regionI == -1)
881  {
882  // edge inbetween uncollapsed faces. Handle these later on.
883  }
884  else
885  {
886  // some faces around edge are collapsed.
887 
888  // Find additional set of points on edge to be used to split
889  // the remaining faces.
890 
891  labelList splitVerts;
892  scalarField splitWeights;
893  getSplitVerts
894  (
895  surf,
896  regionI,
897  spanPoints[regionI],
898  orderedVertices[regionI],
899  orderedWeights[regionI],
900  edgeI,
901 
902  splitVerts,
903  splitWeights
904  );
905 
906  if (splitVerts.size())
907  {
908  // Split edge using splitVerts. All non-collapsed triangles
909  // using edge will get split.
910 
911  //{
912  // const pointField& localPoints = surf.localPoints();
913  // Pout<< "edge " << edgeI << ' ' << e
914  // << " points "
915  // << localPoints[e[0]] << ' ' << localPoints[e[1]]
916  // << " split into edges with extra points:"
917  // << endl;
918  // forAll(splitVerts, i)
919  // {
920  // Pout<< " " << splitVerts[i] << " weight "
921  // << splitWeights[i] << nl;
922  // }
923  //}
924 
925  const labelList& eFaces = surf.edgeFaces()[edgeI];
926 
927  forAll(eFaces, i)
928  {
929  label facei = eFaces[i];
930 
931  if (!faceHandled[facei] && faceToEdge[facei] == -1)
932  {
933  // Split face to use vertices.
934  splitTri
935  (
936  localFaces[facei],
937  e,
938  splitVerts,
939  newTris
940  );
941 
942  faceHandled[facei] = true;
943 
944  nSplit++;
945  }
946  }
947  }
948  }
949  }
950 
951  // Copy all unsplit faces
952  forAll(faceHandled, facei)
953  {
954  if (!faceHandled[facei] && faceToEdge[facei] == -1)
955  {
956  newTris.append(localFaces[facei]);
957  }
958  }
959 
960  Info<< "collapseBase : collapsing " << nSplit
961  << " triangles by splitting their base edge."
962  << endl;
963 
964  nTotalSplit += nSplit;
965 
966  if (nSplit == 0)
967  {
968  break;
969  }
970 
971  // Pack the triangles
972  newTris.shrink();
973 
974  //Pout<< "Resetting surface from " << surf.size() << " to "
975  // << newTris.size() << " triangles" << endl;
976  surf = triSurface(newTris, surf.patches(), surf.localPoints());
977 
978  //{
979  // fileName fName("bla" + name(iter) + ".obj");
980  // Pout<< "Writing surf to " << fName << endl;
981  // surf.write(fName);
982  //}
983 
984  iter++;
985  }
986 
987  // Remove any unused vertices
988  surf = triSurface(surf.localFaces(), surf.patches(), surf.localPoints());
989 
990  return nTotalSplit;
991 }
992 
993 
994 // ************************************************************************* //
linePointRef line(const pointField &) const
Return edge line.
Definition: edgeI.H:169
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
label collapseBase(triSurface &surf, const scalar minLen, const scalar minQuality)
Keep collapsing all triangles whose height is < minLen or quality < minQ.
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
label region() const
Return region label.
Definition: labelledTriI.H:68
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const List< Face > & localFaces() const
Return patch faces addressing into local point list.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
const Field< PointType > & points() const
Return reference to global points.
Routines collapse sliver triangles by splitting the base edge.
static label oppositeVertex(const triSurface &surf, const label facei, const label edgeI)
Get vertex (local numbering) opposite edge.
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
A List obtained as a section of another List.
Definition: SubList.H:53
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
line< point, const point & > linePointRef
Line using referred points.
Definition: linePointRef.H:45
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
const geometricSurfacePatchList & patches() const
Definition: triSurface.H:314
const labelListList & faceEdges() const
Return face-edge addressing.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:292
Triangle with additional region number.
Definition: labelledTri.H:57
label fcIndex(const label i) const
Return the forward circular index, i.e. the next index.
Definition: FixedListI.H:115
errorManip< error > abort(error &err)
Definition: errorManip.H:131
void reverse(UList< T > &, const label n)
Definition: UListI.H:322
static const char nl
Definition: Ostream.H:262
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
void setSize(const label)
Reset size of List.
Definition: List.C:295
#define WarningInFunction
Report a warning using Foam::Warning.
triPointRef tri(const pointField &) const
Return the triangle.
Definition: triFaceI.H:152
const Point & hitPoint() const
Return hit point.
Definition: PointHit.H:126
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
const labelListList & edgeFaces() const
Return edge-face addressing.
PointHit< Point > nearestDist(const Point &p) const
Return nearest distance to line from a given point.
Definition: lineI.H:95
bool hit() const
Is there a hit.
Definition: PointHit.H:120
T & last()
Return the last element of the list.
Definition: UListI.H:128
Triangulated surface description with patch information.
Definition: triSurface.H:65
const Field< PointType > & localPoints() const
Return pointField of points in patch.
Namespace for OpenFOAM.