boundaryCutter.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 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
27\*---------------------------------------------------------------------------*/
28
29#include "boundaryCutter.H"
30#include "polyMesh.H"
31#include "polyTopoChange.H"
32#include "polyAddCell.H"
33#include "polyAddFace.H"
34#include "polyAddPoint.H"
35#include "polyModifyFace.H"
36#include "polyModifyPoint.H"
37#include "mapPolyMesh.H"
38#include "meshTools.H"
39
40// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41
42namespace Foam
43{
45}
46
47
48// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49
50void Foam::boundaryCutter::getFaceInfo
51(
52 const label facei,
53 label& patchID,
54 label& zoneID,
55 label& zoneFlip
56) const
57{
58 patchID = -1;
59
60 if (!mesh_.isInternalFace(facei))
61 {
62 patchID = mesh_.boundaryMesh().whichPatch(facei);
63 }
64
65 zoneID = mesh_.faceZones().whichZone(facei);
66
67 zoneFlip = false;
68
69 if (zoneID >= 0)
70 {
71 const faceZone& fZone = mesh_.faceZones()[zoneID];
72
73 zoneFlip = fZone.flipMap()[fZone.whichFace(facei)];
74 }
75}
76
77
78// Adds additional vertices (from edge cutting) to face. Used for faces which
79// are not split but still might use edge that has been cut.
80Foam::face Foam::boundaryCutter::addEdgeCutsToFace
81(
82 const label facei,
83 const Map<labelList>& edgeToAddedPoints
84) const
85{
86 const edgeList& edges = mesh_.edges();
87 const face& f = mesh_.faces()[facei];
88 const labelList& fEdges = mesh_.faceEdges()[facei];
89
90 // Storage for face
91 DynamicList<label> newFace(2 * f.size());
92
93 forAll(f, fp)
94 {
95 // Duplicate face vertex .
96 newFace.append(f[fp]);
97
98 // Check if edge has been cut.
99 label v1 = f.nextLabel(fp);
100
101 label edgeI = meshTools::findEdge(edges, fEdges, f[fp], v1);
102
103 const auto fnd = edgeToAddedPoints.cfind(edgeI);
104
105 if (fnd.found())
106 {
107 // edge has been cut. Introduce new vertices. Check order.
108 const labelList& addedPoints = fnd();
109
110 if (edges[edgeI].start() == f[fp])
111 {
112 // Introduce in same order.
113 forAll(addedPoints, i)
114 {
115 newFace.append(addedPoints[i]);
116 }
117 }
118 else
119 {
120 // Introduce in opposite order.
121 forAllReverse(addedPoints, i)
122 {
123 newFace.append(addedPoints[i]);
124 }
125 }
126 }
127 }
128
129 face returnFace;
130 returnFace.transfer(newFace);
131
132 if (debug)
133 {
134 Pout<< "addEdgeCutsToFace:" << nl
135 << " from : " << f << nl
136 << " to : " << returnFace << endl;
137 }
138
139 return returnFace;
140}
141
142
143void Foam::boundaryCutter::addFace
144(
145 const label facei,
146 const face& newFace,
147
148 bool& modifiedFace, // have we already 'used' facei
149 polyTopoChange& meshMod
150) const
151{
152 // Information about old face
153 label patchID, zoneID, zoneFlip;
154 getFaceInfo(facei, patchID, zoneID, zoneFlip);
155 label own = mesh_.faceOwner()[facei];
156 label masterPoint = mesh_.faces()[facei][0];
157
158 if (!modifiedFace)
159 {
160 meshMod.setAction
161 (
162 polyModifyFace
163 (
164 newFace, // face
165 facei,
166 own, // owner
167 -1, // neighbour
168 false, // flux flip
169 patchID, // patch for face
170 false, // remove from zone
171 zoneID, // zone for face
172 zoneFlip // face zone flip
173 )
174 );
175
176 modifiedFace = true;
177 }
178 else
179 {
180 meshMod.setAction
181 (
182 polyAddFace
183 (
184 newFace, // face
185 own, // owner
186 -1, // neighbour
187 masterPoint, // master point
188 -1, // master edge
189 -1, // master face for addition
190 false, // flux flip
191 patchID, // patch for face
192 zoneID, // zone for face
193 zoneFlip // face zone flip
194 )
195 );
196 }
197}
198
199
200
201// Splits a face using the cut edges and modified points
202bool Foam::boundaryCutter::splitFace
203(
204 const label facei,
205 const Map<point>& pointToPos,
206 const Map<labelList>& edgeToAddedPoints,
207 polyTopoChange& meshMod
208) const
209{
210 const edgeList& edges = mesh_.edges();
211 const face& f = mesh_.faces()[facei];
212 const labelList& fEdges = mesh_.faceEdges()[facei];
213
214 // Count number of split edges and total number of splits.
215 label nSplitEdges = 0;
216 label nModPoints = 0;
217 label nTotalSplits = 0;
218
219 forAll(f, fp)
220 {
221 if (pointToPos.found(f[fp]))
222 {
223 nModPoints++;
224 nTotalSplits++;
225 }
226
227 // Check if edge has been cut.
228 label nextV = f.nextLabel(fp);
229
230 label edgeI = meshTools::findEdge(edges, fEdges, f[fp], nextV);
231
232 const auto fnd = edgeToAddedPoints.cfind(edgeI);
233
234 if (fnd.found())
235 {
236 nSplitEdges++;
237 nTotalSplits += fnd().size();
238 }
239 }
240
241 if (debug)
242 {
243 Pout<< "Face:" << facei
244 << " nModPoints:" << nModPoints
245 << " nSplitEdges:" << nSplitEdges
246 << " nTotalSplits:" << nTotalSplits << endl;
247 }
248
249 if (nSplitEdges == 0 && nModPoints == 0)
250 {
252 << " nSplitEdges:" << nSplitEdges
253 << " nTotalSplits:" << nTotalSplits
254 << abort(FatalError);
255 return false;
256 }
257 else if (nSplitEdges + nModPoints == 1)
258 {
259 // single or multiple cuts on a single edge or single modified point
260 // Don't cut and let caller handle this.
261 Warning << "Face " << facei << " has only one edge cut " << endl;
262 return false;
263 }
264 else
265 {
266 // So guaranteed to have two edges cut or points modified. Split face:
267 // - find starting cut
268 // - walk to next cut. Make face
269 // - loop until face done.
270
271 // Information about old face
272 label patchID, zoneID, zoneFlip;
273 getFaceInfo(facei, patchID, zoneID, zoneFlip);
274
275 // Get face with new points on cut edges for ease of looping
276 face extendedFace(addEdgeCutsToFace(facei, edgeToAddedPoints));
277
278 // Find first added point. This is the starting vertex for splitting.
279 label startFp = -1;
280
281 forAll(extendedFace, fp)
282 {
283 if (extendedFace[fp] >= mesh_.nPoints())
284 {
285 startFp = fp;
286 break;
287 }
288 }
289
290 if (startFp == -1)
291 {
292 // No added point. Maybe there is a modified point?
293 forAll(extendedFace, fp)
294 {
295 if (pointToPos.found(extendedFace[fp]))
296 {
297 startFp = fp;
298 break;
299 }
300 }
301 }
302
303 if (startFp == -1)
304 {
306 << "Problem" << abort(FatalError);
307 }
308
309 // Have we already modified existing face (first face gets done
310 // as modification; all following ones as polyAddFace)
311 bool modifiedFace = false;
312
313 // Example face:
314 // +--+
315 // / |
316 // / |
317 // + +
318 // \ |
319 // \ |
320 // +--+
321 //
322 // Needs to get split into:
323 // - three 'side' faces a,b,c
324 // - one middle face d
325 // +--+
326 // /|\A|
327 // / | \|
328 // + C|D +
329 // \ | /|
330 // \|/B|
331 // +--+
332
333
334 // Storage for new face
335 DynamicList<label> newFace(extendedFace.size());
336
337 label fp = startFp;
338
339 forAll(extendedFace, i)
340 {
341 label pointi = extendedFace[fp];
342
343 newFace.append(pointi);
344
345 if
346 (
347 newFace.size() > 2
348 && (
349 pointi >= mesh_.nPoints()
350 || pointToPos.found(pointi)
351 )
352 )
353 {
354 // Enough vertices to create a face from.
355 face tmpFace;
356 tmpFace.transfer(newFace);
357
358 // Add face tmpFace
359 addFace(facei, tmpFace, modifiedFace, meshMod);
360
361 // Starting point is also the starting point for the new face
362 newFace.append(extendedFace[startFp]);
363 newFace.append(extendedFace[fp]);
364 }
365
366 fp = (fp+1) % extendedFace.size();
367 }
368
369 // Check final face.
370 if (newFace.size() > 2)
371 {
372 // Enough vertices to create a face from.
373 face tmpFace;
374 tmpFace.transfer(newFace);
375
376 // Add face tmpFace
377 addFace(facei, tmpFace, modifiedFace, meshMod);
378 }
379
380 // Split something
381 return true;
382 }
383}
384
385
386// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
387
389:
390 mesh_(mesh),
391 edgeAddedPoints_(),
392 faceAddedPoint_()
393{}
394
395
396// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
397
399(
400 const Map<point>& pointToPos,
401 const Map<List<point>>& edgeToCuts,
402 const Map<labelPair>& faceToSplit,
403 const Map<point>& faceToFeaturePoint,
404 polyTopoChange& meshMod
405)
406{
407 // Clear and size maps here since mesh size will change.
408 edgeAddedPoints_.clear();
409
410 faceAddedPoint_.clear();
411 faceAddedPoint_.resize(faceToFeaturePoint.size());
412
413
414 //
415 // Points that just need to be moved
416 // Note: could just as well be handled outside of setRefinement.
417 //
418
419 forAllConstIters(pointToPos, iter)
420 {
421 meshMod.setAction
422 (
424 (
425 iter.key(), // point
426 iter.val(), // position
427 false, // no zone
428 -1, // zone for point
429 true // supports a cell
430 )
431 );
432 }
433
434
435 //
436 // Add new points along cut edges.
437 //
438
439 // Map from edge label to sorted list of points
440 Map<labelList> edgeToAddedPoints(edgeToCuts.size());
441
442 forAllConstIters(edgeToCuts, iter)
443 {
444 const label edgeI = iter.key();
445 const List<point>& cuts = iter.val();
446
447 const edge& e = mesh_.edges()[edgeI];
448
449 // Sorted (from start to end) list of cuts on edge
450
451 forAll(cuts, cutI)
452 {
453 // point on feature to move to
454 const point& featurePoint = cuts[cutI];
455
456 label addedPointi =
457 meshMod.setAction
458 (
460 (
461 featurePoint, // point
462 e.start(), // master point
463 -1, // zone for point
464 true // supports a cell
465 )
466 );
467
468 auto fnd = edgeToAddedPoints.find(edgeI);
469
470 if (fnd.found())
471 {
472 labelList& addedPoints = fnd();
473
474 label sz = addedPoints.size();
475 addedPoints.setSize(sz+1);
476 addedPoints[sz] = addedPointi;
477 }
478 else
479 {
480 edgeToAddedPoints.insert(edgeI, labelList(1, addedPointi));
481 }
482
483 if (debug)
484 {
485 Pout<< "Added point " << addedPointi << " for edge " << edgeI
486 << " with cuts:" << edgeToAddedPoints[edgeI] << endl;
487 }
488 }
489 }
490
491
492 //
493 // Introduce feature points.
494 //
495
496 forAllConstIters(faceToFeaturePoint, iter)
497 {
498 const label facei = iter.key();
499
500 const face& f = mesh_.faces()[facei];
501
502 if (faceToSplit.found(facei))
503 {
505 << "Face " << facei << " vertices " << f
506 << " is both marked for face-centre decomposition and"
507 << " diagonal splitting."
508 << abort(FatalError);
509 }
510
511 if (mesh_.isInternalFace(facei))
512 {
514 << "Face " << facei << " vertices " << f
515 << " is not an external face. Cannot split it"
516 << abort(FatalError);
517 }
518
519 label addedPointi =
520 meshMod.setAction
521 (
523 (
524 iter.val(), // point
525 f[0], // master point
526 -1, // zone for point
527 true // supports a cell
528 )
529 );
530
531 faceAddedPoint_.insert(facei, addedPointi);
532
533 if (debug)
534 {
535 Pout<< "Added point " << addedPointi << " for feature point "
536 << iter.val() << " on face " << facei << " with centre "
537 << mesh_.faceCentres()[facei] << endl;
538 }
539 }
540
541
542 //
543 // Split or retriangulate faces
544 //
545
546
547 // Maintain whether face has been updated (for -split edges
548 // -new owner/neighbour)
549 boolList faceUptodate(mesh_.nFaces(), false);
550
551
552 // Triangulate faces containing feature points
553 forAllConstIters(faceAddedPoint_, iter)
554 {
555 const label facei = iter.key();
556 const label addedPointi = iter.val();
557
558 // Get face with new points on cut edges.
559 face newFace(addEdgeCutsToFace(facei, edgeToAddedPoints));
560
561 // Information about old face
562 label patchID, zoneID, zoneFlip;
563 getFaceInfo(facei, patchID, zoneID, zoneFlip);
564 label own = mesh_.faceOwner()[facei];
565 label masterPoint = mesh_.faces()[facei][0];
566
567 // Triangulate face around mid point
568
569 face tri(3);
570
571 forAll(newFace, fp)
572 {
573 label nextV = newFace.nextLabel(fp);
574
575 tri[0] = newFace[fp];
576 tri[1] = nextV;
577 tri[2] = addedPointi;
578
579 if (fp == 0)
580 {
581 // Modify the existing face.
582 meshMod.setAction
583 (
585 (
586 tri, // face
587 facei,
588 own, // owner
589 -1, // neighbour
590 false, // flux flip
591 patchID, // patch for face
592 false, // remove from zone
593 zoneID, // zone for face
594 zoneFlip // face zone flip
595 )
596 );
597 }
598 else
599 {
600 // Add additional faces
601 meshMod.setAction
602 (
604 (
605 tri, // face
606 own, // owner
607 -1, // neighbour
608 masterPoint, // master point
609 -1, // master edge
610 -1, // master face for addition
611 false, // flux flip
612 patchID, // patch for face
613 zoneID, // zone for face
614 zoneFlip // face zone flip
615 )
616 );
617 }
618 }
619
620 faceUptodate[facei] = true;
621 }
622
623
624 // Diagonally split faces
625 forAllConstIters(faceToSplit, iter)
626 {
627 const label facei = iter.key();
628
629 const face& f = mesh_.faces()[facei];
630
631 if (faceAddedPoint_.found(facei))
632 {
634 << "Face " << facei << " vertices " << f
635 << " is both marked for face-centre decomposition and"
636 << " diagonal splitting."
637 << abort(FatalError);
638 }
639
640
641 // Get face with new points on cut edges.
642 face newFace(addEdgeCutsToFace(facei, edgeToAddedPoints));
643
644 // Information about old face
645 label patchID, zoneID, zoneFlip;
646 getFaceInfo(facei, patchID, zoneID, zoneFlip);
647 label own = mesh_.faceOwner()[facei];
648 label masterPoint = mesh_.faces()[facei][0];
649
650 // Split face from one side of diagonal to other.
651 const labelPair& diag = iter.val();
652
653 label fp0 = newFace.find(f[diag[0]]);
654 label fp1 = newFace.find(f[diag[1]]);
655
656 if (fp0 == -1 || fp1 == -1 || fp0 == fp1)
657 {
659 << "Problem : Face " << facei << " vertices " << f
660 << " newFace:" << newFace << " diagonal:" << f[diag[0]]
661 << ' ' << f[diag[1]]
662 << abort(FatalError);
663 }
664
665 // Replace existing face by newFace from fp0 to fp1 and add new one
666 // from fp1 to fp0.
667
668 DynamicList<label> newVerts(newFace.size());
669
670 // Get vertices from fp0 to (and including) fp1
671 label fp = fp0;
672
673 do
674 {
675 newVerts.append(newFace[fp]);
676
677 fp = (fp == newFace.size()-1 ? 0 : fp+1);
678
679 } while (fp != fp1);
680
681 newVerts.append(newFace[fp1]);
682
683
684 // Modify the existing face.
685 meshMod.setAction
686 (
688 (
689 face(newVerts.shrink()), // face
690 facei,
691 own, // owner
692 -1, // neighbour
693 false, // flux flip
694 patchID, // patch for face
695 false, // remove from zone
696 zoneID, // zone for face
697 zoneFlip // face zone flip
698 )
699 );
700
701
702 newVerts.clear();
703
704 // Get vertices from fp1 to (and including) fp0
705
706 do
707 {
708 newVerts.append(newFace[fp]);
709
710 fp = (fp == newFace.size()-1 ? 0 : fp+1);
711
712 } while (fp != fp0);
713
714 newVerts.append(newFace[fp0]);
715
716 // Add additional face
717 meshMod.setAction
718 (
720 (
721 face(newVerts.shrink()), // face
722 own, // owner
723 -1, // neighbour
724 masterPoint, // master point
725 -1, // master edge
726 -1, // master face for addition
727 false, // flux flip
728 patchID, // patch for face
729 zoneID, // zone for face
730 zoneFlip // face zone flip
731 )
732 );
733
734 faceUptodate[facei] = true;
735 }
736
737
738 // Split external faces without feature point but using cut edges.
739 // Does right handed walk but not really.
740 forAllConstIters(edgeToAddedPoints, iter)
741 {
742 const label edgeI = iter.key();
743
744 const labelList& eFaces = mesh_.edgeFaces()[edgeI];
745
746 forAll(eFaces, i)
747 {
748 label facei = eFaces[i];
749
750 if (!faceUptodate[facei] && !mesh_.isInternalFace(facei))
751 {
752 // Is external face so split
753 if (splitFace(facei, pointToPos, edgeToAddedPoints, meshMod))
754 {
755 // Successful split
756 faceUptodate[facei] = true;
757 }
758 }
759 }
760 }
761
762
763 // Add cut edges (but don't split) any other faces using any cut edge.
764 // These can be external faces where splitFace hasn't cut them or
765 // internal faces.
766 forAllConstIters(edgeToAddedPoints, iter)
767 {
768 const label edgeI = iter.key();
769
770 const labelList& eFaces = mesh_.edgeFaces()[edgeI];
771
772 forAll(eFaces, i)
773 {
774 label facei = eFaces[i];
775
776 if (!faceUptodate[facei])
777 {
778 // Renumber face to include split edges.
779 face newFace(addEdgeCutsToFace(facei, edgeToAddedPoints));
780
781 const label own = mesh_.faceOwner()[facei];
782
783 label nei = -1;
784
785 if (mesh_.isInternalFace(facei))
786 {
787 nei = mesh_.faceNeighbour()[facei];
788 }
789
790 label patchID, zoneID, zoneFlip;
791 getFaceInfo(facei, patchID, zoneID, zoneFlip);
792
793 meshMod.setAction
794 (
796 (
797 newFace, // modified face
798 facei, // label of face being modified
799 own, // owner
800 nei, // neighbour
801 false, // face flip
802 patchID, // patch for face
803 false, // remove from zone
804 zoneID, // zone for face
805 zoneFlip // face flip in zone
806 )
807 );
808
809 faceUptodate[facei] = true;
810 }
811 }
812 }
813
814 // Convert edge to points storage from edge labels (not preserved)
815 // to point labels
816 edgeAddedPoints_.resize(edgeToCuts.size());
817
818 forAllConstIters(edgeToAddedPoints, iter)
819 {
820 edgeAddedPoints_.insert(mesh_.edges()[iter.key()], iter.val());
821 }
822}
823
824
826{
827 // Update stored labels for mesh change.
828
829 //
830 // Do faceToAddedPoint
831 //
832
833 {
834 // Create copy since we're deleting entries.
835 Map<label> newAddedPoints(faceAddedPoint_.size());
836
837 forAllConstIters(faceAddedPoint_, iter)
838 {
839 const label oldFacei = iter.key();
840 const label oldPointi = iter.val();
841
842 const label newFacei = morphMap.reverseFaceMap()[oldFacei];
843 const label newPointi = morphMap.reversePointMap()[oldPointi];
844
845 if (newFacei >= 0 && newPointi >= 0)
846 {
847 newAddedPoints.insert(newFacei, newPointi);
848 }
849 }
850
851 // Copy
852 faceAddedPoint_.transfer(newAddedPoints);
853 }
854
855
856 //
857 // Do edgeToAddedPoints
858 //
859
860
861 {
862 // Create copy since we're deleting entries
863 EdgeMap<labelList> newEdgeAddedPoints(edgeAddedPoints_.size());
864
865 forAllConstIters(edgeAddedPoints_, iter)
866 {
867 const edge& e = iter.key();
868 const labelList& addedPoints = iter.val();
869
870 const label newStart = morphMap.reversePointMap()[e.start()];
871 const label newEnd = morphMap.reversePointMap()[e.end()];
872
873 if (newStart >= 0 && newEnd >= 0)
874 {
875 labelList newAddedPoints(addedPoints.size());
876 label newI = 0;
877
878 forAll(addedPoints, i)
879 {
880 label newAddedPointi =
881 morphMap.reversePointMap()[addedPoints[i]];
882
883 if (newAddedPointi >= 0)
884 {
885 newAddedPoints[newI++] = newAddedPointi;
886 }
887 }
888 if (newI > 0)
889 {
890 newAddedPoints.setSize(newI);
891
892 edge newE = edge(newStart, newEnd);
893
894 newEdgeAddedPoints.insert(newE, newAddedPoints);
895 }
896 }
897 }
898
899 // Copy
900 edgeAddedPoints_.transfer(newEdgeAddedPoints);
901 }
902}
903
904
905// ************************************************************************* //
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:72
void clear() noexcept
Clear the addressed list, i.e. set the size to zero.
Definition: DynamicListI.H:391
DynamicList< T, SizeMin > & shrink()
Shrink the allocated space to the number of elements used.
Definition: DynamicListI.H:434
void append(const T &val)
Copy append an element to the end of this list.
Definition: DynamicListI.H:503
Map from edge (expressed as its endpoints) to value. For easier forward declaration it is currently i...
Definition: EdgeMap.H:54
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
iterator find(const Key &key)
Find and return an iterator set at the hashed entry.
Definition: HashTableI.H:114
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
void setSize(const label n)
Alias for resize()
Definition: List.H:218
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
A HashTable to objects of type <T> with a label key.
Definition: Map.H:60
iterator end() noexcept
Return an iterator to end traversing the UList.
Definition: UListI.H:350
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
Definition: UList.C:212
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
label whichZone(const label objectIndex) const
Given a global object index, return the zone it is in.
Definition: ZoneMesh.C:288
Does modifications to boundary faces.
void setRefinement(const Map< point > &pointToPos, const Map< List< point > > &edgeToCuts, const Map< labelPair > &faceToSplit, const Map< point > &faceToFeaturePoint, polyTopoChange &meshMod)
Do actual cutting with cut description. Inserts mesh changes.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
label nextLabel(const label i) const
Next vertex on face.
Definition: faceI.H:175
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:162
const labelList & reversePointMap() const
Reverse point map.
Definition: mapPolyMesh.H:469
const labelList & reverseFaceMap() const
Reverse face map.
Definition: mapPolyMesh.H:501
void updateMesh()
Update for new mesh topology.
A face addition data class. A face can be inflated either from a point or from another face and can e...
Definition: polyAddFace.H:56
Class containing data for point addition.
Definition: polyAddPoint.H:52
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
Definition: polyMesh.H:498
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:456
Class describing modification of a face.
Class describing modification of a point.
Direct mesh changes based on v1.3 polyTopoChange syntax.
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
bool isInternalFace(const label faceIndex) const noexcept
Return true if given face label is internal to the mesh.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
const labelIOList & zoneID
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
Namespace for OpenFOAM.
List< label > labelList
A List of labels.
Definition: List.H:66
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
errorManip< error > abort(error &err)
Definition: errorManip.H:144
error FatalError
List< edge > edgeList
A List of edges.
Definition: edgeList.H:63
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
messageStream Warning
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
label newPointi
Definition: readKivaGrid.H:496
labelList f(nPoints)
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Definition: stdFoam.H:346
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Definition: stdFoam.H:278