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