surfaceSplitNonManifolds.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2019-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Application
28 surfaceSplitNonManifolds
29
30Group
31 grpSurfaceUtilities
32
33Description
34 Takes multiply connected surface and tries to split surface at
35 multiply connected edges by duplicating points.
36
37 Introduces concept of
38 - borderEdge. Edge with 4 faces connected to it.
39 - borderPoint. Point connected to exactly 2 borderEdges.
40 - borderLine. Connected list of borderEdges.
41
42 By duplicating borderPoints this will split 'borderLines'. As a
43 preprocessing step it can detect borderEdges without any borderPoints
44 and explicitly split these triangles.
45
46 The problems in this algorithm are:
47 - determining which two (of the four) faces form a surface. Done by walking
48 face-edge-face while keeping and edge or point on the borderEdge
49 borderPoint.
50 - determining the outwards pointing normal to be used to slightly offset the
51 duplicated point.
52
53 Uses sortedEdgeFaces quite a bit.
54
55 Is tested on simple borderLines resulting from extracting a surface
56 from a hex mesh. Will quite possibly go wrong on more complicated border
57 lines (i.e. ones forming a loop).
58
59 Dumps surface every so often since might take a long time to complete.
60
61\*---------------------------------------------------------------------------*/
62
63#include "argList.H"
64#include "triSurface.H"
65#include "OFstream.H"
66#include "ListOps.H"
67#include "triSurfaceTools.H"
68
69using namespace Foam;
70
71// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
72
73void writeOBJ(Ostream& os, const pointField& pts)
74{
75 forAll(pts, i)
76 {
77 const point& pt = pts[i];
78
79 os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << endl;
80 }
81}
82
83
84void dumpPoints(const triSurface& surf, const labelList& borderPoint)
85{
86 fileName fName("borderPoints.obj");
87
88 Info<< "Dumping borderPoints as obj file: " << fName << endl;
89
90 OFstream os(fName);
91
92 forAll(borderPoint, pointi)
93 {
94 if (borderPoint[pointi] != -1)
95 {
96 const point& pt = surf.localPoints()[pointi];
97
98 os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << endl;
99 }
100 }
101}
102
103
104void dumpEdges(const triSurface& surf, const boolList& borderEdge)
105{
106 fileName fName("borderEdges.obj");
107
108 Info<< "Dumping borderEdges as obj file: " << fName << endl;
109
110 OFstream os(fName);
111
112 writeOBJ(os, surf.localPoints());
113
114 forAll(borderEdge, edgeI)
115 {
116 if (borderEdge[edgeI])
117 {
118 const edge& e = surf.edges()[edgeI];
119
120 os << "l " << e.start()+1 << ' ' << e.end()+1 << endl;
121 }
122 }
123}
124
125
126void dumpFaces
127(
128 const fileName& fName,
129 const triSurface& surf,
130 const Map<label>& connectedFaces
131)
132{
133 Info<< "Dumping connectedFaces as obj file: " << fName << endl;
134
135 OFstream os(fName);
136
137 forAllConstIters(connectedFaces, iter)
138 {
139 point ctr = surf.localFaces()[iter.key()].centre(surf.localPoints());
140
141 os << "v " << ctr.x() << ' ' << ctr.y() << ' ' << ctr.z() << endl;
142 }
143}
144
145
146void testSortedEdgeFaces(const triSurface& surf)
147{
148 const labelListList& edgeFaces = surf.edgeFaces();
149 const labelListList& sortedEdgeFaces = surf.sortedEdgeFaces();
150
151 forAll(edgeFaces, edgeI)
152 {
153 const labelList& myFaces = edgeFaces[edgeI];
154 const labelList& sortMyFaces = sortedEdgeFaces[edgeI];
155
156 forAll(myFaces, i)
157 {
158 if (!sortMyFaces.found(myFaces[i]))
159 {
160 FatalErrorInFunction << abort(FatalError);
161 }
162 }
163 forAll(sortMyFaces, i)
164 {
165 if (!myFaces.found(sortMyFaces[i]))
166 {
167 FatalErrorInFunction << abort(FatalError);
168 }
169 }
170 }
171}
172
173
174// Mark off all border edges. Return number of border edges.
175label markBorderEdges
176(
177 const bool debug,
178 const triSurface& surf,
179 boolList& borderEdge
180)
181{
182 label nBorderEdges = 0;
183
184 const labelListList& edgeFaces = surf.edgeFaces();
185
186 forAll(edgeFaces, edgeI)
187 {
188 if (edgeFaces[edgeI].size() == 4)
189 {
190 borderEdge[edgeI] = true;
191
192 nBorderEdges++;
193 }
194 }
195
196 if (debug)
197 {
198 dumpEdges(surf, borderEdge);
199 }
200
201 return nBorderEdges;
202}
203
204
205// Mark off all border points. Return number of border points. Border points
206// marked by setting value to newly introduced point.
207label markBorderPoints
208(
209 const bool debug,
210 const triSurface& surf,
211 const boolList& borderEdge,
212 labelList& borderPoint
213)
214{
215 label nPoints = surf.nPoints();
216
217 const labelListList& pointEdges = surf.pointEdges();
218
219 forAll(pointEdges, pointi)
220 {
221 const labelList& pEdges = pointEdges[pointi];
222
223 label nBorderEdges = 0;
224
225 forAll(pEdges, i)
226 {
227 if (borderEdge[pEdges[i]])
228 {
229 nBorderEdges++;
230 }
231 }
232
233 if (nBorderEdges == 2 && borderPoint[pointi] == -1)
234 {
235 borderPoint[pointi] = nPoints++;
236 }
237 }
238
239 label nBorderPoints = nPoints - surf.nPoints();
240
241 if (debug)
242 {
243 dumpPoints(surf, borderPoint);
244 }
245
246 return nBorderPoints;
247}
248
249
250// Get minimum length of edges connected to pointi
251// Serves to get some local length scale.
252scalar minEdgeLen(const triSurface& surf, const label pointi)
253{
254 const pointField& points = surf.localPoints();
255
256 const labelList& pEdges = surf.pointEdges()[pointi];
257
258 scalar minLen = GREAT;
259
260 forAll(pEdges, i)
261 {
262 label edgeI = pEdges[i];
263
264 scalar len = surf.edges()[edgeI].mag(points);
265
266 if (len < minLen)
267 {
268 minLen = len;
269 }
270 }
271 return minLen;
272}
273
274
275// Find edge among edgeLabels with endpoints v0,v1
276label findEdge
277(
278 const triSurface& surf,
279 const labelList& edgeLabels,
280 const label v0,
281 const label v1
282)
283{
284 forAll(edgeLabels, i)
285 {
286 label edgeI = edgeLabels[i];
287
288 const edge& e = surf.edges()[edgeI];
289
290 if
291 (
292 (
293 e.start() == v0
294 && e.end() == v1
295 )
296 || (
297 e.start() == v1
298 && e.end() == v0
299 )
300 )
301 {
302 return edgeI;
303 }
304 }
305
306
308 << ' ' << v1 << " in candidates " << edgeLabels
309 << " with vertices:" << UIndirectList<edge>(surf.edges(), edgeLabels)
310 << abort(FatalError);
311
312 return -1;
313}
314
315
316// Get the other edge connected to pointi on facei.
317label otherEdge
318(
319 const triSurface& surf,
320 const label facei,
321 const label otherEdgeI,
322 const label pointi
323)
324{
325 const labelList& fEdges = surf.faceEdges()[facei];
326
327 forAll(fEdges, i)
328 {
329 label edgeI = fEdges[i];
330
331 const edge& e = surf.edges()[edgeI];
332
333 if
334 (
335 edgeI != otherEdgeI
336 && (
337 e.start() == pointi
338 || e.end() == pointi
339 )
340 )
341 {
342 return edgeI;
343 }
344 }
345
347 << " verts:" << surf.localPoints()[facei]
348 << " connected to point " << pointi
349 << " faceEdges:" << UIndirectList<edge>(surf.edges(), fEdges)
350 << abort(FatalError);
351
352 return -1;
353}
354
355
356// Starting from startPoint on startEdge on startFace walk along border
357// and insert faces along the way. Walk keeps always one point or one edge
358// on the border line.
359void walkSplitLine
360(
361 const triSurface& surf,
362 const boolList& borderEdge,
363 const labelList& borderPoint,
364
365 const label startFacei,
366 const label startEdgeI, // is border edge
367 const label startPointi, // is border point
368
369 Map<label>& faceToEdge,
371)
372{
373 label facei = startFacei;
374 label edgeI = startEdgeI;
375 label pointi = startPointi;
376
377 do
378 {
379 //
380 // Stick to pointi and walk face-edge-face until back on border edge.
381 //
382
383 do
384 {
385 // Cross face to next edge.
386 edgeI = otherEdge(surf, facei, edgeI, pointi);
387
388 if (borderEdge[edgeI])
389 {
390 if (!faceToEdge.insert(facei, edgeI))
391 {
392 // Was already visited.
393 return;
394 }
395 else
396 {
397 // First visit to this borderEdge. We're back on borderline.
398 break;
399 }
400 }
401 else if (!faceToPoint.insert(facei, pointi))
402 {
403 // Was already visited.
404 return;
405 }
406
407 // Cross edge to other face
408 const labelList& eFaces = surf.edgeFaces()[edgeI];
409
410 if (eFaces.size() != 2)
411 {
413 << "Can only handle edges with 2 or 4 edges for now."
414 << abort(FatalError);
415 }
416
417 if (eFaces[0] == facei)
418 {
419 facei = eFaces[1];
420 }
421 else if (eFaces[1] == facei)
422 {
423 facei = eFaces[0];
424 }
425 else
426 {
427 FatalErrorInFunction << abort(FatalError);
428 }
429 }
430 while (true);
431
432
433 //
434 // Back on border edge. Cross to other point on edge.
435 //
436
437 pointi = surf.edges()[edgeI].otherVertex(pointi);
438
439 if (borderPoint[pointi] == -1)
440 {
441 return;
442 }
443
444 }
445 while (true);
446}
447
448
449// Find second face which is from same surface i.e. has points on the
450// shared edge in reverse order.
451label sharedFace
452(
453 const triSurface& surf,
454 const label firstFacei,
455 const label sharedEdgeI
456)
457{
458 // Find ordering of face points in edge.
459
460 const edge& e = surf.edges()[sharedEdgeI];
461
462 const triSurface::face_type& f = surf.localFaces()[firstFacei];
463
464 label startIndex = f.find(e.start());
465
466 // points in face in same order as edge
467 const bool edgeOrder = (f[f.fcIndex(startIndex)] == e.end());
468
469 // Get faces using edge in sorted order. (sorted such that walking
470 // around them in anti-clockwise order corresponds to edge vector
471 // acc. to right-hand rule)
472 const labelList& eFaces = surf.sortedEdgeFaces()[sharedEdgeI];
473
474 // Get position of face in sorted edge faces
475 const label faceIndex = eFaces.find(firstFacei);
476
477 if (edgeOrder)
478 {
479 // Get face before firstFacei
480 return eFaces[eFaces.rcIndex(faceIndex)];
481 }
482 else
483 {
484 // Get face after firstFacei
485 return eFaces[eFaces.fcIndex(faceIndex)];
486 }
487}
488
489
490// Calculate (inward pointing) normals on edges shared by faces in
491// faceToEdge and averages them to pointNormals.
492void calcPointVecs
493(
494 const triSurface& surf,
495 const Map<label>& faceToEdge,
496 const Map<label>& faceToPoint,
497 vectorField& borderPointVec
498)
499{
500 const labelListList& sortedEdgeFaces = surf.sortedEdgeFaces();
501 const edgeList& edges = surf.edges();
502 const pointField& points = surf.localPoints();
503
504 boolList edgeDone(surf.nEdges(), false);
505
506 forAllConstIters(faceToEdge, iter)
507 {
508 const label edgeI = iter.val();
509
510 if (!edgeDone[edgeI])
511 {
512 edgeDone[edgeI] = true;
513
514 // Get the two connected faces in sorted order
515 // Note: should have stored this when walking ...
516
517 label face0I = -1;
518 label face1I = -1;
519
520 const labelList& eFaces = sortedEdgeFaces[edgeI];
521
522 forAll(eFaces, i)
523 {
524 label facei = eFaces[i];
525
526 if (faceToEdge.found(facei))
527 {
528 if (face0I == -1)
529 {
530 face0I = facei;
531 }
532 else if (face1I == -1)
533 {
534 face1I = facei;
535
536 break;
537 }
538 }
539 }
540
541 if (face0I == -1 && face1I == -1)
542 {
543 Info<< "Writing surface to errorSurf.obj" << endl;
544
545 surf.write("errorSurf.obj");
546
548 << "Cannot find two faces using border edge " << edgeI
549 << " verts:" << edges[edgeI]
550 << " eFaces:" << eFaces << endl
551 << "face0I:" << face0I
552 << " face1I:" << face1I << nl
553 << "faceToEdge:" << faceToEdge << nl
554 << "faceToPoint:" << faceToPoint
555 << "Written surface to errorSurf.obj"
556 << abort(FatalError);
557 }
558
559 // Now we have edge and the two faces in counter-clockwise order
560 // as seen from edge vector. Calculate normal.
561 const edge& e = edges[edgeI];
562 vector eVec = e.vec(points);
563
564 // Determine vector as average of the vectors in the two faces.
565 // If there is only one face available use only one vector.
566 vector midVec(Zero);
567
568 if (face0I != -1)
569 {
570 label v0 = triSurfaceTools::oppositeVertex(surf, face0I, edgeI);
571 const vector e0 =
573 (
574 (points[v0] - points[e.start()]) ^ eVec
575 );
576
577 midVec = e0;
578 }
579
580 if (face1I != -1)
581 {
582 label v1 = triSurfaceTools::oppositeVertex(surf, face1I, edgeI);
583 const vector e1 =
585 (
586 (points[e.start()] - points[v1]) ^ eVec
587 );
588
589 midVec += e1;
590 }
591
592 scalar magMidVec = mag(midVec);
593
594 if (magMidVec > SMALL)
595 {
596 midVec /= magMidVec;
597
598 // Average to pointVec
599 borderPointVec[e.start()] += midVec;
600 borderPointVec[e.end()] += midVec;
601 }
602 }
603 }
604}
605
606
607// Renumbers vertices (of triangles in faceToEdge) of which the pointMap is
608// not -1.
609void renumberFaces
610(
611 const triSurface& surf,
612 const labelList& pointMap,
613 const Map<label>& faceToEdge,
615)
616{
617 forAllConstIters(faceToEdge, iter)
618 {
619 const label facei = iter.key();
620 const triSurface::face_type& f = surf.localFaces()[facei];
621
622 forAll(f, fp)
623 {
624 if (pointMap[f[fp]] != -1)
625 {
626 newTris[facei][fp] = pointMap[f[fp]];
627 }
628 }
629 }
630}
631
632
633// Split all borderEdges that don't have borderPoint. Return true if split
634// anything.
635bool splitBorderEdges
636(
637 triSurface& surf,
638 const boolList& borderEdge,
639 const labelList& borderPoint
640)
641{
642 labelList edgesToBeSplit(surf.nEdges());
643 label nSplit = 0;
644
645 forAll(borderEdge, edgeI)
646 {
647 if (borderEdge[edgeI])
648 {
649 const edge& e = surf.edges()[edgeI];
650
651 if (borderPoint[e.start()] == -1 && borderPoint[e.end()] == -1)
652 {
653 // None of the points of the edge is borderPoint. Split edge
654 // to introduce border point.
655 edgesToBeSplit[nSplit++] = edgeI;
656 }
657 }
658 }
659 edgesToBeSplit.setSize(nSplit);
660
661 if (nSplit > 0)
662 {
663 Info<< "Splitting surface along " << nSplit << " borderEdges that don't"
664 << " neighbour other borderEdges" << nl << endl;
665
666 surf = triSurfaceTools::greenRefine(surf, edgesToBeSplit);
667
668 return true;
669 }
670
671 Info<< "No edges to be split" <<nl << endl;
672 return false;
673}
674
675
676
677int main(int argc, char *argv[])
678{
679 argList::addNote
680 (
681 "Split multiply connected surface edges by duplicating points"
682 );
683 argList::noParallel();
684 argList::addArgument("input", "The input surface file");
685 argList::addArgument("output", "The output surface file");
686 argList::addBoolOption
687 (
688 "debug",
689 "Add debugging output"
690 );
691
692 argList args(argc, argv);
693
694 const auto inSurfName = args.get<fileName>(1);
695 const auto outSurfName = args.get<fileName>(2);
696 const bool debug = args.found("debug");
697
698 Info<< "Reading surface from " << inSurfName << endl;
699
700 triSurface surf(inSurfName);
701
702 // Make sure sortedEdgeFaces is calculated correctly
703 testSortedEdgeFaces(surf);
704
705 // Get all quad connected edges. These are seen as borders when walking.
706 boolList borderEdge(surf.nEdges(), false);
707 markBorderEdges(debug, surf, borderEdge);
708
709 // Points on two sides connected to borderEdges are called
710 // borderPoints and will be duplicated. borderPoint contains label
711 // of newly introduced vertex.
712 labelList borderPoint(surf.nPoints(), -1);
713 markBorderPoints(debug, surf, borderEdge, borderPoint);
714
715 // Split edges where there would be no borderPoint to duplicate.
716 splitBorderEdges(surf, borderEdge, borderPoint);
717
718
719 Info<< "Writing split surface to " << outSurfName << nl << endl;
720 surf.write(outSurfName);
721 Info<< "Finished writing surface to " << outSurfName << nl << endl;
722
723
724 // Last iteration values.
725 label nOldBorderEdges = -1;
726 label nOldBorderPoints = -1;
727
728 label iteration = 0;
729
730 do
731 {
732 // Redo borderEdge/borderPoint calculation.
733 boolList borderEdge(surf.nEdges(), false);
734
735 label nBorderEdges = markBorderEdges(debug, surf, borderEdge);
736
737 if (nBorderEdges == 0)
738 {
739 Info<< "Found no border edges. Exiting." << nl << nl;
740
741 break;
742 }
743
744 // Label of newly introduced duplicate.
745 labelList borderPoint(surf.nPoints(), -1);
746
747 label nBorderPoints =
748 markBorderPoints
749 (
750 debug,
751 surf,
752 borderEdge,
753 borderPoint
754 );
755
756 if (nBorderPoints == 0)
757 {
758 Info<< "Found no border points. Exiting." << nl << nl;
759
760 break;
761 }
762
763 Info<< "Found:\n"
764 << " border edges :" << nBorderEdges << nl
765 << " border points:" << nBorderPoints << nl
766 << endl;
767
768 if
769 (
770 nBorderPoints == nOldBorderPoints
771 && nBorderEdges == nOldBorderEdges
772 )
773 {
774 Info<< "Stopping since number of border edges and point is same"
775 << " as in previous iteration" << nl << endl;
776
777 break;
778 }
779
780 //
781 // Define splitLine as a series of connected borderEdges. Find start
782 // of one (as edge and point on edge)
783 //
784
785 label startEdgeI = -1;
786 label startPointi = -1;
787
788 forAll(borderEdge, edgeI)
789 {
790 if (borderEdge[edgeI])
791 {
792 const edge& e = surf.edges()[edgeI];
793
794 if ((borderPoint[e[0]] != -1) && (borderPoint[e[1]] == -1))
795 {
796 startEdgeI = edgeI;
797 startPointi = e[0];
798
799 break;
800 }
801 else if ((borderPoint[e[0]] == -1) && (borderPoint[e[1]] != -1))
802 {
803 startEdgeI = edgeI;
804 startPointi = e[1];
805
806 break;
807 }
808 }
809 }
810
811 if (startEdgeI == -1)
812 {
813 Info<< "Cannot find starting point of splitLine\n" << endl;
814
815 break;
816 }
817
818 // Pick any face using edge to start from.
819 const labelList& eFaces = surf.edgeFaces()[startEdgeI];
820
821 label firstFacei = eFaces[0];
822
823 // Find second face which is from same surface i.e. has outwards
824 // pointing normal as well (actually bit more complex than this)
825 label secondFacei = sharedFace(surf, firstFacei, startEdgeI);
826
827 Info<< "Starting local walk from:" << endl
828 << " edge :" << startEdgeI << endl
829 << " point:" << startPointi << endl
830 << " face0:" << firstFacei << endl
831 << " face1:" << secondFacei << endl
832 << endl;
833
834 // From face on border edge to edge.
835 Map<label> faceToEdge(2*nBorderEdges);
836 // From face connected to border point (but not border edge) to point.
837 Map<label> faceToPoint(2*nBorderPoints);
838
839 faceToEdge.insert(firstFacei, startEdgeI);
840
841 walkSplitLine
842 (
843 surf,
844 borderEdge,
845 borderPoint,
846
847 firstFacei,
848 startEdgeI,
849 startPointi,
850
851 faceToEdge,
853 );
854
855 faceToEdge.insert(secondFacei, startEdgeI);
856
857 walkSplitLine
858 (
859 surf,
860 borderEdge,
861 borderPoint,
862
863 secondFacei,
864 startEdgeI,
865 startPointi,
866
867 faceToEdge,
869 );
870
871 Info<< "Finished local walk and visited" << nl
872 << " border edges :" << faceToEdge.size() << nl
873 << " border points(but not edges):" << faceToPoint.size() << nl
874 << endl;
875
876 if (debug)
877 {
878 dumpFaces("faceToEdge.obj", surf, faceToEdge);
879 dumpFaces("faceToPoint.obj", surf, faceToPoint);
880 }
881
882 //
883 // Create coordinates for borderPoints by duplicating the existing
884 // point and then slightly shifting it inwards. To determine the
885 // inwards direction get the average normal of both connectedFaces on
886 // the edge and then interpolate this to the (border)point.
887 //
888
889 vectorField borderPointVec(surf.nPoints(), vector(GREAT, GREAT, GREAT));
890
891 calcPointVecs(surf, faceToEdge, faceToPoint, borderPointVec);
892
893
894 // New position. Start off from copy of old points.
895 pointField newPoints(surf.localPoints());
896 newPoints.setSize(newPoints.size() + nBorderPoints);
897
898 forAll(borderPoint, pointi)
899 {
900 label newPointi = borderPoint[pointi];
901
902 if (newPointi != -1)
903 {
904 scalar minLen = minEdgeLen(surf, pointi);
905
906 const vector n = normalised(borderPointVec[pointi]);
907
908 newPoints[newPointi] = newPoints[pointi] + 0.1 * minLen * n;
909 }
910 }
911
912
913 //
914 // Renumber all faces in connectedFaces
915 //
916
917 // Start off from copy of faces.
918 List<labelledTri> newTris(surf.size());
919
920 forAll(surf, facei)
921 {
922 newTris[facei] = surf.localFaces()[facei];
923 newTris[facei].region() = surf[facei].region();
924 }
925
926 // Renumber all faces in faceToEdge
927 renumberFaces(surf, borderPoint, faceToEdge, newTris);
928 // Renumber all faces in faceToPoint
929 renumberFaces(surf, borderPoint, faceToPoint, newTris);
930
931
932 // Check if faces use unmoved points.
933 forAll(newTris, facei)
934 {
935 const triSurface::face_type& f = newTris[facei];
936
937 forAll(f, fp)
938 {
939 const point& pt = newPoints[f[fp]];
940
941 if (mag(pt) >= GREAT/2)
942 {
943 Info<< "newTri:" << facei << " verts:" << f
944 << " vert:" << f[fp] << " point:" << pt << endl;
945 }
946 }
947 }
948
949
950 surf = triSurface(newTris, surf.patches(), newPoints);
951
952 if (debug || (iteration != 0 && (iteration % 20) == 0))
953 {
954 Info<< "Writing surface to " << outSurfName << nl << endl;
955
956 surf.write(outSurfName);
957
958 Info<< "Finished writing surface to " << outSurfName << nl << endl;
959 }
960
961 // Save prev iteration values.
962 nOldBorderEdges = nBorderEdges;
963 nOldBorderPoints = nBorderPoints;
964
965 iteration++;
966 }
967 while (true);
968
969 Info<< "Writing final surface to " << outSurfName << nl << endl;
970
971 surf.write(outSurfName);
972
973 Info<< "End\n" << endl;
974
975 return 0;
976}
977
978
979// ************************************************************************* //
Various functions to operate on Lists.
label n
bool found(const Key &key) const
Return true if hashed entry is found in table.
Definition: HashTableI.H:100
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
Definition: HashTableI.H:180
label size() const noexcept
The number of elements in table.
Definition: HashTableI.H:52
void setSize(const label n)
Alias for resize()
Definition: List.H:218
A HashTable to objects of type <T> with a label key.
Definition: Map.H:60
Output to file stream, using an OSstream.
Definition: OFstream.H:57
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
label nEdges() const
Number of edges in patch.
label nPoints() const
Number of points supporting patch faces.
const labelListList & pointEdges() const
Return point-edge addressing.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
const Field< point_type > & localPoints() const
Return pointField of points in patch.
const labelListList & edgeFaces() const
Return edge-face addressing.
const labelListList & faceEdges() const
Return face-edge addressing.
const List< face_type > & localFaces() const
Return patch faces addressing into local point list.
A List with indirect addressing. Like IndirectList but does not store addressing.
Definition: IndirectList.H:79
bool found(const T &val, label pos=0) const
True if the value if found in the list.
Definition: UListI.H:265
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
Definition: UList.C:212
label rcIndex(const label i) const noexcept
Definition: UListI.H:67
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
label fcIndex(const label i) const noexcept
Definition: UListI.H:60
const Cmpt & z() const
Access to the vector z component.
Definition: VectorI.H:85
const Cmpt & y() const
Access to the vector y component.
Definition: VectorI.H:79
const Cmpt & x() const
Access to the vector x component.
Definition: VectorI.H:73
Extract command arguments and options from the supplied argc and argv parameters.
Definition: argList.H:124
T get(const label index) const
Get a value from the argument at index.
Definition: argListI.H:278
bool found(const word &optName) const
Return true if the named option is found.
Definition: argListI.H:178
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
A topoSetPointSource to select all points based on usage in given faceSet(s).
Definition: faceToPoint.H:177
A class for handling file names.
Definition: fileName.H:76
A triFace with additional (region) index.
Definition: labelledTri.H:60
Triangulated surface description with patch information.
Definition: triSurface.H:79
void write(Ostream &os) const
Write to Ostream in simple OpenFOAM format.
Definition: triSurfaceIO.C:336
const geometricSurfacePatchList & patches() const noexcept
Definition: triSurface.H:399
const labelListList & sortedEdgeFaces() const
Return edge-face addressing sorted (for edges with more than.
Definition: triSurface.C:590
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
OBJstream os(runTime.globalPath()/outputName)
const pointField & points
label nPoints
int debug
Static debugging option.
label findEdge(const edgeList &edges, const labelList &candidates, const label v0, const label v1)
Return edge among candidates that uses the two vertices.
Definition: meshTools.C:359
label otherEdge(const primitiveMesh &mesh, const labelList &edgeLabels, const label thisEdgeI, const label thisVertI)
Return label of other edge (out of candidates edgeLabels)
Definition: meshTools.C:521
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of a point.
Definition: meshTools.C:203
Namespace for OpenFOAM.
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Definition: quaternionI.H:680
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManip< error > abort(error &err)
Definition: errorManip.H:144
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
label newPointi
Definition: readKivaGrid.H:496
labelList f(nPoints)
Foam::argList args(argc, argv)
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Definition: stdFoam.H:278