combineFaces.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) 2018-2022 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 "combineFaces.H"
30#include "polyMesh.H"
31#include "polyTopoChange.H"
32#include "polyRemoveFace.H"
33#include "polyAddFace.H"
34#include "polyModifyFace.H"
35#include "polyRemovePoint.H"
36#include "polyAddPoint.H"
37#include "syncTools.H"
38#include "meshTools.H"
39
40// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41
42namespace Foam
43{
45}
46
47
48// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49
50// Test face for (almost) convexeness. Allows certain convexity before
51// complaining.
52// For every two consecutive edges calculate the normal. If it differs too
53// much from the average face normal complain.
54bool Foam::combineFaces::convexFace
55(
56 const scalar minConcaveCos,
57 const pointField& points,
58 const face& f
59)
60{
61 // Get outwards pointing normal of f, only the sign matters.
62 const vector areaNorm = f.areaNormal(points);
63
64 // Normalized vector from f[size-1] to f[0];
65 vector ePrev(points[f.first()] - points[f.last()]);
66 scalar magEPrev = mag(ePrev);
67 ePrev /= magEPrev + VSMALL;
68
69 forAll(f, fp0)
70 {
71 // Normalized vector between two consecutive points
72 vector e10(points[f.nextLabel(fp0)] - points[f.thisLabel(fp0)]);
73 scalar magE10 = mag(e10);
74 e10 /= magE10 + VSMALL;
75
76 if (magEPrev > SMALL && magE10 > SMALL)
77 {
78 vector edgeNormal = ePrev ^ e10;
79
80 if ((edgeNormal & areaNorm) < 0)
81 {
82 // Concave. Check angle.
83 if ((ePrev & e10) < minConcaveCos)
84 {
85 return false;
86 }
87 }
88 }
89
90 ePrev = e10;
91 magEPrev = magE10;
92 }
93
94 // Not a single internal angle is concave so face is convex.
95 return true;
96}
97
98
99// Determines if set of faces is valid to collapse into single face.
100bool Foam::combineFaces::validFace
101(
102 const scalar minConcaveCos,
103 const indirectPrimitivePatch& bigFace
104)
105{
106 // Get outside vertices (in local vertex numbering)
107 const labelListList& edgeLoops = bigFace.edgeLoops();
108
109 if (edgeLoops.size() > 1)
110 {
111 return false;
112 }
113
114 bool isNonManifold = bigFace.checkPointManifold(false, nullptr);
115 if (isNonManifold)
116 {
117 return false;
118 }
119
120 // Check for convexness
121 face f(getOutsideFace(bigFace));
122
123 return convexFace(minConcaveCos, bigFace.points(), f);
124}
125
126
127void Foam::combineFaces::regioniseFaces
128(
129 const scalar minCos,
130 const bool mergeAcrossPatches,
131 const label celli,
132 const labelList& cEdges,
133 Map<label>& faceRegion
134) const
135{
136 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
137
138 for (const label edgei : cEdges)
139 {
140 label f0, f1;
141 meshTools::getEdgeFaces(mesh_, celli, edgei, f0, f1);
142
143 const vector& a0 = mesh_.faceAreas()[f0];
144 const vector& a1 = mesh_.faceAreas()[f1];
145
146 const label p0 = patches.whichPatch(f0);
147 const label p1 = patches.whichPatch(f1);
148
149 // Face can be merged if
150 // - small angle
151 // - mergeAcrossPatches=false : same non-constraint patch
152 // - mergeAcrossPatches=true : always (if non-constraint patch)
153 // (this logic could be extended to e.g. merge faces on symm plane
154 // if they have similar normals. But there might be lots of other
155 // constraints which disallow merging so this decision ideally should
156 // be up to patch type)
157 if
158 (
159 p0 != -1
160 && p1 != -1
161 && !(
164 )
165 )
166 {
167 if (!mergeAcrossPatches && (p0 != p1))
168 {
169 continue;
170 }
171
172 const vector f0Normal = normalised(a0);
173 const vector f1Normal = normalised(a1);
174
175 if ((f0Normal & f1Normal) > minCos)
176 {
177 const label region0 = faceRegion.lookup(f0, -1);
178 const label region1 = faceRegion.lookup(f1, -1);
179
180 if (region0 == -1)
181 {
182 if (region1 == -1)
183 {
184 const label useRegion = faceRegion.size();
185 faceRegion.insert(f0, useRegion);
186 faceRegion.insert(f1, useRegion);
187 }
188 else
189 {
190 faceRegion.insert(f0, region1);
191 }
192 }
193 else
194 {
195 if (region1 == -1)
196 {
197 faceRegion.insert(f1, region0);
198 }
199 else if (region0 != region1)
200 {
201 // Merge the two regions
202 const label useRegion = min(region0, region1);
203 const label freeRegion = max(region0, region1);
204
205 forAllIters(faceRegion, iter)
206 {
207 if (iter.val() == freeRegion)
208 {
209 iter.val() = useRegion;
210 }
211 }
212 }
213 }
214 }
215 }
216 }
217}
218
219
220bool Foam::combineFaces::faceNeighboursValid
221(
222 const label celli,
223 const Map<label>& faceRegion
224) const
225{
226 if (faceRegion.size() <= 1)
227 {
228 return true;
229 }
230
231 const cell& cFaces = mesh_.cells()[celli];
232
233 DynamicList<label> storage;
234
235 // Test for face collapsing to edge since too many neighbours merged.
236 forAll(cFaces, cFacei)
237 {
238 label facei = cFaces[cFacei];
239
240 if (!faceRegion.found(facei))
241 {
242 const labelList& fEdges = mesh_.faceEdges(facei, storage);
243
244 // Count number of remaining faces neighbouring facei. This has
245 // to be 3 or more.
246
247 // Unregioned neighbouring faces
248 DynamicList<label> neighbourFaces(cFaces.size());
249 // Regioned neighbouring faces
250 labelHashSet neighbourRegions;
251
252 forAll(fEdges, i)
253 {
254 const label edgeI = fEdges[i];
255 label nbrI = meshTools::otherFace(mesh_, celli, facei, edgeI);
256
257 const auto iter = faceRegion.cfind(nbrI);
258
259 if (iter.found())
260 {
261 neighbourRegions.insert(iter.val());
262 }
263 else
264 {
265 neighbourFaces.appendUniq(nbrI);
266 }
267 }
268
269 if ((neighbourFaces.size()+neighbourRegions.size()) < 3)
270 {
271 return false;
272 }
273 }
274 }
275 return true;
276}
277
278
279// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
280
281// Construct from mesh
283(
284 const polyMesh& mesh,
285 const bool undoable
286)
287:
288 mesh_(mesh),
289 undoable_(undoable),
290 masterFace_(0),
291 faceSetsVertices_(0),
292 savedPointLabels_(0),
293 savedPoints_(0)
294{}
295
296
297// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
298
300(
301 const scalar featureCos,
302 const scalar minConcaveCos,
303 const labelHashSet& boundaryCells,
304 const bool mergeAcrossPatches
305) const
306{
307 // Lists of faces that can be merged.
308 DynamicList<labelList> allFaceSets(boundaryCells.size() / 10);
309 // Storage for on-the-fly cell-edge addressing.
310 labelHashSet set;
311 DynamicList<label> storage;
312
313 // On all cells regionise the faces
314 for (const label celli : boundaryCells)
315 {
316 const cell& cFaces = mesh_.cells()[celli];
317
318 const labelList& cEdges = mesh_.cellEdges(celli, set, storage);
319
320 // Region per face
321 Map<label> faceRegion(cFaces.size());
322 regioniseFaces
323 (
324 featureCos,
325 mergeAcrossPatches,
326 celli,
327 cEdges,
328 faceRegion
329 );
330
331 // Now we have in faceRegion for every face the region with planar
332 // face sharing the same region. We now check whether the resulting
333 // sets cause a face
334 // - to become a set of edges since too many faces are merged.
335 // - to become convex
336
337 if (faceNeighboursValid(celli, faceRegion))
338 {
339 // Create region-to-faces addressing
340 Map<labelList> regionToFaces(faceRegion.size());
341
342 forAllConstIters(faceRegion, iter)
343 {
344 const label facei = iter.key();
345 const label region = iter.val();
346
347 auto regionFnd = regionToFaces.find(region);
348
349 if (regionFnd.found())
350 {
351 labelList& setFaces = regionFnd();
352 label sz = setFaces.size();
353 setFaces.setSize(sz+1);
354 setFaces[sz] = facei;
355 }
356 else
357 {
358 regionToFaces.insert(region, labelList(1, facei));
359 }
360 }
361
362 // For every set check if it forms a valid convex face
363
364 forAllIters(regionToFaces, iter)
365 {
366 // Make face out of setFaces
368 (
370 (
371 mesh_.faces(),
372 iter.val()
373 ),
374 mesh_.points()
375 );
376
377 // Only store if -only one outside loop -which forms convex face
378 if (validFace(minConcaveCos, bigFace))
379 {
380 labelList& faceIDs = iter.val();
381
382 // For cross-patch merging we want to make the
383 // largest face the one to decide the final patch
384 // (i.e. master face)
385 if (mergeAcrossPatches)
386 {
387 const vectorField& areas = mesh_.faceAreas();
388
389 label maxIndex = 0;
390 scalar maxMagSqr = magSqr(areas[faceIDs[0]]);
391 for (label i = 1; i < faceIDs.size(); ++i)
392 {
393 const scalar a2 = magSqr(areas[faceIDs[i]]);
394 if (a2 > maxMagSqr)
395 {
396 maxMagSqr = a2;
397 maxIndex = i;
398 }
399 }
400 if (maxIndex != 0)
401 {
402 std::swap(faceIDs[0], faceIDs[maxIndex]);
403 }
404 }
405
406 allFaceSets.append(faceIDs);
407 }
408 }
409 }
410 }
411
412 return allFaceSets.shrink();
413}
414
415
417(
418 const scalar featureCos,
419 const scalar minConcaveCos,
420 const bool mergeAcrossPatches
421) const
422{
423 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
424
425 // Pick up all cells on boundary
426 labelHashSet boundaryCells(mesh_.nBoundaryFaces());
427
428 forAll(patches, patchi)
429 {
430 const polyPatch& patch = patches[patchi];
431
432 if (!patch.coupled())
433 {
434 forAll(patch, i)
435 {
436 boundaryCells.insert(mesh_.faceOwner()[patch.start()+i]);
437 }
438 }
439 }
440
441 return getMergeSets
442 (
443 featureCos,
444 minConcaveCos,
445 boundaryCells,
446 mergeAcrossPatches
447 );
448}
449
450
451// Gets outside edgeloop as a face
452// - in same order as faces
453// - in mesh vertex labels
455(
456 const indirectPrimitivePatch& fp
457)
458{
459 if (fp.edgeLoops().size() != 1)
460 {
462 << "Multiple outside loops:" << fp.edgeLoops()
463 << abort(FatalError);
464 }
465
466 // Get first boundary edge. Since guaranteed one edgeLoop when in here this
467 // edge must be on it.
468 label bEdgeI = fp.nInternalEdges();
469
470 const edge& e = fp.edges()[bEdgeI];
471
472 const labelList& eFaces = fp.edgeFaces()[bEdgeI];
473
474 if (eFaces.size() != 1)
475 {
477 << "boundary edge:" << bEdgeI
478 << " points:" << fp.meshPoints()[e[0]]
479 << ' ' << fp.meshPoints()[e[1]]
480 << " on indirectPrimitivePatch has " << eFaces.size()
481 << " faces using it" << abort(FatalError);
482 }
483
484
485 // Outside loop
486 const labelList& outsideLoop = fp.edgeLoops()[0];
487
488
489 // Get order of edge e in outside loop
490 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
491
492 // edgeLoopConsistent : true = edge in same order as outsideloop
493 // false = opposite order
494 bool edgeLoopConsistent = false;
495
496 {
497 label index0 = outsideLoop.find(e[0]);
498 label index1 = outsideLoop.find(e[1]);
499
500 if (index0 == -1 || index1 == -1)
501 {
503 << "Cannot find boundary edge:" << e
504 << " points:" << fp.meshPoints()[e[0]]
505 << ' ' << fp.meshPoints()[e[1]]
506 << " in edgeLoop:" << outsideLoop << abort(FatalError);
507 }
508 else if (index1 == outsideLoop.fcIndex(index0))
509 {
510 edgeLoopConsistent = true;
511 }
512 else if (index0 == outsideLoop.fcIndex(index1))
513 {
514 edgeLoopConsistent = false;
515 }
516 else
517 {
519 << "Cannot find boundary edge:" << e
520 << " points:" << fp.meshPoints()[e[0]]
521 << ' ' << fp.meshPoints()[e[1]]
522 << " on consecutive points in edgeLoop:"
523 << outsideLoop << abort(FatalError);
524 }
525 }
526
527
528 // Get face in local vertices
529 const face& localF = fp.localFaces()[eFaces[0]];
530
531 // Get order of edge in localF
532 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
533
534 // faceEdgeConsistent : true = edge in same order as localF
535 // false = opposite order
536 bool faceEdgeConsistent = false;
537
538 {
539 // Find edge in face.
540 label index = fp.faceEdges()[eFaces[0]].find(bEdgeI);
541
542 if (index == -1)
543 {
545 << "Cannot find boundary edge:" << e
546 << " points:" << fp.meshPoints()[e[0]]
547 << ' ' << fp.meshPoints()[e[1]]
548 << " in face:" << eFaces[0]
549 << " edges:" << fp.faceEdges()[eFaces[0]]
550 << abort(FatalError);
551 }
552 else if (localF[index] == e[0] && localF.nextLabel(index) == e[1])
553 {
554 faceEdgeConsistent = true;
555 }
556 else if (localF[index] == e[1] && localF.nextLabel(index) == e[0])
557 {
558 faceEdgeConsistent = false;
559 }
560 else
561 {
563 << "Cannot find boundary edge:" << e
564 << " points:" << fp.meshPoints()[e[0]]
565 << ' ' << fp.meshPoints()[e[1]]
566 << " in face:" << eFaces[0] << " verts:" << localF
567 << abort(FatalError);
568 }
569 }
570
571 // Get face in mesh points.
572 face meshFace(renumber(fp.meshPoints(), outsideLoop));
573
574 if (faceEdgeConsistent != edgeLoopConsistent)
575 {
576 reverse(meshFace);
577 }
578 return meshFace;
579}
580
581
583(
584 const labelListList& faceSets,
585 polyTopoChange& meshMod
586)
587{
588 if (undoable_)
589 {
590 masterFace_.setSize(faceSets.size());
591 faceSetsVertices_.setSize(faceSets.size());
592 forAll(faceSets, setI)
593 {
594 const labelList& setFaces = faceSets[setI];
595
596 masterFace_[setI] = setFaces[0];
597 faceSetsVertices_[setI] = UIndirectList<face>
598 (
599 mesh_.faces(),
600 setFaces
601 );
602 }
603 }
604
605 // Running count of number of faces using a point
606 labelList nPointFaces(mesh_.nPoints(), Zero);
607
608 const labelListList& pointFaces = mesh_.pointFaces();
609
610 forAll(pointFaces, pointi)
611 {
612 nPointFaces[pointi] = pointFaces[pointi].size();
613 }
614
615 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
616
617 forAll(faceSets, setI)
618 {
619 const labelList& setFaces = faceSets[setI];
620
621 // Check
622 if (debug)
623 {
624 forAll(setFaces, i)
625 {
626 label patchi = patches.whichPatch(setFaces[i]);
627
628 if (patchi == -1 || patches[patchi].coupled())
629 {
631 << "Can only merge non-coupled boundary faces"
632 << " but found internal or coupled face:"
633 << setFaces[i] << " in set " << setI
634 << abort(FatalError);
635 }
636 }
637 }
638
639 // Make face out of setFaces
641 (
643 (
644 mesh_.faces(),
645 setFaces
646 ),
647 mesh_.points()
648 );
649
650 // Get outside vertices (in local vertex numbering)
651 const labelListList& edgeLoops = bigFace.edgeLoops();
652
653 if (edgeLoops.size() != 1)
654 {
656 << "Faces to-be-merged " << setFaces
657 << " do not form a single big face." << nl
658 << abort(FatalError);
659 }
660
661
662 // Delete all faces except master, whose face gets modified.
663
664 // Modify master face
665 // ~~~~~~~~~~~~~~~~~~
666
667 label masterFacei = setFaces[0];
668
669 // Get outside face in mesh vertex labels
670 face outsideFace(getOutsideFace(bigFace));
671
672 label zoneID = mesh_.faceZones().whichZone(masterFacei);
673
674 bool zoneFlip = false;
675
676 if (zoneID >= 0)
677 {
678 const faceZone& fZone = mesh_.faceZones()[zoneID];
679
680 zoneFlip = fZone.flipMap()[fZone.whichFace(masterFacei)];
681 }
682
683 label patchi = mesh_.boundaryMesh().whichPatch(masterFacei);
684
685 meshMod.setAction
686 (
688 (
689 outsideFace, // modified face
690 masterFacei, // label of face being modified
691 mesh_.faceOwner()[masterFacei], // owner
692 -1, // neighbour
693 false, // face flip
694 patchi, // patch for face
695 false, // remove from zone
696 zoneID, // zone for face
697 zoneFlip // face flip in zone
698 )
699 );
700
701
702 // Delete all non-master faces
703 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
704
705 for (label i = 1; i < setFaces.size(); i++)
706 {
707 meshMod.setAction(polyRemoveFace(setFaces[i]));
708 }
709
710
711 // Mark unused points for deletion
712 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
713
714 // Remove count of all faces combined
715 forAll(setFaces, i)
716 {
717 const face& f = mesh_.faces()[setFaces[i]];
718
719 forAll(f, fp)
720 {
721 nPointFaces[f[fp]]--;
722 }
723 }
724 // But keep count on new outside face
725 forAll(outsideFace, fp)
726 {
727 nPointFaces[outsideFace[fp]]++;
728 }
729 }
730
731
732 // If one or more processors use the point it needs to be kept.
734 (
735 mesh_,
736 nPointFaces,
738 label(0) // null value
739 );
740
741 // Remove all unused points. Store position if undoable.
742 if (!undoable_)
743 {
744 forAll(nPointFaces, pointi)
745 {
746 if (nPointFaces[pointi] == 0)
747 {
748 meshMod.setAction(polyRemovePoint(pointi));
749 }
750 }
751 }
752 else
753 {
754 // Count removed points
755 label n = 0;
756 forAll(nPointFaces, pointi)
757 {
758 if (nPointFaces[pointi] == 0)
759 {
760 n++;
761 }
762 }
763
764 savedPointLabels_.setSize(n);
765 savedPoints_.setSize(n);
766 Map<label> meshToSaved(2*n);
767
768 // Remove points and store position
769 n = 0;
770 forAll(nPointFaces, pointi)
771 {
772 if (nPointFaces[pointi] == 0)
773 {
774 meshMod.setAction(polyRemovePoint(pointi));
775
776 savedPointLabels_[n] = pointi;
777 savedPoints_[n] = mesh_.points()[pointi];
778
779 meshToSaved.insert(pointi, n);
780 n++;
781 }
782 }
783
784 // Update stored vertex labels. Negative indices index into local points
785 forAll(faceSetsVertices_, setI)
786 {
787 faceList& setFaces = faceSetsVertices_[setI];
788
789 forAll(setFaces, i)
790 {
791 face& f = setFaces[i];
792
793 forAll(f, fp)
794 {
795 label pointi = f[fp];
796
797 if (nPointFaces[pointi] == 0)
798 {
799 f[fp] = -meshToSaved[pointi]-1;
800 }
801 }
802 }
803 }
804 }
805}
806
807
809{
810 if (undoable_)
811 {
812 // Master face just renumbering of point labels
813 inplaceRenumber(map.reverseFaceMap(), masterFace_);
814
815 // Stored faces refer to backed-up vertices (not changed)
816 // and normal vertices (need to be renumbered)
817 forAll(faceSetsVertices_, setI)
818 {
819 faceList& faces = faceSetsVertices_[setI];
820
821 forAll(faces, i)
822 {
823 // Note: faces[i] can have negative elements.
824 face& f = faces[i];
825
826 forAll(f, fp)
827 {
828 label pointi = f[fp];
829
830 if (pointi >= 0)
831 {
832 f[fp] = map.reversePointMap()[pointi];
833
834 if (f[fp] < 0)
835 {
837 << "In set " << setI << " at position " << i
838 << " with master face "
839 << masterFace_[setI] << nl
840 << "the points of the slave face " << faces[i]
841 << " don't exist anymore."
842 << abort(FatalError);
843 }
844 }
845 }
846 }
847 }
848 }
849}
850
851
852// Note that faces on coupled patches are never combined (or at least
853// getMergeSets never picks them up. Thus the points on them will never get
854// deleted since still used by the coupled faces. So never a risk of undoing
855// things (faces/points) on coupled patches.
857(
858 const labelList& masterFaces,
859 polyTopoChange& meshMod,
860 Map<label>& restoredPoints,
861 Map<label>& restoredFaces,
862 Map<label>& restoredCells
863)
864{
865 if (!undoable_)
866 {
868 << "Can only call setUnrefinement if constructed with"
869 << " unrefinement capability." << exit(FatalError);
870 }
871
872
873 // Restored points
874 labelList addedPoints(savedPoints_.size(), -1);
875
876 // Invert set-to-master-face
877 Map<label> masterToSet(masterFace_.size());
878
879 forAll(masterFace_, setI)
880 {
881 if (masterFace_[setI] >= 0)
882 {
883 masterToSet.insert(masterFace_[setI], setI);
884 }
885 }
886
887 forAll(masterFaces, i)
888 {
889 const label masterFacei = masterFaces[i];
890
891 const auto iter = masterToSet.cfind(masterFacei);
892
893 if (!iter.found())
894 {
896 << "Master face " << masterFacei
897 << " is not the master of one of the merge sets"
898 << " or has already been merged"
899 << abort(FatalError);
900 }
901
902 const label setI = iter.val();
903
904
905 // Update faces of the merge setI for reintroduced vertices
906 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
907
908 faceList& faces = faceSetsVertices_[setI];
909
910 if (faces.empty())
911 {
913 << "Set " << setI << " with master face " << masterFacei
914 << " has already been merged." << abort(FatalError);
915 }
916
917 forAll(faces, i)
918 {
919 face& f = faces[i];
920
921 forAll(f, fp)
922 {
923 label pointi = f[fp];
924
925 if (pointi < 0)
926 {
927 label localI = -pointi-1;
928
929 if (addedPoints[localI] == -1)
930 {
931 // First occurrence of saved point. Reintroduce point
932 addedPoints[localI] = meshMod.setAction
933 (
935 (
936 savedPoints_[localI], // point
937 -1, // master point
938 -1, // zone for point
939 true // supports a cell
940 )
941 );
942 restoredPoints.insert
943 (
944 addedPoints[localI], // current point label
945 savedPointLabels_[localI] // point label when it
946 // was stored
947 );
948 }
949 f[fp] = addedPoints[localI];
950 }
951 }
952 }
953
954
955 // Restore
956 // ~~~~~~~
957
958 label own = mesh_.faceOwner()[masterFacei];
959 label zoneID = mesh_.faceZones().whichZone(masterFacei);
960 bool zoneFlip = false;
961 if (zoneID >= 0)
962 {
963 const faceZone& fZone = mesh_.faceZones()[zoneID];
964 zoneFlip = fZone.flipMap()[fZone.whichFace(masterFacei)];
965 }
966 label patchi = mesh_.boundaryMesh().whichPatch(masterFacei);
967
968 if (mesh_.boundaryMesh()[patchi].coupled())
969 {
971 << "Master face " << masterFacei << " is on coupled patch "
972 << mesh_.boundaryMesh()[patchi].name()
973 << abort(FatalError);
974 }
975
976 //Pout<< "Restoring new master face " << masterFacei
977 // << " to vertices " << faces[0] << endl;
978
979 // Modify the master face.
980 meshMod.setAction
981 (
983 (
984 faces[0], // original face
985 masterFacei, // label of face
986 own, // owner
987 -1, // neighbour
988 false, // face flip
989 patchi, // patch for face
990 false, // remove from zone
991 zoneID, // zone for face
992 zoneFlip // face flip in zone
993 )
994 );
995 restoredFaces.insert(masterFacei, masterFacei);
996
997 // Add the previously removed faces
998 for (label i = 1; i < faces.size(); i++)
999 {
1000 //Pout<< "Restoring removed face with vertices " << faces[i]
1001 // << endl;
1002
1003 label facei = meshMod.setAction
1004 (
1006 (
1007 faces[i], // vertices
1008 own, // owner,
1009 -1, // neighbour,
1010 -1, // masterPointID,
1011 -1, // masterEdgeID,
1012 masterFacei, // masterFaceID,
1013 false, // flipFaceFlux,
1014 patchi, // patchID,
1015 zoneID, // zoneID,
1016 zoneFlip // zoneFlip
1017 )
1018 );
1019 restoredFaces.insert(facei, masterFacei);
1020 }
1021
1022 // Clear out restored set
1023 faceSetsVertices_[setI].clear();
1024 masterFace_[setI] = -1;
1025 }
1026}
1027
1028
1029// ************************************************************************* //
label n
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:72
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
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition: HashSet.H:191
const_iterator cfind(const Key &key) const
Find and return an const_iterator set at the hashed entry.
Definition: HashTableI.H:141
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 List with indirect addressing.
Definition: IndirectList.H:119
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
A list of faces which address into the list of points.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
label nInternalEdges() const
Number of internal edges.
const labelList & meshPoints() const
Return labelList of mesh points in patch.
const labelListList & edgeLoops() const
Return list of closed loops of boundary vertices.
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
T & first()
Return the first element of the list.
Definition: UListI.H:202
bool empty() const noexcept
True if the UList is empty (ie, size() is zero)
Definition: UListI.H:427
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
T & last()
Return the last element of the list.
Definition: UListI.H:216
label fcIndex(const label i) const noexcept
Definition: UListI.H:60
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:57
Combines boundary faces into single face. The faces get the patch of the first face ('the master')
Definition: combineFaces.H:59
labelListList getMergeSets(const scalar featureCos, const scalar minConcaveCos, const labelHashSet &boundaryCells, const bool mergeAcrossPatches=false) const
Extract lists of all (non-coupled) boundary faces on selected.
Definition: combineFaces.C:300
void setUnrefinement(const labelList &masterFaces, polyTopoChange &meshMod, Map< label > &restoredPoints, Map< label > &restoredFaces, Map< label > &restoredCells)
Play commands into polyTopoChange to reinsert original faces.
Definition: combineFaces.C:857
void setRefinement(const labelListList &, polyTopoChange &)
Play commands into polyTopoChange to combine faces. Gets.
Definition: combineFaces.C:583
static face getOutsideFace(const indirectPrimitivePatch &)
Gets outside of patch as a face (in mesh point labels)
Definition: combineFaces.C:455
virtual const word & constraintType() const
Return the constraint type this pointPatch implements.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:67
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
Definition: faceZone.C:372
const boolList & flipMap() const noexcept
Return face flip map.
Definition: faceZone.H:272
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
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
Class describing modification of a face.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
Class containing data for face removal.
Class containing data for point removal.
Direct mesh changes based on v1.3 polyTopoChange syntax.
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
static void syncPointList(const polyMesh &mesh, List< T > &pointValues, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh points.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
const volScalarField & p0
Definition: EEqn.H:36
const polyBoundaryMesh & patches
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
const pointField & points
const labelIOList & zoneID
const dimensionedScalar a0
Bohr radius: default SI units: [m].
label otherFace(const primitiveMesh &mesh, const label celli, const label facei, const label edgeI)
Return face on cell using edgeI but not facei. Throws error.
Definition: meshTools.C:555
void getEdgeFaces(const primitiveMesh &mesh, const label celli, const label edgeI, label &face0, label &face1)
Get faces on cell using edgeI. Throws error if no two found.
Definition: meshTools.C:479
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
List< label > labelList
A List of labels.
Definition: List.H:66
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Definition: quaternionI.H:680
void reverse(UList< T > &list, const label n)
Reverse the first n elements of the list.
Definition: UListI.H:449
IntListType renumber(const labelUList &oldToNew, const IntListType &input)
Renumber the values (not the indices) of a list.
Type maxMagSqr(const UList< Type > &f)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:56
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
errorManip< error > abort(error &err)
Definition: errorManip.H:144
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
Definition: HashSet.H:85
error FatalError
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
A PrimitivePatch with an IndirectList for the faces, const reference for the point field.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
labelList f(nPoints)
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllIters(container, iter)
Iterate across all elements in the container object.
Definition: stdFoam.H:260
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Definition: stdFoam.H:278