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 | Copyright (C) 2011-2016 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "boundaryCutter.H"
27 #include "polyMesh.H"
28 #include "polyTopoChange.H"
29 #include "polyAddCell.H"
30 #include "polyAddFace.H"
31 #include "polyAddPoint.H"
32 #include "polyModifyFace.H"
33 #include "polyModifyPoint.H"
34 #include "mapPolyMesh.H"
35 #include "meshTools.H"
36 
37 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41 defineTypeNameAndDebug(boundaryCutter, 0);
42 }
43 
44 
45 // * * * * * * * * * * * * * Private Static Functions * * * * * * * * * * * //
46 
47 
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 
50 void 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.
80 Foam::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  Map<labelList>::const_iterator fnd = edgeToAddedPoints.find(edgeI);
104 
105  if (fnd != edgeToAddedPoints.end())
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 
143 void 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
202 bool 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  Map<labelList>::const_iterator fnd = edgeToAddedPoints.find(edgeI);
233 
234  if (fnd != edgeToAddedPoints.end())
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  // Dont 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 
388 // Construct from components
389 Foam::boundaryCutter::boundaryCutter(const polyMesh& mesh)
390 :
391  mesh_(mesh),
392  edgeAddedPoints_(),
393  faceAddedPoint_()
394 {}
395 
396 
397 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
398 
400 {}
401 
402 
403 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
404 
406 (
407  const Map<point>& pointToPos,
408  const Map<List<point>>& edgeToCuts,
409  const Map<labelPair>& faceToSplit,
410  const Map<point>& faceToFeaturePoint,
411  polyTopoChange& meshMod
412 )
413 {
414  // Clear and size maps here since mesh size will change.
415  edgeAddedPoints_.clear();
416 
417  faceAddedPoint_.clear();
418  faceAddedPoint_.resize(faceToFeaturePoint.size());
419 
420 
421  //
422  // Points that just need to be moved
423  // Note: could just as well be handled outside of setRefinement.
424  //
425 
426  forAllConstIter(Map<point>, pointToPos, iter)
427  {
428  meshMod.setAction
429  (
431  (
432  iter.key(), // point
433  iter(), // position
434  false, // no zone
435  -1, // zone for point
436  true // supports a cell
437  )
438  );
439  }
440 
441 
442  //
443  // Add new points along cut edges.
444  //
445 
446  // Map from edge label to sorted list of points
447  Map<labelList> edgeToAddedPoints(edgeToCuts.size());
448 
449  forAllConstIter(Map<List<point>>, edgeToCuts, iter)
450  {
451  label edgeI = iter.key();
452 
453  const edge& e = mesh_.edges()[edgeI];
454 
455  // Sorted (from start to end) list of cuts on edge
456  const List<point>& cuts = iter();
457 
458  forAll(cuts, cutI)
459  {
460  // point on feature to move to
461  const point& featurePoint = cuts[cutI];
462 
463  label addedPointi =
464  meshMod.setAction
465  (
467  (
468  featurePoint, // point
469  e.start(), // master point
470  -1, // zone for point
471  true // supports a cell
472  )
473  );
474 
475  Map<labelList>::iterator fnd = edgeToAddedPoints.find(edgeI);
476 
477  if (fnd != edgeToAddedPoints.end())
478  {
479  labelList& addedPoints = fnd();
480 
481  label sz = addedPoints.size();
482  addedPoints.setSize(sz+1);
483  addedPoints[sz] = addedPointi;
484  }
485  else
486  {
487  edgeToAddedPoints.insert(edgeI, labelList(1, addedPointi));
488  }
489 
490  if (debug)
491  {
492  Pout<< "Added point " << addedPointi << " for edge " << edgeI
493  << " with cuts:" << edgeToAddedPoints[edgeI] << endl;
494  }
495  }
496  }
497 
498 
499  //
500  // Introduce feature points.
501  //
502 
503  forAllConstIter(Map<point>, faceToFeaturePoint, iter)
504  {
505  label facei = iter.key();
506 
507  const face& f = mesh_.faces()[facei];
508 
509  if (faceToSplit.found(facei))
510  {
512  << "Face " << facei << " vertices " << f
513  << " is both marked for face-centre decomposition and"
514  << " diagonal splitting."
515  << abort(FatalError);
516  }
517 
518  if (mesh_.isInternalFace(facei))
519  {
521  << "Face " << facei << " vertices " << f
522  << " is not an external face. Cannot split it"
523  << abort(FatalError);
524  }
525 
526  label addedPointi =
527  meshMod.setAction
528  (
530  (
531  iter(), // point
532  f[0], // master point
533  -1, // zone for point
534  true // supports a cell
535  )
536  );
537  faceAddedPoint_.insert(facei, addedPointi);
538 
539  if (debug)
540  {
541  Pout<< "Added point " << addedPointi << " for feature point "
542  << iter() << " on face " << facei << " with centre "
543  << mesh_.faceCentres()[facei] << endl;
544  }
545  }
546 
547 
548  //
549  // Split or retriangulate faces
550  //
551 
552 
553  // Maintain whether face has been updated (for -split edges
554  // -new owner/neighbour)
555  boolList faceUptodate(mesh_.nFaces(), false);
556 
557 
558  // Triangulate faces containing feature points
559  forAllConstIter(Map<label>, faceAddedPoint_, iter)
560  {
561  label facei = iter.key();
562 
563  // Get face with new points on cut edges.
564  face newFace(addEdgeCutsToFace(facei, edgeToAddedPoints));
565 
566  label addedPointi = iter();
567 
568  // Information about old face
569  label patchID, zoneID, zoneFlip;
570  getFaceInfo(facei, patchID, zoneID, zoneFlip);
571  label own = mesh_.faceOwner()[facei];
572  label masterPoint = mesh_.faces()[facei][0];
573 
574  // Triangulate face around mid point
575 
576  face tri(3);
577 
578  forAll(newFace, fp)
579  {
580  label nextV = newFace.nextLabel(fp);
581 
582  tri[0] = newFace[fp];
583  tri[1] = nextV;
584  tri[2] = addedPointi;
585 
586  if (fp == 0)
587  {
588  // Modify the existing face.
589  meshMod.setAction
590  (
592  (
593  tri, // face
594  facei,
595  own, // owner
596  -1, // neighbour
597  false, // flux flip
598  patchID, // patch for face
599  false, // remove from zone
600  zoneID, // zone for face
601  zoneFlip // face zone flip
602  )
603  );
604  }
605  else
606  {
607  // Add additional faces
608  meshMod.setAction
609  (
611  (
612  tri, // face
613  own, // owner
614  -1, // neighbour
615  masterPoint, // master point
616  -1, // master edge
617  -1, // master face for addition
618  false, // flux flip
619  patchID, // patch for face
620  zoneID, // zone for face
621  zoneFlip // face zone flip
622  )
623  );
624  }
625  }
626 
627  faceUptodate[facei] = true;
628  }
629 
630 
631  // Diagonally split faces
632  forAllConstIter(Map<labelPair>, faceToSplit, iter)
633  {
634  label facei = iter.key();
635 
636  const face& f = mesh_.faces()[facei];
637 
638  if (faceAddedPoint_.found(facei))
639  {
641  << "Face " << facei << " vertices " << f
642  << " is both marked for face-centre decomposition and"
643  << " diagonal splitting."
644  << abort(FatalError);
645  }
646 
647 
648  // Get face with new points on cut edges.
649  face newFace(addEdgeCutsToFace(facei, edgeToAddedPoints));
650 
651  // Information about old face
652  label patchID, zoneID, zoneFlip;
653  getFaceInfo(facei, patchID, zoneID, zoneFlip);
654  label own = mesh_.faceOwner()[facei];
655  label masterPoint = mesh_.faces()[facei][0];
656 
657  // Split face from one side of diagonal to other.
658  const labelPair& diag = iter();
659 
660  label fp0 = findIndex(newFace, f[diag[0]]);
661  label fp1 = findIndex(newFace, f[diag[1]]);
662 
663  if (fp0 == -1 || fp1 == -1 || fp0 == fp1)
664  {
666  << "Problem : Face " << facei << " vertices " << f
667  << " newFace:" << newFace << " diagonal:" << f[diag[0]]
668  << ' ' << f[diag[1]]
669  << abort(FatalError);
670  }
671 
672  // Replace existing face by newFace from fp0 to fp1 and add new one
673  // from fp1 to fp0.
674 
675  DynamicList<label> newVerts(newFace.size());
676 
677  // Get vertices from fp0 to (and including) fp1
678  label fp = fp0;
679 
680  do
681  {
682  newVerts.append(newFace[fp]);
683 
684  fp = (fp == newFace.size()-1 ? 0 : fp+1);
685 
686  } while (fp != fp1);
687 
688  newVerts.append(newFace[fp1]);
689 
690 
691  // Modify the existing face.
692  meshMod.setAction
693  (
695  (
696  face(newVerts.shrink()), // face
697  facei,
698  own, // owner
699  -1, // neighbour
700  false, // flux flip
701  patchID, // patch for face
702  false, // remove from zone
703  zoneID, // zone for face
704  zoneFlip // face zone flip
705  )
706  );
707 
708 
709  newVerts.clear();
710 
711  // Get vertices from fp1 to (and including) fp0
712 
713  do
714  {
715  newVerts.append(newFace[fp]);
716 
717  fp = (fp == newFace.size()-1 ? 0 : fp+1);
718 
719  } while (fp != fp0);
720 
721  newVerts.append(newFace[fp0]);
722 
723  // Add additional face
724  meshMod.setAction
725  (
727  (
728  face(newVerts.shrink()), // face
729  own, // owner
730  -1, // neighbour
731  masterPoint, // master point
732  -1, // master edge
733  -1, // master face for addition
734  false, // flux flip
735  patchID, // patch for face
736  zoneID, // zone for face
737  zoneFlip // face zone flip
738  )
739  );
740 
741  faceUptodate[facei] = true;
742  }
743 
744 
745  // Split external faces without feature point but using cut edges.
746  // Does right handed walk but not really.
747  forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
748  {
749  label edgeI = iter.key();
750 
751  const labelList& eFaces = mesh_.edgeFaces()[edgeI];
752 
753  forAll(eFaces, i)
754  {
755  label facei = eFaces[i];
756 
757  if (!faceUptodate[facei] && !mesh_.isInternalFace(facei))
758  {
759  // Is external face so split
760  if (splitFace(facei, pointToPos, edgeToAddedPoints, meshMod))
761  {
762  // Successfull split
763  faceUptodate[facei] = true;
764  }
765  }
766  }
767  }
768 
769 
770  // Add cut edges (but don't split) any other faces using any cut edge.
771  // These can be external faces where splitFace hasn't cut them or
772  // internal faces.
773  forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
774  {
775  label edgeI = iter.key();
776 
777  const labelList& eFaces = mesh_.edgeFaces()[edgeI];
778 
779  forAll(eFaces, i)
780  {
781  label facei = eFaces[i];
782 
783  if (!faceUptodate[facei])
784  {
785  // Renumber face to include split edges.
786  face newFace(addEdgeCutsToFace(facei, edgeToAddedPoints));
787 
788  label own = mesh_.faceOwner()[facei];
789 
790  label nei = -1;
791 
792  if (mesh_.isInternalFace(facei))
793  {
794  nei = mesh_.faceNeighbour()[facei];
795  }
796 
797  label patchID, zoneID, zoneFlip;
798  getFaceInfo(facei, patchID, zoneID, zoneFlip);
799 
800  meshMod.setAction
801  (
803  (
804  newFace, // modified face
805  facei, // label of face being modified
806  own, // owner
807  nei, // neighbour
808  false, // face flip
809  patchID, // patch for face
810  false, // remove from zone
811  zoneID, // zone for face
812  zoneFlip // face flip in zone
813  )
814  );
815 
816  faceUptodate[facei] = true;
817  }
818  }
819  }
820 
821  // Convert edge to points storage from edge labels (not preserved)
822  // to point labels
823  edgeAddedPoints_.resize(edgeToCuts.size());
824 
825  forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
826  {
827  edgeAddedPoints_.insert(mesh_.edges()[iter.key()], iter());
828  }
829 }
830 
831 
833 {
834  // Update stored labels for mesh change.
835 
836  //
837  // Do faceToAddedPoint
838  //
839 
840  {
841  // Create copy since we're deleting entries.
842  Map<label> newAddedPoints(faceAddedPoint_.size());
843 
844  forAllConstIter(Map<label>, faceAddedPoint_, iter)
845  {
846  label oldFacei = iter.key();
847 
848  label newFacei = morphMap.reverseFaceMap()[oldFacei];
849 
850  label oldPointi = iter();
851 
852  label newPointi = morphMap.reversePointMap()[oldPointi];
853 
854  if (newFacei >= 0 && newPointi >= 0)
855  {
856  newAddedPoints.insert(newFacei, newPointi);
857  }
858  }
859 
860  // Copy
861  faceAddedPoint_.transfer(newAddedPoints);
862  }
863 
864 
865  //
866  // Do edgeToAddedPoints
867  //
868 
869 
870  {
871  // Create copy since we're deleting entries
873  newEdgeAddedPoints(edgeAddedPoints_.size());
874 
875  for
876  (
877  HashTable<labelList, edge, Hash<edge>>::const_iterator iter =
878  edgeAddedPoints_.begin();
879  iter != edgeAddedPoints_.end();
880  ++iter
881  )
882  {
883  const edge& e = iter.key();
884 
885  label newStart = morphMap.reversePointMap()[e.start()];
886 
887  label newEnd = morphMap.reversePointMap()[e.end()];
888 
889  if (newStart >= 0 && newEnd >= 0)
890  {
891  const labelList& addedPoints = iter();
892 
893  labelList newAddedPoints(addedPoints.size());
894  label newI = 0;
895 
896  forAll(addedPoints, i)
897  {
898  label newAddedPointi =
899  morphMap.reversePointMap()[addedPoints[i]];
900 
901  if (newAddedPointi >= 0)
902  {
903  newAddedPoints[newI++] = newAddedPointi;
904  }
905  }
906  if (newI > 0)
907  {
908  newAddedPoints.setSize(newI);
909 
910  edge newE = edge(newStart, newEnd);
911 
912  newEdgeAddedPoints.insert(newE, newAddedPoints);
913  }
914  }
915  }
916 
917  // Copy
918  edgeAddedPoints_.transfer(newEdgeAddedPoints);
919  }
920 }
921 
922 
923 // ************************************************************************* //
label end() const
Return end vertex label.
Definition: edgeI.H:92
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
Class describing modification of a face.
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
const double e
Elementary charge.
Definition: doubleFloat.H:78
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Definition: UList.H:440
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:337
HashTable< T, label, Hash< label > >::const_iterator const_iterator
Definition: Map.H:59
Class containing data for point addition.
Definition: polyAddPoint.H:47
label size() const
Return number of elements in table.
void transfer(HashTable< T, Key, Hash > &)
Transfer the contents of the argument table into this table.
Definition: HashTable.C:509
const vectorField & faceCentres() const
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
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 labelListList & edgeFaces() const
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.
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:138
An STL-conforming iterator.
Definition: HashTable.H:415
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
dynamicFvMesh & mesh
const labelList & reverseFaceMap() const
Reverse face map.
Definition: mapPolyMesh.H:495
List< edge > edgeList
Definition: edgeList.H:38
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
void clear()
Clear all entries from table.
Definition: HashTable.C:464
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:292
List< label > labelList
A List of labels.
Definition: labelList.H:56
iterator begin()
Iterator set to the beginning of the HashTable.
Definition: HashTableI.H:419
An STL-conforming hash table.
Definition: HashTable.H:61
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
bool found(const label &) const
Return true if hashedEntry is found in table.
Class describing modification of a point.
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
label start() const
Return start vertex label.
Definition: edgeI.H:81
static const char nl
Definition: Ostream.H:262
defineTypeNameAndDebug(combustionModel, 0)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
messageStream Warning
label nFaces() const
void setSize(const label)
Reset size of List.
Definition: List.C:295
~boundaryCutter()
Destructor.
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1023
label newPointi
Definition: readKivaGrid.H:501
Direct mesh changes based on v1.3 polyTopoChange syntax.
Hash function class for primitives. All non-primitives used to hash entries on hash tables likely nee...
Definition: Hash.H:54
void resize(const label newSize)
Resize the hash table for efficiency.
Definition: HashTable.C:428
label nextLabel(const label i) const
Next vertex on face.
Definition: faceI.H:117
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1017
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1004
Namespace for OpenFOAM.
void updateMesh(const mapPolyMesh &)
Force recalculation of locally stored data on topological change.
const labelList & reversePointMap() const
Reverse point map.
Definition: mapPolyMesh.H:463
A HashTable to objects of type <T> with a label key.
Definition: Map.H:49
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.