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