fvMeshSubset.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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 Description
25  Post-processing mesh subset tool. Given the original mesh and the
26  list of selected cells, it creates the mesh consisting only of the
27  desired cells, with the mapping list for points, faces, and cells.
28 
29 \*---------------------------------------------------------------------------*/
30 
31 #include "fvMeshSubset.H"
32 #include "boolList.H"
33 #include "Pstream.H"
34 #include "emptyPolyPatch.H"
35 #include "demandDrivenData.H"
36 #include "cyclicPolyPatch.H"
37 #include "removeCells.H"
38 #include "polyTopoChange.H"
39 #include "mapPolyMesh.H"
40 
41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47 
48 bool Foam::fvMeshSubset::checkCellSubset() const
49 {
50  if (fvMeshSubsetPtr_.empty())
51  {
53  << "void setCellSubset(const labelHashSet& cellsToSubset)" << endl
54  << "before attempting to access subset data"
55  << abort(FatalError);
56 
57  return false;
58  }
59  else
60  {
61  return true;
62  }
63 }
64 
65 
66 void Foam::fvMeshSubset::markPoints
67 (
68  const labelList& curPoints,
69  Map<label>& pointMap
70 )
71 {
72  forAll(curPoints, pointi)
73  {
74  // Note: insert will only insert if not yet there.
75  pointMap.insert(curPoints[pointi], 0);
76  }
77 }
78 
79 
80 void Foam::fvMeshSubset::markPoints
81 (
82  const labelList& curPoints,
83  labelList& pointMap
84 )
85 {
86  forAll(curPoints, pointi)
87  {
88  pointMap[curPoints[pointi]] = 0;
89  }
90 }
91 
92 
93 void Foam::fvMeshSubset::doCoupledPatches
94 (
95  const bool syncPar,
96  Map<label>& facesToSubset,
97  labelList& nCellsUsingFace
98 ) const
99 {
100  // Synchronize facesToSubset on both sides of coupled patches.
101  // Marks faces that become 'uncoupled' with 3.
102 
103  const polyBoundaryMesh& oldPatches = baseMesh().boundaryMesh();
104 
105  label nUncoupled = 0;
106 
107  if (syncPar && Pstream::parRun())
108  {
109  PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
110 
111  // Send face usage across processor patches
112  forAll(oldPatches, oldPatchi)
113  {
114  const polyPatch& pp = oldPatches[oldPatchi];
115 
116  if (isA<processorPolyPatch>(pp))
117  {
118  const processorPolyPatch& procPatch =
119  refCast<const processorPolyPatch>(pp);
120 
121  UOPstream toNeighbour(procPatch.neighbProcNo(), pBufs);
122 
123  if (!facesToSubset.empty())
124  {
125  DynamicList<label> patchFacesToSubset;
126  forAll(pp, i)
127  {
128  if
129  (
130  facesToSubset.found(pp.start()+i)
131  && facesToSubset[pp.start()+i] == 1
132  )
133  {
134  patchFacesToSubset.append(i);
135  }
136  }
137  toNeighbour << patchFacesToSubset;
138  }
139  else if (!nCellsUsingFace.empty())
140  {
141  toNeighbour <<
142  SubList<label>(nCellsUsingFace, pp.size(), pp.start());
143  }
144  else
145  {
146  toNeighbour << labelList();
147  }
148  }
149  }
150 
151  pBufs.finishedSends();
152 
153  // Receive face usage count and check for faces that become uncoupled.
154  forAll(oldPatches, oldPatchi)
155  {
156  const polyPatch& pp = oldPatches[oldPatchi];
157 
158  if (isA<processorPolyPatch>(pp))
159  {
160  const processorPolyPatch& procPatch =
161  refCast<const processorPolyPatch>(pp);
162 
163  UIPstream fromNeighbour(procPatch.neighbProcNo(), pBufs);
164 
165  const labelList nbrList(fromNeighbour);
166 
167  // Combine with this side.
168 
169  if (!facesToSubset.empty())
170  {
171  const labelHashSet nbrPatchFacesToSubset(nbrList);
172 
173  forAll(pp, i)
174  {
175  if
176  (
177  facesToSubset.found(pp.start()+i)
178  && facesToSubset[pp.start()+i] == 1
179  && !nbrPatchFacesToSubset.found(i)
180  )
181  {
182  // Face's neighbour is no longer there. Mark face
183  // off as coupled
184  facesToSubset[pp.start()+i] = 3;
185  nUncoupled++;
186  }
187  }
188  }
189  else if (!nCellsUsingFace.empty())
190  {
191  const labelList& nbrCellsUsingFace(nbrList);
192 
193  // Combine with this side.
194 
195  forAll(pp, i)
196  {
197  if
198  (
199  nCellsUsingFace[pp.start()+i] == 1
200  && nbrCellsUsingFace[i] == 0
201  )
202  {
203  // Face's neighbour is no longer there. Mark face
204  // off as coupled
205  nCellsUsingFace[pp.start()+i] = 3;
206  nUncoupled++;
207  }
208  }
209  }
210  }
211  }
212  }
213 
214  // Do same for cyclics.
215  forAll(oldPatches, oldPatchi)
216  {
217  const polyPatch& pp = oldPatches[oldPatchi];
218 
219  if (isA<cyclicPolyPatch>(pp))
220  {
221  const cyclicPolyPatch& cycPatch =
222  refCast<const cyclicPolyPatch>(pp);
223 
224  if (!facesToSubset.empty())
225  {
226  forAll(cycPatch, i)
227  {
228  label thisFacei = cycPatch.start() + i;
229  label otherFacei = cycPatch.transformGlobalFace(thisFacei);
230 
231  if
232  (
233  facesToSubset.found(thisFacei)
234  && facesToSubset[thisFacei] == 1
235  && !facesToSubset.found(otherFacei)
236  )
237  {
238  facesToSubset[thisFacei] = 3;
239  nUncoupled++;
240  }
241  }
242  }
243  else if (!nCellsUsingFace.empty())
244  {
245  forAll(cycPatch, i)
246  {
247  label thisFacei = cycPatch.start() + i;
248  label otherFacei = cycPatch.transformGlobalFace(thisFacei);
249 
250  if
251  (
252  nCellsUsingFace[thisFacei] == 1
253  && nCellsUsingFace[otherFacei] == 0
254  )
255  {
256  nCellsUsingFace[thisFacei] = 3;
257  nUncoupled++;
258  }
259  }
260  }
261  }
262  }
263 
264  if (syncPar)
265  {
266  reduce(nUncoupled, sumOp<label>());
267  }
268 
269  if (nUncoupled > 0)
270  {
271  Info<< "Uncoupled " << nUncoupled << " faces on coupled patches. "
272  << "(processorPolyPatch, cyclicPolyPatch)" << endl;
273  }
274 }
275 
276 
277 labelList Foam::fvMeshSubset::subset
278 (
279  const label nElems,
280  const labelList& selectedElements,
281  const labelList& subsetMap
282 )
283 {
284  // Mark selected elements.
285  boolList selected(nElems, false);
286  forAll(selectedElements, i)
287  {
288  selected[selectedElements[i]] = true;
289  }
290 
291  // Count subset of selected elements
292  label n = 0;
293  forAll(subsetMap, i)
294  {
295  if (selected[subsetMap[i]])
296  {
297  n++;
298  }
299  }
300 
301  // Collect selected elements
302  labelList subsettedElements(n);
303  n = 0;
304 
305  forAll(subsetMap, i)
306  {
307  if (selected[subsetMap[i]])
308  {
309  subsettedElements[n++] = i;
310  }
311  }
312 
313  return subsettedElements;
314 }
315 
316 
317 void Foam::fvMeshSubset::subsetZones()
318 {
319  // Keep all zones, even if zero size.
320 
321  const pointZoneMesh& pointZones = baseMesh().pointZones();
322 
323  // PointZones
324  List<pointZone*> pZonePtrs(pointZones.size());
325 
326  forAll(pointZones, i)
327  {
328  const pointZone& pz = pointZones[i];
329 
330  pZonePtrs[i] = new pointZone
331  (
332  pz.name(),
333  subset(baseMesh().nPoints(), pz, pointMap()),
334  i,
335  fvMeshSubsetPtr_().pointZones()
336  );
337  }
338 
339 
340  // FaceZones
341 
342  const faceZoneMesh& faceZones = baseMesh().faceZones();
343 
344 
345  // Do we need to remove zones where the side we're interested in
346  // no longer exists? Guess not.
347  List<faceZone*> fZonePtrs(faceZones.size());
348 
349  forAll(faceZones, i)
350  {
351  const faceZone& fz = faceZones[i];
352 
353  // Expand faceZone to full mesh
354  // +1 : part of faceZone, flipped
355  // -1 : ,, , unflipped
356  // 0 : not part of faceZone
357  labelList zone(baseMesh().nFaces(), 0);
358  forAll(fz, j)
359  {
360  if (fz.flipMap()[j])
361  {
362  zone[fz[j]] = 1;
363  }
364  else
365  {
366  zone[fz[j]] = -1;
367  }
368  }
369 
370  // Select faces
371  label nSub = 0;
372  forAll(faceMap(), j)
373  {
374  if (zone[faceMap()[j]] != 0)
375  {
376  nSub++;
377  }
378  }
379  labelList subAddressing(nSub);
380  boolList subFlipStatus(nSub);
381  nSub = 0;
382  forAll(faceMap(), subFacei)
383  {
384  label meshFacei = faceMap()[subFacei];
385  if (zone[meshFacei] != 0)
386  {
387  subAddressing[nSub] = subFacei;
388  label subOwner = subMesh().faceOwner()[subFacei];
389  label baseOwner = baseMesh().faceOwner()[meshFacei];
390  // If subowner is the same cell as the base keep the flip status
391  bool sameOwner = (cellMap()[subOwner] == baseOwner);
392  bool flip = (zone[meshFacei] == 1);
393  subFlipStatus[nSub] = (sameOwner == flip);
394 
395  nSub++;
396  }
397  }
398 
399  fZonePtrs[i] = new faceZone
400  (
401  fz.name(),
402  subAddressing,
403  subFlipStatus,
404  i,
405  fvMeshSubsetPtr_().faceZones()
406  );
407  }
408 
409 
410  const cellZoneMesh& cellZones = baseMesh().cellZones();
411 
412  List<cellZone*> cZonePtrs(cellZones.size());
413 
414  forAll(cellZones, i)
415  {
416  const cellZone& cz = cellZones[i];
417 
418  cZonePtrs[i] = new cellZone
419  (
420  cz.name(),
421  subset(baseMesh().nCells(), cz, cellMap()),
422  i,
423  fvMeshSubsetPtr_().cellZones()
424  );
425  }
426 
427 
428  // Add the zones
429  fvMeshSubsetPtr_().addZones(pZonePtrs, fZonePtrs, cZonePtrs);
430 }
431 
432 
433 Foam::labelList Foam::fvMeshSubset::getCellsToRemove
434 (
435  const labelList& region,
436  const label currentRegion
437 ) const
438 {
439  // Count
440  label nKeep = 0;
441  forAll(region, cellI)
442  {
443  if (region[cellI] == currentRegion)
444  {
445  nKeep++;
446  }
447  }
448 
449  // Collect cells to remove
450  label nRemove = baseMesh().nCells() - nKeep;
451  labelList cellsToRemove(nRemove);
452 
453  nRemove = 0;
454  forAll(region, cellI)
455  {
456  if (region[cellI] != currentRegion)
457  {
458  cellsToRemove[nRemove++] = cellI;
459  }
460  }
461 
462  return cellsToRemove;
463 }
464 
465 
466 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
467 
468 Foam::fvMeshSubset::fvMeshSubset(const fvMesh& baseMesh)
469 :
470  baseMesh_(baseMesh),
471  fvMeshSubsetPtr_(nullptr),
472  pointMap_(0),
473  faceMap_(0),
474  cellMap_(0),
475  patchMap_(0),
476  faceFlipMapPtr_()
477 {}
478 
479 
480 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
481 
483 (
484  const labelHashSet& globalCellMap,
485  const label patchID,
486  const bool syncPar
487 )
488 {
489  // Initial check on patches before doing anything time consuming.
490  const polyBoundaryMesh& oldPatches = baseMesh().boundaryMesh();
491  const cellList& oldCells = baseMesh().cells();
492  const faceList& oldFaces = baseMesh().faces();
493  const pointField& oldPoints = baseMesh().points();
494  const labelList& oldOwner = baseMesh().faceOwner();
495  const labelList& oldNeighbour = baseMesh().faceNeighbour();
496 
497  label wantedPatchID = patchID;
498 
499  if (wantedPatchID == -1)
500  {
501  // No explicit patch specified. Put in oldInternalFaces patch.
502  // Check if patch with this name already exists.
503  wantedPatchID = oldPatches.findPatchID("oldInternalFaces");
504  }
505  else if (wantedPatchID < 0 || wantedPatchID >= oldPatches.size())
506  {
508  << "Non-existing patch index " << wantedPatchID << endl
509  << "Should be between 0 and " << oldPatches.size()-1
510  << abort(FatalError);
511  }
512 
513 
514  // Clear demand driven data
515  faceFlipMapPtr_.clear();
516 
517 
518  cellMap_ = globalCellMap.toc();
519 
520  // Sort the cell map in the ascending order
521  sort(cellMap_);
522 
523  // Approximate sizing parameters for face and point lists
524  const label avgNFacesPerCell = 6;
525  const label avgNPointsPerFace = 4;
526 
527 
528  label nCellsInSet = cellMap_.size();
529 
530  // Mark all used faces
531 
532  Map<label> facesToSubset(avgNFacesPerCell*nCellsInSet);
533 
534  forAll(cellMap_, celli)
535  {
536  // Mark all faces from the cell
537  const labelList& curFaces = oldCells[cellMap_[celli]];
538 
539  forAll(curFaces, facei)
540  {
541  if (!facesToSubset.found(curFaces[facei]))
542  {
543  facesToSubset.insert(curFaces[facei], 1);
544  }
545  else
546  {
547  facesToSubset[curFaces[facei]]++;
548  }
549  }
550  }
551 
552  // Handle coupled faces. Modifies patch faces to be uncoupled to 3.
553  labelList empty;
554  doCoupledPatches(syncPar, facesToSubset, empty);
555 
556  // Mark all used points and make a global-to-local face map
557  Map<label> globalFaceMap(facesToSubset.size());
558 
559  // Make a global-to-local point map
560  Map<label> globalPointMap(avgNPointsPerFace*facesToSubset.size());
561 
562  // This is done in two goes, so that the boundary faces are last
563  // in the list. Because of this, I need to create the face map
564  // along the way rather than just grab the table of contents.
565  labelList facesToc = facesToSubset.toc();
566  sort(facesToc);
567  faceMap_.setSize(facesToc.size());
568 
569  // 1. Get all faces that will be internal to the submesh.
570  forAll(facesToc, facei)
571  {
572  if (facesToSubset[facesToc[facei]] == 2)
573  {
574  // Mark face and increment number of points in set
575  faceMap_[globalFaceMap.size()] = facesToc[facei];
576  globalFaceMap.insert(facesToc[facei], globalFaceMap.size());
577 
578  // Mark all points from the face
579  markPoints(oldFaces[facesToc[facei]], globalPointMap);
580  }
581  }
582 
583  // These are all the internal faces in the mesh.
584  label nInternalFaces = globalFaceMap.size();
585 
586 
587  // Where to insert old internal faces.
588  label oldPatchStart = labelMax;
589  if (wantedPatchID != -1)
590  {
591  oldPatchStart = oldPatches[wantedPatchID].start();
592  }
593 
594 
595  label facei = 0;
596 
597  // 2. Boundary faces up to where we want to insert old internal faces
598  for (; facei< facesToc.size(); facei++)
599  {
600  if (facesToc[facei] >= oldPatchStart)
601  {
602  break;
603  }
604  if
605  (
606  !baseMesh().isInternalFace(facesToc[facei])
607  && facesToSubset[facesToc[facei]] == 1
608  )
609  {
610  // Mark face and increment number of points in set
611  faceMap_[globalFaceMap.size()] = facesToc[facei];
612  globalFaceMap.insert(facesToc[facei], globalFaceMap.size());
613 
614  // Mark all points from the face
615  markPoints(oldFaces[facesToc[facei]], globalPointMap);
616  }
617  }
618 
619  // 3. old internal faces and uncoupled faces
620  forAll(facesToc, intFacei)
621  {
622  if
623  (
624  (
625  baseMesh().isInternalFace(facesToc[intFacei])
626  && facesToSubset[facesToc[intFacei]] == 1
627  )
628  || (
629  !baseMesh().isInternalFace(facesToc[intFacei])
630  && facesToSubset[facesToc[intFacei]] == 3
631  )
632  )
633  {
634  // Mark face and increment number of points in set
635  faceMap_[globalFaceMap.size()] = facesToc[intFacei];
636  globalFaceMap.insert(facesToc[intFacei], globalFaceMap.size());
637 
638  // Mark all points from the face
639  markPoints(oldFaces[facesToc[intFacei]], globalPointMap);
640  }
641  }
642 
643  // 4. Remaining boundary faces
644  for (; facei< facesToc.size(); facei++)
645  {
646  if
647  (
648  !baseMesh().isInternalFace(facesToc[facei])
649  && facesToSubset[facesToc[facei]] == 1
650  )
651  {
652  // Mark face and increment number of points in set
653  faceMap_[globalFaceMap.size()] = facesToc[facei];
654  globalFaceMap.insert(facesToc[facei], globalFaceMap.size());
655 
656  // Mark all points from the face
657  markPoints(oldFaces[facesToc[facei]], globalPointMap);
658  }
659  }
660 
661 
662 
663  // Grab the points map
664  pointMap_ = globalPointMap.toc();
665  sort(pointMap_);
666 
667  forAll(pointMap_, pointi)
668  {
669  globalPointMap[pointMap_[pointi]] = pointi;
670  }
671 
672  // Pout<< "Number of cells in new mesh: " << nCellsInSet << endl;
673  // Pout<< "Number of faces in new mesh: " << globalFaceMap.size() << endl;
674  // Pout<< "Number of points in new mesh: " << globalPointMap.size() << endl;
675 
676  // Make a new mesh
677  pointField newPoints(globalPointMap.size());
678 
679  label nNewPoints = 0;
680 
681  forAll(pointMap_, pointi)
682  {
683  newPoints[nNewPoints] = oldPoints[pointMap_[pointi]];
684  nNewPoints++;
685  }
686 
687  faceList newFaces(globalFaceMap.size());
688 
689  label nNewFaces = 0;
690 
691  // Make internal faces
692  for (label facei = 0; facei < nInternalFaces; facei++)
693  {
694  const face& oldF = oldFaces[faceMap_[facei]];
695 
696  face newF(oldF.size());
697 
698  forAll(newF, i)
699  {
700  newF[i] = globalPointMap[oldF[i]];
701  }
702 
703  newFaces[nNewFaces] = newF;
704  nNewFaces++;
705  }
706 
707  // Make boundary faces
708 
709  label nbSize = oldPatches.size();
710  label oldInternalPatchID = -1;
711 
712  if (wantedPatchID == -1)
713  {
714  // Create 'oldInternalFaces' patch at the end
715  // and put all exposed internal faces in there.
716  oldInternalPatchID = nbSize;
717  nbSize++;
718 
719  }
720  else
721  {
722  oldInternalPatchID = wantedPatchID;
723  }
724 
725 
726  // Grad size and start of each patch on the fly. Because of the
727  // structure of the underlying mesh, the patches will appear in the
728  // ascending order
729  labelList boundaryPatchSizes(nbSize, 0);
730 
731  // Assign boundary faces. Visited in order of faceMap_.
732  for (label facei = nInternalFaces; facei < faceMap_.size(); facei++)
733  {
734  label oldFacei = faceMap_[facei];
735 
736  face oldF = oldFaces[oldFacei];
737 
738  // Turn the faces as necessary to point outwards
739  if (baseMesh().isInternalFace(oldFacei))
740  {
741  // Internal face. Possibly turned the wrong way round
742  if
743  (
744  !globalCellMap.found(oldOwner[oldFacei])
745  && globalCellMap.found(oldNeighbour[oldFacei])
746  )
747  {
748  oldF = oldFaces[oldFacei].reverseFace();
749  }
750 
751  // Update count for patch
752  boundaryPatchSizes[oldInternalPatchID]++;
753  }
754  else if (facesToSubset[oldFacei] == 3)
755  {
756  // Uncoupled face. Increment the old patch.
757  boundaryPatchSizes[oldInternalPatchID]++;
758  }
759  else
760  {
761  // Boundary face. Increment the appropriate patch
762  label patchOfFace = oldPatches.whichPatch(oldFacei);
763 
764  // Update count for patch
765  boundaryPatchSizes[patchOfFace]++;
766  }
767 
768  face newF(oldF.size());
769 
770  forAll(newF, i)
771  {
772  newF[i] = globalPointMap[oldF[i]];
773  }
774 
775  newFaces[nNewFaces] = newF;
776  nNewFaces++;
777  }
778 
779 
780 
781  // Create cells
782  cellList newCells(nCellsInSet);
783 
784  label nNewCells = 0;
785 
786  forAll(cellMap_, celli)
787  {
788  const labelList& oldC = oldCells[cellMap_[celli]];
789 
790  labelList newC(oldC.size());
791 
792  forAll(newC, i)
793  {
794  newC[i] = globalFaceMap[oldC[i]];
795  }
796 
797  newCells[nNewCells] = cell(newC);
798  nNewCells++;
799  }
800 
801 
802  // Delete any old mesh
803  fvMeshSubsetPtr_.clear();
804  // Make a new mesh
805  fvMeshSubsetPtr_.reset
806  (
807  new fvMesh
808  (
809  IOobject
810  (
811  baseMesh().name() + "SubSet",
812  baseMesh().time().timeName(),
813  baseMesh().time(),
816  ),
817  xferMove(newPoints),
818  xferMove(newFaces),
819  xferMove(newCells)
820  )
821  );
822 
823 
824  // Add old patches
825  List<polyPatch*> newBoundary(nbSize);
826  patchMap_.setSize(nbSize);
827  label nNewPatches = 0;
828  label patchStart = nInternalFaces;
829 
830 
831  forAll(oldPatches, patchi)
832  {
833  if (boundaryPatchSizes[patchi] > 0)
834  {
835  // Patch still exists. Add it
836  newBoundary[nNewPatches] = oldPatches[patchi].clone
837  (
838  fvMeshSubsetPtr_().boundaryMesh(),
839  nNewPatches,
840  boundaryPatchSizes[patchi],
841  patchStart
842  ).ptr();
843 
844  patchStart += boundaryPatchSizes[patchi];
845  patchMap_[nNewPatches] = patchi;
846  nNewPatches++;
847  }
848  }
849 
850  if (wantedPatchID == -1)
851  {
852  // Newly created patch so is at end. Check if any faces in it.
853  if (boundaryPatchSizes[oldInternalPatchID] > 0)
854  {
855  newBoundary[nNewPatches] = new emptyPolyPatch
856  (
857  "oldInternalFaces",
858  boundaryPatchSizes[oldInternalPatchID],
859  patchStart,
860  nNewPatches,
861  fvMeshSubsetPtr_().boundaryMesh(),
862  emptyPolyPatch::typeName
863  );
864 
865  // The index for the first patch is -1 as it originates from
866  // the internal faces
867  patchMap_[nNewPatches] = -1;
868  nNewPatches++;
869  }
870  }
871 
872  // Reset the patch lists
873  newBoundary.setSize(nNewPatches);
874  patchMap_.setSize(nNewPatches);
875 
876  // Add the fvPatches
877  fvMeshSubsetPtr_().addFvPatches(newBoundary);
878 
879  // Subset and add any zones
880  subsetZones();
881 }
882 
883 
885 (
886  const labelList& region,
887  const label currentRegion,
888  const label patchID,
889  const bool syncPar
890 )
891 {
892  const cellList& oldCells = baseMesh().cells();
893  const faceList& oldFaces = baseMesh().faces();
894  const pointField& oldPoints = baseMesh().points();
895  const labelList& oldOwner = baseMesh().faceOwner();
896  const labelList& oldNeighbour = baseMesh().faceNeighbour();
897  const polyBoundaryMesh& oldPatches = baseMesh().boundaryMesh();
898  const label oldNInternalFaces = baseMesh().nInternalFaces();
899 
900  // Initial checks
901 
902  if (region.size() != oldCells.size())
903  {
905  << "Size of region " << region.size()
906  << " is not equal to number of cells in mesh " << oldCells.size()
907  << abort(FatalError);
908  }
909 
910 
911  label wantedPatchID = patchID;
912 
913  if (wantedPatchID == -1)
914  {
915  // No explicit patch specified. Put in oldInternalFaces patch.
916  // Check if patch with this name already exists.
917  wantedPatchID = oldPatches.findPatchID("oldInternalFaces");
918  }
919  else if (wantedPatchID < 0 || wantedPatchID >= oldPatches.size())
920  {
922  << "Non-existing patch index " << wantedPatchID << endl
923  << "Should be between 0 and " << oldPatches.size()-1
924  << abort(FatalError);
925  }
926 
927  // Clear demand driven data
928  faceFlipMapPtr_.clear();
929 
930  // Get the cells for the current region.
931  cellMap_.setSize(oldCells.size());
932  label nCellsInSet = 0;
933 
934  forAll(region, oldCelli)
935  {
936  if (region[oldCelli] == currentRegion)
937  {
938  cellMap_[nCellsInSet++] = oldCelli;
939  }
940  }
941  cellMap_.setSize(nCellsInSet);
942 
943 
944  // Mark all used faces. Count number of cells using them
945  // 0: face not used anymore
946  // 1: face used by one cell, face becomes/stays boundary face
947  // 2: face still used and remains internal face
948  // 3: face coupled and used by one cell only (so should become normal,
949  // non-coupled patch face)
950  //
951  // Note that this is not really necessary - but means we can size things
952  // correctly. Also makes handling coupled faces much easier.
953 
954  labelList nCellsUsingFace(oldFaces.size(), 0);
955 
956  label nFacesInSet = 0;
957  forAll(oldFaces, oldFacei)
958  {
959  bool faceUsed = false;
960 
961  if (region[oldOwner[oldFacei]] == currentRegion)
962  {
963  nCellsUsingFace[oldFacei]++;
964  faceUsed = true;
965  }
966 
967  if
968  (
969  baseMesh().isInternalFace(oldFacei)
970  && (region[oldNeighbour[oldFacei]] == currentRegion)
971  )
972  {
973  nCellsUsingFace[oldFacei]++;
974  faceUsed = true;
975  }
976 
977  if (faceUsed)
978  {
979  nFacesInSet++;
980  }
981  }
982  faceMap_.setSize(nFacesInSet);
983 
984  // Handle coupled faces. Modifies patch faces to be uncoupled to 3.
985  Map<label> empty;
986  doCoupledPatches(syncPar, empty, nCellsUsingFace);
987 
988 
989  // See which patch to use for exposed internal faces.
990  label oldInternalPatchID = 0;
991 
992  // Insert faces before which patch
993  label nextPatchID = oldPatches.size();
994 
995  // old to new patches
996  labelList globalPatchMap(oldPatches.size());
997 
998  // New patch size
999  label nbSize = oldPatches.size();
1000 
1001  if (wantedPatchID == -1)
1002  {
1003  // Create 'oldInternalFaces' patch at the end (or before
1004  // processorPatches)
1005  // and put all exposed internal faces in there.
1006 
1007  forAll(oldPatches, patchi)
1008  {
1009  if (isA<processorPolyPatch>(oldPatches[patchi]))
1010  {
1011  nextPatchID = patchi;
1012  break;
1013  }
1014  oldInternalPatchID++;
1015  }
1016 
1017  nbSize++;
1018 
1019  // adapt old to new patches for inserted patch
1020  for (label oldPatchi = 0; oldPatchi < nextPatchID; oldPatchi++)
1021  {
1022  globalPatchMap[oldPatchi] = oldPatchi;
1023  }
1024  for
1025  (
1026  label oldPatchi = nextPatchID;
1027  oldPatchi < oldPatches.size();
1028  oldPatchi++
1029  )
1030  {
1031  globalPatchMap[oldPatchi] = oldPatchi+1;
1032  }
1033  }
1034  else
1035  {
1036  oldInternalPatchID = wantedPatchID;
1037  nextPatchID = wantedPatchID+1;
1038 
1039  // old to new patches
1040  globalPatchMap = identity(oldPatches.size());
1041  }
1042 
1043  labelList boundaryPatchSizes(nbSize, 0);
1044 
1045 
1046  // Make a global-to-local point map
1047  labelList globalPointMap(oldPoints.size(), -1);
1048 
1049  labelList globalFaceMap(oldFaces.size(), -1);
1050  label facei = 0;
1051 
1052  // 1. Pick up all preserved internal faces.
1053  for (label oldFacei = 0; oldFacei < oldNInternalFaces; oldFacei++)
1054  {
1055  if (nCellsUsingFace[oldFacei] == 2)
1056  {
1057  globalFaceMap[oldFacei] = facei;
1058  faceMap_[facei++] = oldFacei;
1059 
1060  // Mark all points from the face
1061  markPoints(oldFaces[oldFacei], globalPointMap);
1062  }
1063  }
1064 
1065  // These are all the internal faces in the mesh.
1066  label nInternalFaces = facei;
1067 
1068  // 2. Boundary faces up to where we want to insert old internal faces
1069  for
1070  (
1071  label oldPatchi = 0;
1072  oldPatchi < oldPatches.size()
1073  && oldPatchi < nextPatchID;
1074  oldPatchi++
1075  )
1076  {
1077  const polyPatch& oldPatch = oldPatches[oldPatchi];
1078 
1079  label oldFacei = oldPatch.start();
1080 
1081  forAll(oldPatch, i)
1082  {
1083  if (nCellsUsingFace[oldFacei] == 1)
1084  {
1085  // Boundary face is kept.
1086 
1087  // Mark face and increment number of points in set
1088  globalFaceMap[oldFacei] = facei;
1089  faceMap_[facei++] = oldFacei;
1090 
1091  // Mark all points from the face
1092  markPoints(oldFaces[oldFacei], globalPointMap);
1093 
1094  // Increment number of patch faces
1095  boundaryPatchSizes[globalPatchMap[oldPatchi]]++;
1096  }
1097  oldFacei++;
1098  }
1099  }
1100 
1101  // 3a. old internal faces that have become exposed.
1102  for (label oldFacei = 0; oldFacei < oldNInternalFaces; oldFacei++)
1103  {
1104  if (nCellsUsingFace[oldFacei] == 1)
1105  {
1106  globalFaceMap[oldFacei] = facei;
1107  faceMap_[facei++] = oldFacei;
1108 
1109  // Mark all points from the face
1110  markPoints(oldFaces[oldFacei], globalPointMap);
1111 
1112  // Increment number of patch faces
1113  boundaryPatchSizes[oldInternalPatchID]++;
1114  }
1115  }
1116 
1117  // 3b. coupled patch faces that have become uncoupled.
1118  for
1119  (
1120  label oldFacei = oldNInternalFaces;
1121  oldFacei < oldFaces.size();
1122  oldFacei++
1123  )
1124  {
1125  if (nCellsUsingFace[oldFacei] == 3)
1126  {
1127  globalFaceMap[oldFacei] = facei;
1128  faceMap_[facei++] = oldFacei;
1129 
1130  // Mark all points from the face
1131  markPoints(oldFaces[oldFacei], globalPointMap);
1132 
1133  // Increment number of patch faces
1134  boundaryPatchSizes[oldInternalPatchID]++;
1135  }
1136  }
1137 
1138  // 4. Remaining boundary faces
1139  for
1140  (
1141  label oldPatchi = nextPatchID;
1142  oldPatchi < oldPatches.size();
1143  oldPatchi++
1144  )
1145  {
1146  const polyPatch& oldPatch = oldPatches[oldPatchi];
1147 
1148  label oldFacei = oldPatch.start();
1149 
1150  forAll(oldPatch, i)
1151  {
1152  if (nCellsUsingFace[oldFacei] == 1)
1153  {
1154  // Boundary face is kept.
1155 
1156  // Mark face and increment number of points in set
1157  globalFaceMap[oldFacei] = facei;
1158  faceMap_[facei++] = oldFacei;
1159 
1160  // Mark all points from the face
1161  markPoints(oldFaces[oldFacei], globalPointMap);
1162 
1163  // Increment number of patch faces
1164  boundaryPatchSizes[globalPatchMap[oldPatchi]]++;
1165  }
1166  oldFacei++;
1167  }
1168  }
1169 
1170  if (facei != nFacesInSet)
1171  {
1173  << "Problem" << abort(FatalError);
1174  }
1175 
1176 
1177  // Grab the points map
1178  label nPointsInSet = 0;
1179 
1180  forAll(globalPointMap, pointi)
1181  {
1182  if (globalPointMap[pointi] != -1)
1183  {
1184  nPointsInSet++;
1185  }
1186  }
1187  pointMap_.setSize(nPointsInSet);
1188 
1189  nPointsInSet = 0;
1190 
1191  forAll(globalPointMap, pointi)
1192  {
1193  if (globalPointMap[pointi] != -1)
1194  {
1195  pointMap_[nPointsInSet] = pointi;
1196  globalPointMap[pointi] = nPointsInSet;
1197  nPointsInSet++;
1198  }
1199  }
1200 
1201  // Pout<< "Number of cells in new mesh : " << cellMap_.size() << endl;
1202  // Pout<< "Number of faces in new mesh : " << faceMap_.size() << endl;
1203  // Pout<< "Number of points in new mesh: " << pointMap_.size() << endl;
1204 
1205  // Make a new mesh
1206  pointField newPoints(pointMap_.size());
1207 
1208  label nNewPoints = 0;
1209 
1210  forAll(pointMap_, pointi)
1211  {
1212  newPoints[nNewPoints] = oldPoints[pointMap_[pointi]];
1213  nNewPoints++;
1214  }
1215 
1216  faceList newFaces(faceMap_.size());
1217 
1218  label nNewFaces = 0;
1219 
1220  // Make internal faces
1221  for (label facei = 0; facei < nInternalFaces; facei++)
1222  {
1223  const face& oldF = oldFaces[faceMap_[facei]];
1224 
1225  face newF(oldF.size());
1226 
1227  forAll(newF, i)
1228  {
1229  newF[i] = globalPointMap[oldF[i]];
1230  }
1231 
1232  newFaces[nNewFaces] = newF;
1233  nNewFaces++;
1234  }
1235 
1236 
1237  // Make boundary faces. (different from internal since might need to be
1238  // flipped)
1239  for (label facei = nInternalFaces; facei < faceMap_.size(); facei++)
1240  {
1241  label oldFacei = faceMap_[facei];
1242 
1243  face oldF = oldFaces[oldFacei];
1244 
1245  // Turn the faces as necessary to point outwards
1246  if (baseMesh().isInternalFace(oldFacei))
1247  {
1248  // Was internal face. Possibly turned the wrong way round
1249  if
1250  (
1251  region[oldOwner[oldFacei]] != currentRegion
1252  && region[oldNeighbour[oldFacei]] == currentRegion
1253  )
1254  {
1255  oldF = oldFaces[oldFacei].reverseFace();
1256  }
1257  }
1258 
1259  // Relabel vertices of the (possibly turned) face.
1260  face newF(oldF.size());
1261 
1262  forAll(newF, i)
1263  {
1264  newF[i] = globalPointMap[oldF[i]];
1265  }
1266 
1267  newFaces[nNewFaces] = newF;
1268  nNewFaces++;
1269  }
1270 
1271 
1272 
1273  // Create cells
1274  cellList newCells(nCellsInSet);
1275 
1276  label nNewCells = 0;
1277 
1278  forAll(cellMap_, celli)
1279  {
1280  const labelList& oldC = oldCells[cellMap_[celli]];
1281 
1282  labelList newC(oldC.size());
1283 
1284  forAll(newC, i)
1285  {
1286  newC[i] = globalFaceMap[oldC[i]];
1287  }
1288 
1289  newCells[nNewCells] = cell(newC);
1290  nNewCells++;
1291  }
1292 
1293 
1294  // Delete any old one
1295  fvMeshSubsetPtr_.clear();
1296 
1297  // Make a new mesh
1298  // Note that mesh gets registered with same name as original mesh. This is
1299  // not proper but cannot be avoided since otherwise surfaceInterpolation
1300  // cannot find its fvSchemes (it will try to read e.g.
1301  // system/region0SubSet/fvSchemes)
1302  // Make a new mesh
1303  fvMeshSubsetPtr_.reset
1304  (
1305  new fvMesh
1306  (
1307  IOobject
1308  (
1309  baseMesh().name(),
1310  baseMesh().time().timeName(),
1311  baseMesh().time(),
1314  ),
1315  xferMove(newPoints),
1316  xferMove(newFaces),
1317  xferMove(newCells),
1318  syncPar // parallel synchronisation
1319  )
1320  );
1321 
1322  // Add old patches
1323  List<polyPatch*> newBoundary(nbSize);
1324  patchMap_.setSize(nbSize);
1325  label nNewPatches = 0;
1326  label patchStart = nInternalFaces;
1327 
1328 
1329  // For parallel: only remove patch if none of the processors has it.
1330  // This only gets done for patches before the one being inserted
1331  // (so patches < nextPatchID)
1332 
1333  // Get sum of patch sizes. Zero if patch can be deleted.
1334  labelList globalPatchSizes(boundaryPatchSizes);
1335  globalPatchSizes.setSize(nextPatchID);
1336 
1337  if (syncPar && Pstream::parRun())
1338  {
1339  // Get patch names (up to nextPatchID)
1341  patchNames[Pstream::myProcNo()] = oldPatches.names();
1342  patchNames[Pstream::myProcNo()].setSize(nextPatchID);
1343  Pstream::gatherList(patchNames);
1344  Pstream::scatterList(patchNames);
1345 
1346  // Get patch sizes (up to nextPatchID).
1347  // Note that up to nextPatchID the globalPatchMap is an identity so
1348  // no need to index through that.
1349  Pstream::listCombineGather(globalPatchSizes, plusEqOp<label>());
1350  Pstream::listCombineScatter(globalPatchSizes);
1351 
1352  // Now all processors have all the patchnames.
1353  // Decide: if all processors have the same patch names and size is zero
1354  // everywhere remove the patch.
1355  bool samePatches = true;
1356 
1357  for (label proci = 1; proci < patchNames.size(); proci++)
1358  {
1359  if (patchNames[proci] != patchNames[0])
1360  {
1361  samePatches = false;
1362  break;
1363  }
1364  }
1365 
1366  if (!samePatches)
1367  {
1368  // Patchnames not sync on all processors so disable removal of
1369  // zero sized patches.
1370  globalPatchSizes = labelMax;
1371  }
1372  }
1373 
1374 
1375  // Old patches
1376 
1377  for
1378  (
1379  label oldPatchi = 0;
1380  oldPatchi < oldPatches.size()
1381  && oldPatchi < nextPatchID;
1382  oldPatchi++
1383  )
1384  {
1385  label newSize = boundaryPatchSizes[globalPatchMap[oldPatchi]];
1386 
1387  // Clone (even if 0 size)
1388  newBoundary[nNewPatches] = oldPatches[oldPatchi].clone
1389  (
1390  fvMeshSubsetPtr_().boundaryMesh(),
1391  nNewPatches,
1392  newSize,
1393  patchStart
1394  ).ptr();
1395 
1396  patchStart += newSize;
1397  patchMap_[nNewPatches] = oldPatchi; // compact patchMap
1398  nNewPatches++;
1399  }
1400 
1401  // Inserted patch
1402 
1403  if (wantedPatchID == -1)
1404  {
1405  label oldInternalSize = boundaryPatchSizes[oldInternalPatchID];
1406 
1407  if (syncPar)
1408  {
1409  reduce(oldInternalSize, sumOp<label>());
1410  }
1411 
1412  // Newly created patch so is at end. Check if any faces in it.
1413  if (oldInternalSize > 0)
1414  {
1415  newBoundary[nNewPatches] = new emptyPolyPatch
1416  (
1417  "oldInternalFaces",
1418  boundaryPatchSizes[oldInternalPatchID],
1419  patchStart,
1420  nNewPatches,
1421  fvMeshSubsetPtr_().boundaryMesh(),
1422  emptyPolyPatch::typeName
1423  );
1424 
1425  // Pout<< " oldInternalFaces : "
1426  // << boundaryPatchSizes[oldInternalPatchID] << endl;
1427 
1428  // The index for the first patch is -1 as it originates from
1429  // the internal faces
1430  patchStart += boundaryPatchSizes[oldInternalPatchID];
1431  patchMap_[nNewPatches] = -1;
1432  nNewPatches++;
1433  }
1434  }
1435 
1436  // Old patches
1437 
1438  for
1439  (
1440  label oldPatchi = nextPatchID;
1441  oldPatchi < oldPatches.size();
1442  oldPatchi++
1443  )
1444  {
1445  label newSize = boundaryPatchSizes[globalPatchMap[oldPatchi]];
1446 
1447  // Patch still exists. Add it
1448  newBoundary[nNewPatches] = oldPatches[oldPatchi].clone
1449  (
1450  fvMeshSubsetPtr_().boundaryMesh(),
1451  nNewPatches,
1452  newSize,
1453  patchStart
1454  ).ptr();
1455 
1456  // Pout<< " " << oldPatches[oldPatchi].name() << " : "
1457  // << newSize << endl;
1458 
1459  patchStart += newSize;
1460  patchMap_[nNewPatches] = oldPatchi; // compact patchMap
1461  nNewPatches++;
1462  }
1463 
1464 
1465  // Reset the patch lists
1466  newBoundary.setSize(nNewPatches);
1467  patchMap_.setSize(nNewPatches);
1468 
1469 
1470  // Add the fvPatches
1471  fvMeshSubsetPtr_().addFvPatches(newBoundary, syncPar);
1472 
1473  // Subset and add any zones
1474  subsetZones();
1475 }
1476 
1477 
1480  const labelHashSet& globalCellMap,
1481  const label patchID,
1482  const bool syncPar
1483 )
1484 {
1485  labelList region(baseMesh().nCells(), 0);
1486 
1487  forAllConstIter(labelHashSet, globalCellMap, iter)
1488  {
1489  region[iter.key()] = 1;
1490  }
1491  setLargeCellSubset(region, 1, patchID, syncPar);
1492 }
1493 
1494 
1497  const labelList& region,
1498  const label currentRegion,
1499  const bool syncCouples
1500 ) const
1501 {
1502  // Collect cells to remove
1503  labelList cellsToRemove(getCellsToRemove(region, currentRegion));
1504 
1505  return removeCells(baseMesh(), syncCouples).getExposedFaces(cellsToRemove);
1506 }
1507 
1508 
1511  const labelList& region,
1512  const label currentRegion,
1513  const labelList& exposedFaces,
1514  const labelList& patchIDs,
1515  const bool syncCouples
1516 )
1517 {
1518  // Collect cells to remove
1519  labelList cellsToRemove(getCellsToRemove(region, currentRegion));
1520 
1521  // Mesh changing engine.
1522  polyTopoChange meshMod(baseMesh());
1523 
1524  removeCells cellRemover(baseMesh(), syncCouples);
1525 
1526  cellRemover.setRefinement
1527  (
1528  cellsToRemove,
1529  exposedFaces,
1530  patchIDs,
1531  meshMod
1532  );
1533 
1534  // Create mesh, return map from old to new mesh.
1535  autoPtr<mapPolyMesh> map = meshMod.makeMesh
1536  (
1537  fvMeshSubsetPtr_,
1538  IOobject
1539  (
1540  baseMesh().name(),
1541  baseMesh().time().timeName(),
1542  baseMesh().time(),
1545  ),
1546  baseMesh(),
1547  syncCouples
1548  );
1549 
1550  pointMap_ = map().pointMap();
1551  faceMap_ = map().faceMap();
1552  cellMap_ = map().cellMap();
1553  patchMap_ = identity(baseMesh().boundaryMesh().size());
1554 }
1555 
1556 
1558 {
1559  return fvMeshSubsetPtr_.valid();
1560 }
1561 
1562 
1564 {
1565  checkCellSubset();
1566 
1567  return fvMeshSubsetPtr_();
1568 }
1569 
1570 
1572 {
1573  checkCellSubset();
1574 
1575  return fvMeshSubsetPtr_();
1576 }
1577 
1578 
1580 {
1581  checkCellSubset();
1582 
1583  return pointMap_;
1584 }
1585 
1586 
1588 {
1589  checkCellSubset();
1590 
1591  return faceMap_;
1592 }
1593 
1594 
1596 {
1597  if (!faceFlipMapPtr_.valid())
1598  {
1599  const labelList& subToBaseFace = faceMap();
1600  const labelList& subToBaseCell = cellMap();
1601 
1602  faceFlipMapPtr_.reset(new labelList(subToBaseFace.size()));
1603  labelList& faceFlipMap = faceFlipMapPtr_();
1604 
1605  // Only exposed internal faces might be flipped (since we don't do
1606  // any cell renumbering, just compacting)
1607  label subInt = subMesh().nInternalFaces();
1608  const labelList& subOwn = subMesh().faceOwner();
1609  const labelList& own = baseMesh_.faceOwner();
1610 
1611  for (label subFaceI = 0; subFaceI < subInt; subFaceI++)
1612  {
1613  faceFlipMap[subFaceI] = subToBaseFace[subFaceI]+1;
1614  }
1615  for (label subFaceI = subInt; subFaceI < subOwn.size(); subFaceI++)
1616  {
1617  label faceI = subToBaseFace[subFaceI];
1618  if (subToBaseCell[subOwn[subFaceI]] == own[faceI])
1619  {
1620  faceFlipMap[subFaceI] = faceI+1;
1621  }
1622  else
1623  {
1624  faceFlipMap[subFaceI] = -faceI-1;
1625  }
1626  }
1627  }
1628 
1629  return faceFlipMapPtr_();
1630 }
1631 
1632 
1634 {
1635  checkCellSubset();
1636 
1637  return cellMap_;
1638 }
1639 
1640 
1642 {
1643  checkCellSubset();
1644 
1645  return patchMap_;
1646 }
1647 
1648 
1649 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1650 
1651 } // End namespace Foam
1652 
1653 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
const labelList & patchMap() const
Return patch map.
autoPtr< IOobject > clone() const
Clone.
Definition: IOobject.H:270
ZoneMesh< faceZone, polyMesh > faceZoneMesh
A ZoneMesh with the type faceZone.
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#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
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
void setRefinement(const labelList &cellsToRemove, const labelList &facesToExpose, const labelList &patchIDs, polyTopoChange &) const
Play commands into polyTopoChange to remove cells.
Definition: removeCells.C:179
error FatalError
autoPtr< mapPolyMesh > makeMesh(autoPtr< fvMesh > &newMesh, const IOobject &io, const polyMesh &mesh, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Create new mesh with old mesh patches.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
label nInternalFaces() const
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1047
const labelList & faceFlipMap() const
Return face map with sign to encode flipped faces.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:427
Given list of cells to remove insert all the topology changes.
Definition: removeCells.H:59
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
label findPatchID(const word &patchName) const
Find patch index given a name.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
const cellList & cells() const
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1003
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:210
labelList getExposedFaces(const labelList &region, const label currentRegion, const bool syncCouples=true) const
Two step subsetting.
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:113
Xfer< T > xferMove(T &)
Construct by transferring the contents of the arg.
label nPoints
void sort(UList< T > &)
Definition: UList.C:115
wordList names() const
Return a list of patch names.
wordList patchNames(nPatches)
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1041
const labelList & faceMap() const
Return face map.
static const label labelMax
Definition: label.H:62
word timeName
Definition: getTimeIndex.H:3
List< label > labelList
A List of labels.
Definition: labelList.H:56
ZoneMesh< pointZone, polyMesh > pointZoneMesh
A ZoneMesh with the type pointZone.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1028
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Foam::polyBoundaryMesh.
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
void setCellSubset(const labelHashSet &globalCellMap, const label patchID=-1, const bool syncPar=true)
Set the subset. Create "oldInternalFaces" patch for exposed.
Definition: fvMeshSubset.C:483
Empty front and back plane patch. Used for 2-D geometries.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
void setSize(const label)
Reset size of List.
Definition: List.C:281
Template functions to aid in the implementation of demand driven data.
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:397
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:409
label patchi
const fvMesh & subMesh() const
Return reference to subset mesh.
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Direct mesh changes based on v1.3 polyTopoChange syntax.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
ListType subset(const UList< T > &select, const T &value, const ListType &)
Extract elements of List when select is a certain value.
messageStream Info
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
labelList getExposedFaces(const labelList &cellsToRemove) const
Get labels of exposed faces.
Definition: removeCells.C:73
List< Key > toc() const
Return the table of contents.
Definition: HashTable.C:202
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
void setLargeCellSubset(const labelList &region, const label currentRegion, const label patchID=-1, const bool syncCouples=true)
Set the subset from all cells with region == currentRegion.
Definition: fvMeshSubset.C:885
const labelList & cellMap() const
Return cell map.
const fvMesh & baseMesh() const
Original mesh.
Definition: fvMeshSubset.H:219
ZoneMesh< cellZone, polyMesh > cellZoneMesh
A ZoneMesh with the type cellZone.
const labelList & pointMap() const
Return point map.
Namespace for OpenFOAM.
bool hasSubMesh() const
Have subMesh?