fvMeshDistribute.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-2014 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "fvMeshDistribute.H"
28 #include "fvMeshAdder.H"
29 #include "faceCoupleInfo.H"
30 #include "processorFvPatchField.H"
31 #include "processorFvsPatchField.H"
34 #include "polyTopoChange.H"
35 #include "removeCells.H"
36 #include "polyModifyFace.H"
37 #include "polyRemovePoint.H"
38 #include "mapDistributePolyMesh.H"
39 #include "surfaceFields.H"
40 #include "syncTools.H"
41 #include "CompactListList.H"
42 #include "fvMeshTools.H"
43 
44 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 defineTypeNameAndDebug(fvMeshDistribute, 0);
49 }
50 
51 
52 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
53 
54 Foam::labelList Foam::fvMeshDistribute::select
55 (
56  const bool selectEqual,
57  const labelList& values,
58  const label value
59 )
60 {
61  label n = 0;
62 
63  forAll(values, i)
64  {
65  if (selectEqual == (values[i] == value))
66  {
67  n++;
68  }
69  }
70 
71  labelList indices(n);
72  n = 0;
73 
74  forAll(values, i)
75  {
76  if (selectEqual == (values[i] == value))
77  {
78  indices[n++] = i;
79  }
80  }
81  return indices;
82 }
83 
84 
85 // Check all procs have same names and in exactly same order.
86 void Foam::fvMeshDistribute::checkEqualWordList
87 (
88  const string& msg,
89  const wordList& lst
90 )
91 {
92  List<wordList> allNames(Pstream::nProcs());
93  allNames[Pstream::myProcNo()] = lst;
94  Pstream::gatherList(allNames);
95  Pstream::scatterList(allNames);
96 
97  for (label procI = 1; procI < Pstream::nProcs(); procI++)
98  {
99  if (allNames[procI] != allNames[0])
100  {
101  FatalErrorIn("fvMeshDistribute::checkEqualWordList(..)")
102  << "When checking for equal " << msg.c_str() << " :" << endl
103  << "processor0 has:" << allNames[0] << endl
104  << "processor" << procI << " has:" << allNames[procI] << endl
105  << msg.c_str() << " need to be synchronised on all processors."
106  << exit(FatalError);
107  }
108  }
109 }
110 
111 
112 Foam::wordList Foam::fvMeshDistribute::mergeWordList(const wordList& procNames)
113 {
114  List<wordList> allNames(Pstream::nProcs());
115  allNames[Pstream::myProcNo()] = procNames;
116  Pstream::gatherList(allNames);
117  Pstream::scatterList(allNames);
118 
119  HashSet<word> mergedNames;
120  forAll(allNames, procI)
121  {
122  forAll(allNames[procI], i)
123  {
124  mergedNames.insert(allNames[procI][i]);
125  }
126  }
127  return mergedNames.toc();
128 }
129 
130 
131 // Print some info on mesh.
133 {
134  Pout<< "Primitives:" << nl
135  << " points :" << mesh.nPoints() << nl
136  << " bb :" << boundBox(mesh.points(), false) << nl
137  << " internalFaces:" << mesh.nInternalFaces() << nl
138  << " faces :" << mesh.nFaces() << nl
139  << " cells :" << mesh.nCells() << nl;
140 
141  const fvBoundaryMesh& patches = mesh.boundary();
142 
143  Pout<< "Patches:" << endl;
144  forAll(patches, patchI)
145  {
146  const polyPatch& pp = patches[patchI].patch();
147 
148  Pout<< " " << patchI << " name:" << pp.name()
149  << " size:" << pp.size()
150  << " start:" << pp.start()
151  << " type:" << pp.type()
152  << endl;
153  }
154 
155  if (mesh.pointZones().size())
156  {
157  Pout<< "PointZones:" << endl;
158  forAll(mesh.pointZones(), zoneI)
159  {
160  const pointZone& pz = mesh.pointZones()[zoneI];
161  Pout<< " " << zoneI << " name:" << pz.name()
162  << " size:" << pz.size()
163  << endl;
164  }
165  }
166  if (mesh.faceZones().size())
167  {
168  Pout<< "FaceZones:" << endl;
169  forAll(mesh.faceZones(), zoneI)
170  {
171  const faceZone& fz = mesh.faceZones()[zoneI];
172  Pout<< " " << zoneI << " name:" << fz.name()
173  << " size:" << fz.size()
174  << endl;
175  }
176  }
177  if (mesh.cellZones().size())
178  {
179  Pout<< "CellZones:" << endl;
180  forAll(mesh.cellZones(), zoneI)
181  {
182  const cellZone& cz = mesh.cellZones()[zoneI];
183  Pout<< " " << zoneI << " name:" << cz.name()
184  << " size:" << cz.size()
185  << endl;
186  }
187  }
188 }
189 
190 
192 (
193  const primitiveMesh& mesh,
194  const labelList& sourceFace,
195  const labelList& sourceProc,
196  const labelList& sourcePatch,
197  const labelList& sourceNewNbrProc
198 )
199 {
200  Pout<< nl
201  << "Current coupling info:"
202  << endl;
203 
204  forAll(sourceFace, bFaceI)
205  {
206  label meshFaceI = mesh.nInternalFaces() + bFaceI;
207 
208  Pout<< " meshFace:" << meshFaceI
209  << " fc:" << mesh.faceCentres()[meshFaceI]
210  << " connects to proc:" << sourceProc[bFaceI]
211  << "/face:" << sourceFace[bFaceI]
212  << " which will move to proc:" << sourceNewNbrProc[bFaceI]
213  << endl;
214  }
215 }
216 
217 
218 // Finds (non-empty) patch that exposed internal and proc faces can be put into.
219 Foam::label Foam::fvMeshDistribute::findNonEmptyPatch() const
220 {
221  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
222 
223  label nonEmptyPatchI = -1;
224 
225  forAllReverse(patches, patchI)
226  {
227  const polyPatch& pp = patches[patchI];
228 
229  if (!isA<emptyPolyPatch>(pp) && !pp.coupled())
230  {
231  nonEmptyPatchI = patchI;
232  break;
233  }
234  }
235 
236  if (nonEmptyPatchI == -1)
237  {
238  FatalErrorIn("fvMeshDistribute::findNonEmptyPatch() const")
239  << "Cannot find a patch which is neither of type empty nor"
240  << " coupled in patches " << patches.names() << endl
241  << "There has to be at least one such patch for"
242  << " distribution to work" << abort(FatalError);
243  }
244 
245  if (debug)
246  {
247  Pout<< "findNonEmptyPatch : using patch " << nonEmptyPatchI
248  << " name:" << patches[nonEmptyPatchI].name()
249  << " type:" << patches[nonEmptyPatchI].type()
250  << " to put exposed faces into." << endl;
251  }
252 
253 
254  // Do additional test for processor patches intermingled with non-proc
255  // patches.
256  label procPatchI = -1;
257 
258  forAll(patches, patchI)
259  {
260  if (isA<processorPolyPatch>(patches[patchI]))
261  {
262  procPatchI = patchI;
263  }
264  else if (procPatchI != -1)
265  {
266  FatalErrorIn("fvMeshDistribute::findNonEmptyPatch() const")
267  << "Processor patches should be at end of patch list."
268  << endl
269  << "Have processor patch " << procPatchI
270  << " followed by non-processor patch " << patchI
271  << " in patches " << patches.names()
272  << abort(FatalError);
273  }
274  }
275 
276  return nonEmptyPatchI;
277 }
278 
279 
280 // Delete all processor patches. Move any processor faces into the last
281 // non-processor patch.
282 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::deleteProcPatches
283 (
284  const label destinationPatch
285 )
286 {
287  // New patchID per boundary faces to be repatched. Is -1 (no change)
288  // or new patchID
289  labelList newPatchID(mesh_.nFaces() - mesh_.nInternalFaces(), -1);
290 
291  label nProcPatches = 0;
292 
293  forAll(mesh_.boundaryMesh(), patchI)
294  {
295  const polyPatch& pp = mesh_.boundaryMesh()[patchI];
296 
297  if (isA<processorPolyPatch>(pp))
298  {
299  if (debug)
300  {
301  Pout<< "Moving all faces of patch " << pp.name()
302  << " into patch " << destinationPatch
303  << endl;
304  }
305 
306  label offset = pp.start() - mesh_.nInternalFaces();
307 
308  forAll(pp, i)
309  {
310  newPatchID[offset+i] = destinationPatch;
311  }
312 
313  nProcPatches++;
314  }
315  }
316 
317  // Note: order of boundary faces been kept the same since the
318  // destinationPatch is at the end and we have visited the patches in
319  // incremental order.
320  labelListList dummyFaceMaps;
321  autoPtr<mapPolyMesh> map = repatch(newPatchID, dummyFaceMaps);
322 
323 
324  // Delete (now empty) processor patches.
325  {
326  labelList oldToNew(identity(mesh_.boundaryMesh().size()));
327  label newI = 0;
328  // Non processor patches first
329  forAll(mesh_.boundaryMesh(), patchI)
330  {
331  if (!isA<processorPolyPatch>(mesh_.boundaryMesh()[patchI]))
332  {
333  oldToNew[patchI] = newI++;
334  }
335  }
336  label nNonProcPatches = newI;
337 
338  // Processor patches as last
339  forAll(mesh_.boundaryMesh(), patchI)
340  {
341  if (isA<processorPolyPatch>(mesh_.boundaryMesh()[patchI]))
342  {
343  oldToNew[patchI] = newI++;
344  }
345  }
346  fvMeshTools::reorderPatches(mesh_, oldToNew, nNonProcPatches, false);
347  }
348 
349  return map;
350 }
351 
352 
353 // Repatch the mesh.
354 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::repatch
355 (
356  const labelList& newPatchID, // per boundary face -1 or new patchID
357  labelListList& constructFaceMap
358 )
359 {
360  polyTopoChange meshMod(mesh_);
361 
362  forAll(newPatchID, bFaceI)
363  {
364  if (newPatchID[bFaceI] != -1)
365  {
366  label faceI = mesh_.nInternalFaces() + bFaceI;
367 
368  label zoneID = mesh_.faceZones().whichZone(faceI);
369  bool zoneFlip = false;
370 
371  if (zoneID >= 0)
372  {
373  const faceZone& fZone = mesh_.faceZones()[zoneID];
374  zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)];
375  }
376 
377  meshMod.setAction
378  (
380  (
381  mesh_.faces()[faceI], // modified face
382  faceI, // label of face
383  mesh_.faceOwner()[faceI], // owner
384  -1, // neighbour
385  false, // face flip
386  newPatchID[bFaceI], // patch for face
387  false, // remove from zone
388  zoneID, // zone for face
389  zoneFlip // face flip in zone
390  )
391  );
392  }
393  }
394 
395 
396  // Do mapping of fields from one patchField to the other ourselves since
397  // is currently not supported by updateMesh.
398 
399  // Store boundary fields (we only do this for surfaceFields)
401  saveBoundaryFields<scalar, surfaceMesh>(sFlds);
403  saveBoundaryFields<vector, surfaceMesh>(vFlds);
405  saveBoundaryFields<sphericalTensor, surfaceMesh>(sptFlds);
407  saveBoundaryFields<symmTensor, surfaceMesh>(sytFlds);
409  saveBoundaryFields<tensor, surfaceMesh>(tFlds);
410 
411  // Change the mesh (no inflation). Note: parallel comms allowed.
412  //
413  // NOTE: there is one very particular problem with this ordering.
414  // We first create the processor patches and use these to merge out
415  // shared points (see mergeSharedPoints below). So temporarily points
416  // and edges do not match!
417 
418  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, true);
419 
420  // Update fields. No inflation, parallel sync.
421  mesh_.updateMesh(map);
422 
423  // Map patch fields using stored boundary fields. Note: assumes order
424  // of fields has not changed in object registry!
425  mapBoundaryFields<scalar, surfaceMesh>(map, sFlds);
426  mapBoundaryFields<vector, surfaceMesh>(map, vFlds);
427  mapBoundaryFields<sphericalTensor, surfaceMesh>(map, sptFlds);
428  mapBoundaryFields<symmTensor, surfaceMesh>(map, sytFlds);
429  mapBoundaryFields<tensor, surfaceMesh>(map, tFlds);
430 
431 
432  // Move mesh (since morphing does not do this)
433  if (map().hasMotionPoints())
434  {
435  mesh_.movePoints(map().preMotionPoints());
436  }
437 
438  // Adapt constructMaps.
439 
440  if (debug)
441  {
442  label index = findIndex(map().reverseFaceMap(), -1);
443 
444  if (index != -1)
445  {
447  (
448  "fvMeshDistribute::repatch(const labelList&, labelListList&)"
449  ) << "reverseFaceMap contains -1 at index:"
450  << index << endl
451  << "This means that the repatch operation was not just"
452  << " a shuffle?" << abort(FatalError);
453  }
454  }
455 
456  forAll(constructFaceMap, procI)
457  {
458  inplaceRenumber(map().reverseFaceMap(), constructFaceMap[procI]);
459  }
460 
461 
462  return map;
463 }
464 
465 
466 // Detect shared points. Need processor patches to be present.
467 // Background: when adding bits of mesh one can get points which
468 // share the same position but are only detectable to be topologically
469 // the same point when doing parallel analysis. This routine will
470 // merge those points.
471 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::mergeSharedPoints
472 (
473  labelListList& constructPointMap
474 )
475 {
476  // Find out which sets of points get merged and create a map from
477  // mesh point to unique point.
478  Map<label> pointToMaster
479  (
481  (
482  mesh_,
483  mergeTol_
484  )
485  );
486 
487  if (returnReduce(pointToMaster.size(), sumOp<label>()) == 0)
488  {
489  return autoPtr<mapPolyMesh>(NULL);
490  }
491 
492  polyTopoChange meshMod(mesh_);
493 
494  fvMeshAdder::mergePoints(mesh_, pointToMaster, meshMod);
495 
496  // Change the mesh (no inflation). Note: parallel comms allowed.
497  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, true);
498 
499  // Update fields. No inflation, parallel sync.
500  mesh_.updateMesh(map);
501 
502  // Adapt constructMaps for merged points.
503  forAll(constructPointMap, procI)
504  {
505  labelList& constructMap = constructPointMap[procI];
506 
507  forAll(constructMap, i)
508  {
509  label oldPointI = constructMap[i];
510 
511  label newPointI = map().reversePointMap()[oldPointI];
512 
513  if (newPointI < -1)
514  {
515  constructMap[i] = -newPointI-2;
516  }
517  else if (newPointI >= 0)
518  {
519  constructMap[i] = newPointI;
520  }
521  else
522  {
523  FatalErrorIn("fvMeshDistribute::mergeSharedPoints()")
524  << "Problem. oldPointI:" << oldPointI
525  << " newPointI:" << newPointI << abort(FatalError);
526  }
527  }
528  }
529  return map;
530 }
531 
532 
533 // Construct the local environment of all boundary faces.
534 void Foam::fvMeshDistribute::getNeighbourData
535 (
536  const labelList& distribution,
537  labelList& sourceFace,
538  labelList& sourceProc,
539  labelList& sourcePatch,
540  labelList& sourceNewNbrProc
541 ) const
542 {
543  label nBnd = mesh_.nFaces() - mesh_.nInternalFaces();
544  sourceFace.setSize(nBnd);
545  sourceProc.setSize(nBnd);
546  sourcePatch.setSize(nBnd);
547  sourceNewNbrProc.setSize(nBnd);
548 
549  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
550 
551  // Get neighbouring meshFace labels and new processor of coupled boundaries.
552  labelList nbrFaces(nBnd, -1);
553  labelList nbrNewNbrProc(nBnd, -1);
554 
555  forAll(patches, patchI)
556  {
557  const polyPatch& pp = patches[patchI];
558 
559  if (pp.coupled())
560  {
561  label offset = pp.start() - mesh_.nInternalFaces();
562 
563  // Mesh labels of faces on this side
564  forAll(pp, i)
565  {
566  label bndI = offset + i;
567  nbrFaces[bndI] = pp.start()+i;
568  }
569 
570  // Which processor they will end up on
571  SubList<label>(nbrNewNbrProc, pp.size(), offset).assign
572  (
573  UIndirectList<label>(distribution, pp.faceCells())()
574  );
575  }
576  }
577 
578 
579  // Exchange the boundary data
580  syncTools::swapBoundaryFaceList(mesh_, nbrFaces);
581  syncTools::swapBoundaryFaceList(mesh_, nbrNewNbrProc);
582 
583 
584  forAll(patches, patchI)
585  {
586  const polyPatch& pp = patches[patchI];
587  label offset = pp.start() - mesh_.nInternalFaces();
588 
589  if (isA<processorPolyPatch>(pp))
590  {
591  const processorPolyPatch& procPatch =
592  refCast<const processorPolyPatch>(pp);
593 
594  // Check which of the two faces we store.
595 
596  if (procPatch.owner())
597  {
598  // Use my local face labels
599  forAll(pp, i)
600  {
601  label bndI = offset + i;
602  sourceFace[bndI] = pp.start()+i;
603  sourceProc[bndI] = Pstream::myProcNo();
604  sourceNewNbrProc[bndI] = nbrNewNbrProc[bndI];
605  }
606  }
607  else
608  {
609  // Use my neighbours face labels
610  forAll(pp, i)
611  {
612  label bndI = offset + i;
613  sourceFace[bndI] = nbrFaces[bndI];
614  sourceProc[bndI] = procPatch.neighbProcNo();
615  sourceNewNbrProc[bndI] = nbrNewNbrProc[bndI];
616  }
617  }
618 
619 
620  label patchI = -1;
621  if (isA<processorCyclicPolyPatch>(pp))
622  {
623  patchI = refCast<const processorCyclicPolyPatch>
624  (
625  pp
626  ).referPatchID();
627  }
628 
629  forAll(pp, i)
630  {
631  label bndI = offset + i;
632  sourcePatch[bndI] = patchI;
633  }
634  }
635  else if (isA<cyclicPolyPatch>(pp))
636  {
637  const cyclicPolyPatch& cpp = refCast<const cyclicPolyPatch>(pp);
638 
639  if (cpp.owner())
640  {
641  forAll(pp, i)
642  {
643  label bndI = offset + i;
644  sourceFace[bndI] = pp.start()+i;
645  sourceProc[bndI] = Pstream::myProcNo();
646  sourcePatch[bndI] = patchI;
647  sourceNewNbrProc[bndI] = nbrNewNbrProc[bndI];
648  }
649  }
650  else
651  {
652  forAll(pp, i)
653  {
654  label bndI = offset + i;
655  sourceFace[bndI] = nbrFaces[bndI];
656  sourceProc[bndI] = Pstream::myProcNo();
657  sourcePatch[bndI] = patchI;
658  sourceNewNbrProc[bndI] = nbrNewNbrProc[bndI];
659  }
660  }
661  }
662  else
663  {
664  // Normal (physical) boundary
665  forAll(pp, i)
666  {
667  label bndI = offset + i;
668  sourceFace[bndI] = -1;
669  sourceProc[bndI] = -1;
670  sourcePatch[bndI] = patchI;
671  sourceNewNbrProc[bndI] = -1;
672  }
673  }
674  }
675 }
676 
677 
678 // Subset the neighbourCell/neighbourProc fields
679 void Foam::fvMeshDistribute::subsetBoundaryData
680 (
681  const fvMesh& mesh,
682  const labelList& faceMap,
683  const labelList& cellMap,
684 
685  const labelList& oldDistribution,
686  const labelList& oldFaceOwner,
687  const labelList& oldFaceNeighbour,
688  const label oldInternalFaces,
689 
690  const labelList& sourceFace,
691  const labelList& sourceProc,
692  const labelList& sourcePatch,
693  const labelList& sourceNewNbrProc,
694 
695  labelList& subFace,
696  labelList& subProc,
697  labelList& subPatch,
698  labelList& subNewNbrProc
699 )
700 {
701  subFace.setSize(mesh.nFaces() - mesh.nInternalFaces());
702  subProc.setSize(mesh.nFaces() - mesh.nInternalFaces());
703  subPatch.setSize(mesh.nFaces() - mesh.nInternalFaces());
704  subNewNbrProc.setSize(mesh.nFaces() - mesh.nInternalFaces());
705 
706  forAll(subFace, newBFaceI)
707  {
708  label newFaceI = newBFaceI + mesh.nInternalFaces();
709 
710  label oldFaceI = faceMap[newFaceI];
711 
712  // Was oldFaceI internal face? If so which side did we get.
713  if (oldFaceI < oldInternalFaces)
714  {
715  subFace[newBFaceI] = oldFaceI;
716  subProc[newBFaceI] = Pstream::myProcNo();
717  subPatch[newBFaceI] = -1;
718 
719  label oldOwn = oldFaceOwner[oldFaceI];
720  label oldNei = oldFaceNeighbour[oldFaceI];
721 
722  if (oldOwn == cellMap[mesh.faceOwner()[newFaceI]])
723  {
724  // We kept the owner side. Where does the neighbour move to?
725  subNewNbrProc[newBFaceI] = oldDistribution[oldNei];
726  }
727  else
728  {
729  // We kept the neighbour side.
730  subNewNbrProc[newBFaceI] = oldDistribution[oldOwn];
731  }
732  }
733  else
734  {
735  // Was boundary face. Take over boundary information
736  label oldBFaceI = oldFaceI - oldInternalFaces;
737 
738  subFace[newBFaceI] = sourceFace[oldBFaceI];
739  subProc[newBFaceI] = sourceProc[oldBFaceI];
740  subPatch[newBFaceI] = sourcePatch[oldBFaceI];
741  subNewNbrProc[newBFaceI] = sourceNewNbrProc[oldBFaceI];
742  }
743  }
744 }
745 
746 
747 // Find cells on mesh whose faceID/procID match the neighbour cell/proc of
748 // domainMesh. Store the matching face.
749 void Foam::fvMeshDistribute::findCouples
750 (
751  const primitiveMesh& mesh,
752  const labelList& sourceFace,
753  const labelList& sourceProc,
754  const labelList& sourcePatch,
755 
756  const label domain,
757  const primitiveMesh& domainMesh,
758  const labelList& domainFace,
759  const labelList& domainProc,
760  const labelList& domainPatch,
761 
762  labelList& masterCoupledFaces,
763  labelList& slaveCoupledFaces
764 )
765 {
766  // Store domain neighbour as map so we can easily look for pair
767  // with same face+proc.
769 
770  forAll(domainProc, bFaceI)
771  {
772  if (domainProc[bFaceI] != -1 && domainPatch[bFaceI] == -1)
773  {
774  map.insert
775  (
776  labelPair(domainFace[bFaceI], domainProc[bFaceI]),
777  bFaceI
778  );
779  }
780  }
781 
782 
783  // Try to match mesh data.
784 
785  masterCoupledFaces.setSize(domainFace.size());
786  slaveCoupledFaces.setSize(domainFace.size());
787  label coupledI = 0;
788 
789  forAll(sourceFace, bFaceI)
790  {
791  if (sourceProc[bFaceI] != -1 && sourcePatch[bFaceI] == -1)
792  {
793  labelPair myData(sourceFace[bFaceI], sourceProc[bFaceI]);
794 
796  iter = map.find(myData);
797 
798  if (iter != map.end())
799  {
800  label nbrBFaceI = iter();
801 
802  masterCoupledFaces[coupledI] = mesh.nInternalFaces() + bFaceI;
803  slaveCoupledFaces[coupledI] =
804  domainMesh.nInternalFaces()
805  + nbrBFaceI;
806 
807  coupledI++;
808  }
809  }
810  }
811 
812  masterCoupledFaces.setSize(coupledI);
813  slaveCoupledFaces.setSize(coupledI);
814 
815  if (debug)
816  {
817  Pout<< "findCouples : found " << coupledI
818  << " faces that will be stitched" << nl << endl;
819  }
820 }
821 
822 
823 // Map data on boundary faces to new mesh (resulting from adding two meshes)
824 Foam::labelList Foam::fvMeshDistribute::mapBoundaryData
825 (
826  const primitiveMesh& mesh, // mesh after adding
827  const mapAddedPolyMesh& map,
828  const labelList& boundaryData0, // mesh before adding
829  const label nInternalFaces1,
830  const labelList& boundaryData1 // added mesh
831 )
832 {
833  labelList newBoundaryData(mesh.nFaces() - mesh.nInternalFaces());
834 
835  forAll(boundaryData0, oldBFaceI)
836  {
837  label newFaceI = map.oldFaceMap()[oldBFaceI + map.nOldInternalFaces()];
838 
839  // Face still exists (is necessary?) and still boundary face
840  if (newFaceI >= 0 && newFaceI >= mesh.nInternalFaces())
841  {
842  newBoundaryData[newFaceI - mesh.nInternalFaces()] =
843  boundaryData0[oldBFaceI];
844  }
845  }
846 
847  forAll(boundaryData1, addedBFaceI)
848  {
849  label newFaceI = map.addedFaceMap()[addedBFaceI + nInternalFaces1];
850 
851  if (newFaceI >= 0 && newFaceI >= mesh.nInternalFaces())
852  {
853  newBoundaryData[newFaceI - mesh.nInternalFaces()] =
854  boundaryData1[addedBFaceI];
855  }
856  }
857 
858  return newBoundaryData;
859 }
860 
861 
862 // Remove cells. Add all exposed faces to patch oldInternalPatchI
863 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::doRemoveCells
864 (
865  const labelList& cellsToRemove,
866  const label oldInternalPatchI
867 )
868 {
869  // Mesh change engine
870  polyTopoChange meshMod(mesh_);
871 
872  // Cell removal topo engine. Do NOT synchronize parallel since
873  // we are doing a local cell removal.
874  removeCells cellRemover(mesh_, false);
875 
876  // Get all exposed faces
877  labelList exposedFaces(cellRemover.getExposedFaces(cellsToRemove));
878 
879  // Insert the topo changes
880  cellRemover.setRefinement
881  (
882  cellsToRemove,
883  exposedFaces,
884  labelList(exposedFaces.size(), oldInternalPatchI), // patch for exposed
885  // faces.
886  meshMod
887  );
888 
889  // Change the mesh. No inflation. Note: no parallel comms allowed.
890  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, false);
891 
892  // Update fields
893  mesh_.updateMesh(map);
894 
895  // Move mesh (since morphing does not do this)
896  if (map().hasMotionPoints())
897  {
898  mesh_.movePoints(map().preMotionPoints());
899  }
900 
901  return map;
902 }
903 
904 
905 // Delete and add processor patches. Changes mesh and returns per neighbour proc
906 // the processor patchID.
907 void Foam::fvMeshDistribute::addProcPatches
908 (
909  const labelList& nbrProc, // processor that neighbour is now on
910  const labelList& referPatchID, // patchID (or -1) I originated from
911  List<Map<label> >& procPatchID
912 )
913 {
914  // Now use the neighbourFace/Proc to repatch the mesh. These lists
915  // contain for all current boundary faces the global patchID (for non-proc
916  // patch) or the processor.
917 
918  procPatchID.setSize(Pstream::nProcs());
919 
920  forAll(nbrProc, bFaceI)
921  {
922  label procI = nbrProc[bFaceI];
923 
924  if (procI != -1 && procI != Pstream::myProcNo())
925  {
926  if (!procPatchID[procI].found(referPatchID[bFaceI]))
927  {
928  // No patch for neighbour yet. Is either a normal processor
929  // patch or a processorCyclic patch.
930 
931  if (referPatchID[bFaceI] == -1)
932  {
933  // Ordinary processor boundary
934 
935  const word patchName =
936  "procBoundary"
938  + "to"
939  + name(procI);
940 
942  (
943  patchName,
944  0, // size
945  mesh_.nFaces(),
946  mesh_.boundaryMesh().size(),
947  mesh_.boundaryMesh(),
949  nbrProc[bFaceI]
950  );
951 
952  procPatchID[procI].insert
953  (
954  referPatchID[bFaceI],
956  (
957  mesh_,
958  pp,
959  dictionary(), // optional per field patchField
961  false // not parallel sync
962  )
963  );
964  }
965  else
966  {
967  const coupledPolyPatch& pcPatch
968  = refCast<const coupledPolyPatch>
969  (
970  mesh_.boundaryMesh()[referPatchID[bFaceI]]
971  );
972 
973  // Processor boundary originating from cyclic
974  const word& cycName = pcPatch.name();
975 
976  const word patchName =
977  "procBoundary"
979  + "to"
980  + name(procI)
981  + "through"
982  + cycName;
983 
985  (
986  patchName,
987  0, // size
988  mesh_.nFaces(),
989  mesh_.boundaryMesh().size(),
990  mesh_.boundaryMesh(),
992  nbrProc[bFaceI],
993  cycName,
994  pcPatch.transform()
995  );
996 
997  procPatchID[procI].insert
998  (
999  referPatchID[bFaceI],
1001  (
1002  mesh_,
1003  pp,
1004  dictionary(), // optional per field patchField
1006  false // not parallel sync
1007  )
1008  );
1009  }
1010  }
1011  }
1012  }
1013 }
1014 
1015 
1016 // Get boundary faces to be repatched. Is -1 or new patchID
1017 Foam::labelList Foam::fvMeshDistribute::getBoundaryPatch
1018 (
1019  const labelList& nbrProc, // new processor per boundary face
1020  const labelList& referPatchID, // patchID (or -1) I originated from
1021  const List<Map<label> >& procPatchID // per proc the new procPatches
1022 )
1023 {
1024  labelList patchIDs(nbrProc);
1025 
1026  forAll(nbrProc, bFaceI)
1027  {
1028  if (nbrProc[bFaceI] == Pstream::myProcNo())
1029  {
1030  label origPatchI = referPatchID[bFaceI];
1031  patchIDs[bFaceI] = origPatchI;
1032  }
1033  else if (nbrProc[bFaceI] != -1)
1034  {
1035  label origPatchI = referPatchID[bFaceI];
1036  patchIDs[bFaceI] = procPatchID[nbrProc[bFaceI]][origPatchI];
1037  }
1038  else
1039  {
1040  patchIDs[bFaceI] = -1;
1041  }
1042  }
1043  return patchIDs;
1044 }
1045 
1046 
1047 // Send mesh and coupling data.
1048 void Foam::fvMeshDistribute::sendMesh
1049 (
1050  const label domain,
1051  const fvMesh& mesh,
1052 
1053  const wordList& pointZoneNames,
1054  const wordList& faceZoneNames,
1055  const wordList& cellZoneNames,
1056 
1057  const labelList& sourceFace,
1058  const labelList& sourceProc,
1059  const labelList& sourcePatch,
1060  const labelList& sourceNewNbrProc,
1061  Ostream& toDomain
1062 )
1063 {
1064  if (debug)
1065  {
1066  Pout<< "Sending to domain " << domain << nl
1067  << " nPoints:" << mesh.nPoints() << nl
1068  << " nFaces:" << mesh.nFaces() << nl
1069  << " nCells:" << mesh.nCells() << nl
1070  << " nPatches:" << mesh.boundaryMesh().size() << nl
1071  << endl;
1072  }
1073 
1074  // Assume sparse, possibly overlapping point zones. Get contents
1075  // in merged-zone indices.
1076  CompactListList<label> zonePoints;
1077  {
1078  const pointZoneMesh& pointZones = mesh.pointZones();
1079 
1080  labelList rowSizes(pointZoneNames.size(), 0);
1081 
1082  forAll(pointZoneNames, nameI)
1083  {
1084  label myZoneID = pointZones.findZoneID(pointZoneNames[nameI]);
1085 
1086  if (myZoneID != -1)
1087  {
1088  rowSizes[nameI] = pointZones[myZoneID].size();
1089  }
1090  }
1091  zonePoints.setSize(rowSizes);
1092 
1093  forAll(pointZoneNames, nameI)
1094  {
1095  label myZoneID = pointZones.findZoneID(pointZoneNames[nameI]);
1096 
1097  if (myZoneID != -1)
1098  {
1099  zonePoints[nameI].assign(pointZones[myZoneID]);
1100  }
1101  }
1102  }
1103 
1104  // Assume sparse, possibly overlapping face zones
1105  CompactListList<label> zoneFaces;
1106  CompactListList<bool> zoneFaceFlip;
1107  {
1108  const faceZoneMesh& faceZones = mesh.faceZones();
1109 
1110  labelList rowSizes(faceZoneNames.size(), 0);
1111 
1112  forAll(faceZoneNames, nameI)
1113  {
1114  label myZoneID = faceZones.findZoneID(faceZoneNames[nameI]);
1115 
1116  if (myZoneID != -1)
1117  {
1118  rowSizes[nameI] = faceZones[myZoneID].size();
1119  }
1120  }
1121 
1122  zoneFaces.setSize(rowSizes);
1123  zoneFaceFlip.setSize(rowSizes);
1124 
1125  forAll(faceZoneNames, nameI)
1126  {
1127  label myZoneID = faceZones.findZoneID(faceZoneNames[nameI]);
1128 
1129  if (myZoneID != -1)
1130  {
1131  zoneFaces[nameI].assign(faceZones[myZoneID]);
1132  zoneFaceFlip[nameI].assign(faceZones[myZoneID].flipMap());
1133  }
1134  }
1135  }
1136 
1137  // Assume sparse, possibly overlapping cell zones
1138  CompactListList<label> zoneCells;
1139  {
1140  const cellZoneMesh& cellZones = mesh.cellZones();
1141 
1142  labelList rowSizes(cellZoneNames.size(), 0);
1143 
1144  forAll(cellZoneNames, nameI)
1145  {
1146  label myZoneID = cellZones.findZoneID(cellZoneNames[nameI]);
1147 
1148  if (myZoneID != -1)
1149  {
1150  rowSizes[nameI] = cellZones[myZoneID].size();
1151  }
1152  }
1153 
1154  zoneCells.setSize(rowSizes);
1155 
1156  forAll(cellZoneNames, nameI)
1157  {
1158  label myZoneID = cellZones.findZoneID(cellZoneNames[nameI]);
1159 
1160  if (myZoneID != -1)
1161  {
1162  zoneCells[nameI].assign(cellZones[myZoneID]);
1163  }
1164  }
1165  }
1167  //labelList cellZoneID;
1168  //if (hasCellZones)
1169  //{
1170  // cellZoneID.setSize(mesh.nCells());
1171  // cellZoneID = -1;
1172  //
1173  // const cellZoneMesh& cellZones = mesh.cellZones();
1174  //
1175  // forAll(cellZones, zoneI)
1176  // {
1177  // UIndirectList<label>(cellZoneID, cellZones[zoneI]) = zoneI;
1178  // }
1179  //}
1180 
1181  // Send
1182  toDomain
1183  << mesh.points()
1185  << mesh.faceOwner()
1186  << mesh.faceNeighbour()
1187  << mesh.boundaryMesh()
1188 
1189  << zonePoints
1190  << zoneFaces
1191  << zoneFaceFlip
1192  << zoneCells
1193 
1194  << sourceFace
1195  << sourceProc
1196  << sourcePatch
1197  << sourceNewNbrProc;
1198 
1199 
1200  if (debug)
1201  {
1202  Pout<< "Started sending mesh to domain " << domain
1203  << endl;
1204  }
1205 }
1206 
1207 
1208 // Receive mesh. Opposite of sendMesh
1209 Foam::autoPtr<Foam::fvMesh> Foam::fvMeshDistribute::receiveMesh
1210 (
1211  const label domain,
1212  const wordList& pointZoneNames,
1213  const wordList& faceZoneNames,
1214  const wordList& cellZoneNames,
1215  const Time& runTime,
1216  labelList& domainSourceFace,
1217  labelList& domainSourceProc,
1218  labelList& domainSourcePatch,
1219  labelList& domainSourceNewNbrProc,
1220  Istream& fromNbr
1221 )
1222 {
1223  pointField domainPoints(fromNbr);
1224  faceList domainFaces = CompactListList<label, face>(fromNbr)();
1225  labelList domainAllOwner(fromNbr);
1226  labelList domainAllNeighbour(fromNbr);
1227  PtrList<entry> patchEntries(fromNbr);
1228 
1229  CompactListList<label> zonePoints(fromNbr);
1230  CompactListList<label> zoneFaces(fromNbr);
1231  CompactListList<bool> zoneFaceFlip(fromNbr);
1232  CompactListList<label> zoneCells(fromNbr);
1233 
1234  fromNbr
1235  >> domainSourceFace
1236  >> domainSourceProc
1237  >> domainSourcePatch
1238  >> domainSourceNewNbrProc;
1239 
1240  // Construct fvMesh
1241  autoPtr<fvMesh> domainMeshPtr
1242  (
1243  new fvMesh
1244  (
1245  IOobject
1246  (
1248  runTime.timeName(),
1249  runTime,
1251  ),
1252  xferMove(domainPoints),
1253  xferMove(domainFaces),
1254  xferMove(domainAllOwner),
1255  xferMove(domainAllNeighbour),
1256  false // no parallel comms
1257  )
1258  );
1259  fvMesh& domainMesh = domainMeshPtr();
1260 
1261  List<polyPatch*> patches(patchEntries.size());
1262 
1263  forAll(patchEntries, patchI)
1264  {
1265  patches[patchI] = polyPatch::New
1266  (
1267  patchEntries[patchI].keyword(),
1268  patchEntries[patchI].dict(),
1269  patchI,
1270  domainMesh.boundaryMesh()
1271  ).ptr();
1272  }
1273  // Add patches; no parallel comms
1274  domainMesh.addFvPatches(patches, false);
1275 
1276  // Construct zones
1277  List<pointZone*> pZonePtrs(pointZoneNames.size());
1278  forAll(pZonePtrs, i)
1279  {
1280  pZonePtrs[i] = new pointZone
1281  (
1282  pointZoneNames[i],
1283  zonePoints[i],
1284  i,
1285  domainMesh.pointZones()
1286  );
1287  }
1288 
1289  List<faceZone*> fZonePtrs(faceZoneNames.size());
1290  forAll(fZonePtrs, i)
1291  {
1292  fZonePtrs[i] = new faceZone
1293  (
1294  faceZoneNames[i],
1295  zoneFaces[i],
1296  zoneFaceFlip[i],
1297  i,
1298  domainMesh.faceZones()
1299  );
1300  }
1301 
1302  List<cellZone*> cZonePtrs(cellZoneNames.size());
1303  forAll(cZonePtrs, i)
1304  {
1305  cZonePtrs[i] = new cellZone
1306  (
1307  cellZoneNames[i],
1308  zoneCells[i],
1309  i,
1310  domainMesh.cellZones()
1311  );
1312  }
1313  domainMesh.addZones(pZonePtrs, fZonePtrs, cZonePtrs);
1314 
1315  return domainMeshPtr;
1316 }
1317 
1318 
1319 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
1320 
1321 // Construct from components
1322 Foam::fvMeshDistribute::fvMeshDistribute(fvMesh& mesh, const scalar mergeTol)
1323 :
1324  mesh_(mesh),
1325  mergeTol_(mergeTol)
1326 {}
1327 
1328 
1329 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1330 
1333  const labelList& distribution
1334 )
1335 {
1336  labelList nCells(Pstream::nProcs(), 0);
1337  forAll(distribution, cellI)
1338  {
1339  label newProc = distribution[cellI];
1340 
1341  if (newProc < 0 || newProc >= Pstream::nProcs())
1342  {
1343  FatalErrorIn("fvMeshDistribute::distribute(const labelList&)")
1344  << "Distribution should be in range 0.." << Pstream::nProcs()-1
1345  << endl
1346  << "At index " << cellI << " distribution:" << newProc
1347  << abort(FatalError);
1348  }
1349  nCells[newProc]++;
1350  }
1351  return nCells;
1352 }
1353 
1354 
1357  const labelList& distribution
1358 )
1359 {
1360  // Some checks on distribution
1361  if (distribution.size() != mesh_.nCells())
1362  {
1363  FatalErrorIn("fvMeshDistribute::distribute(const labelList&)")
1364  << "Size of distribution:"
1365  << distribution.size() << " mesh nCells:" << mesh_.nCells()
1366  << abort(FatalError);
1367  }
1368 
1369 
1370  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
1371 
1372  // Check all processors have same non-proc patches in same order.
1373  if (patches.checkParallelSync(true))
1374  {
1375  FatalErrorIn("fvMeshDistribute::distribute(const labelList&)")
1376  << "This application requires all non-processor patches"
1377  << " to be present in the same order on all patches" << nl
1378  << "followed by the processor patches (which of course are unique)."
1379  << nl
1380  << "Local patches:" << mesh_.boundaryMesh().names()
1381  << abort(FatalError);
1382  }
1383 
1384  // Save some data for mapping later on
1385  const label nOldPoints(mesh_.nPoints());
1386  const label nOldFaces(mesh_.nFaces());
1387  const label nOldCells(mesh_.nCells());
1388  labelList oldPatchStarts(patches.size());
1389  labelList oldPatchNMeshPoints(patches.size());
1390  forAll(patches, patchI)
1391  {
1392  oldPatchStarts[patchI] = patches[patchI].start();
1393  oldPatchNMeshPoints[patchI] = patches[patchI].nPoints();
1394  }
1395 
1396 
1397 
1398  // Short circuit trivial case.
1399  if (!Pstream::parRun())
1400  {
1401  // Collect all maps and return
1403  (
1405  (
1406  mesh_,
1407 
1408  nOldPoints,
1409  nOldFaces,
1410  nOldCells,
1411  oldPatchStarts.xfer(),
1412  oldPatchNMeshPoints.xfer(),
1413 
1414  labelListList(1, identity(mesh_.nPoints())).xfer(),//subPointMap
1415  labelListList(1, identity(mesh_.nFaces())).xfer(), //subFaceMap
1416  labelListList(1, identity(mesh_.nCells())).xfer(), //subCellMap
1417  labelListList(1, identity(patches.size())).xfer(), //subPatchMap
1418 
1419  labelListList(1, identity(mesh_.nPoints())).xfer(),//pointMap
1420  labelListList(1, identity(mesh_.nFaces())).xfer(), //faceMap
1421  labelListList(1, identity(mesh_.nCells())).xfer(), //cellMap
1422  labelListList(1, identity(patches.size())).xfer() //patchMap
1423  )
1424  );
1425  }
1426 
1427 
1428  // Collect any zone names
1429  const wordList pointZoneNames(mergeWordList(mesh_.pointZones().names()));
1430  const wordList faceZoneNames(mergeWordList(mesh_.faceZones().names()));
1431  const wordList cellZoneNames(mergeWordList(mesh_.cellZones().names()));
1432 
1433 
1434 
1435  // Local environment of all boundary faces
1436  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1437 
1438  // A face is uniquely defined by
1439  // - proc
1440  // - local face no
1441  //
1442  // To glue the parts of meshes which can get sent from anywhere we
1443  // need to know on boundary faces what the above tuple on both sides is.
1444  // So we need to maintain:
1445  // - original face
1446  // - original processor id (= trivial)
1447  // For coupled boundaries (where the faces are 'duplicate') we take the
1448  // lowest numbered processor as the data to store.
1449  //
1450  // Additionally to create the procboundaries we need to know where the owner
1451  // cell on the other side now is: newNeighbourProc.
1452  //
1453 
1454  // physical boundary:
1455  // sourceProc = -1
1456  // sourceNewNbrProc = -1
1457  // sourceFace = -1
1458  // sourcePatch = patchID
1459  // processor boundary:
1460  // sourceProc = proc (on owner side)
1461  // sourceNewNbrProc = distribution of coupled cell
1462  // sourceFace = face (on owner side)
1463  // sourcePatch = -1
1464  // ?cyclic:
1465  // ? sourceProc = proc
1466  // ? sourceNewNbrProc = distribution of coupled cell
1467  // ? sourceFace = face (on owner side)
1468  // ? sourcePatch = patchID
1469  // processor-cyclic boundary:
1470  // sourceProc = proc (on owner side)
1471  // sourceNewNbrProc = distribution of coupled cell
1472  // sourceFace = face (on owner side)
1473  // sourcePatch = patchID
1474 
1475  labelList sourcePatch;
1476  labelList sourceFace;
1477  labelList sourceProc;
1478  labelList sourceNewNbrProc;
1479  getNeighbourData
1480  (
1481  distribution,
1482  sourceFace,
1483  sourceProc,
1484  sourcePatch,
1485  sourceNewNbrProc
1486  );
1487 
1488 
1489  // Remove meshPhi. Since this would otherwise disappear anyway
1490  // during topo changes and we have to guarantee that all the fields
1491  // can be sent.
1492  mesh_.clearOut();
1493  mesh_.resetMotion();
1494 
1495  // Get data to send. Make sure is synchronised
1496  const wordList volScalars(mesh_.names(volScalarField::typeName));
1497  checkEqualWordList("volScalarFields", volScalars);
1498  const wordList volVectors(mesh_.names(volVectorField::typeName));
1499  checkEqualWordList("volVectorFields", volVectors);
1500  const wordList volSphereTensors
1501  (
1503  );
1504  checkEqualWordList("volSphericalTensorFields", volSphereTensors);
1505  const wordList volSymmTensors(mesh_.names(volSymmTensorField::typeName));
1506  checkEqualWordList("volSymmTensorFields", volSymmTensors);
1507  const wordList volTensors(mesh_.names(volTensorField::typeName));
1508  checkEqualWordList("volTensorField", volTensors);
1509 
1510  const wordList surfScalars(mesh_.names(surfaceScalarField::typeName));
1511  checkEqualWordList("surfaceScalarFields", surfScalars);
1512  const wordList surfVectors(mesh_.names(surfaceVectorField::typeName));
1513  checkEqualWordList("surfaceVectorFields", surfVectors);
1514  const wordList surfSphereTensors
1515  (
1517  );
1518  checkEqualWordList("surfaceSphericalTensorFields", surfSphereTensors);
1519  const wordList surfSymmTensors
1520  (
1522  );
1523  checkEqualWordList("surfaceSymmTensorFields", surfSymmTensors);
1524  const wordList surfTensors(mesh_.names(surfaceTensorField::typeName));
1525  checkEqualWordList("surfaceTensorFields", surfTensors);
1526 
1527 
1528 
1529 
1530  // Find patch to temporarily put exposed and processor faces into.
1531  label oldInternalPatchI = findNonEmptyPatch();
1532 
1533 
1534 
1535  // Delete processor patches, starting from the back. Move all faces into
1536  // oldInternalPatchI.
1537  labelList repatchFaceMap;
1538  {
1539  autoPtr<mapPolyMesh> repatchMap = deleteProcPatches(oldInternalPatchI);
1540 
1541  // Store face map (only face ordering that changed)
1542  repatchFaceMap = repatchMap().faceMap();
1543 
1544  // Reorder all boundary face data (sourceProc, sourceFace etc.)
1545  labelList bFaceMap
1546  (
1548  (
1549  repatchMap().reverseFaceMap(),
1550  mesh_.nFaces() - mesh_.nInternalFaces(),
1551  mesh_.nInternalFaces()
1552  )
1553  - mesh_.nInternalFaces()
1554  );
1555 
1556  inplaceReorder(bFaceMap, sourceFace);
1557  inplaceReorder(bFaceMap, sourceProc);
1558  inplaceReorder(bFaceMap, sourcePatch);
1559  inplaceReorder(bFaceMap, sourceNewNbrProc);
1560  }
1561 
1562 
1563 
1564  // Print a bit.
1565  if (debug)
1566  {
1567  Pout<< nl << "MESH WITH PROC PATCHES DELETED:" << endl;
1568  printMeshInfo(mesh_);
1569  printFieldInfo<volScalarField>(mesh_);
1570  printFieldInfo<volVectorField>(mesh_);
1571  printFieldInfo<volSphericalTensorField>(mesh_);
1572  printFieldInfo<volSymmTensorField>(mesh_);
1573  printFieldInfo<volTensorField>(mesh_);
1574  printFieldInfo<surfaceScalarField>(mesh_);
1575  printFieldInfo<surfaceVectorField>(mesh_);
1576  printFieldInfo<surfaceSphericalTensorField>(mesh_);
1577  printFieldInfo<surfaceSymmTensorField>(mesh_);
1578  printFieldInfo<surfaceTensorField>(mesh_);
1579  Pout<< nl << endl;
1580  }
1581 
1582 
1583 
1584  // Maps from subsetted mesh (that is sent) back to original maps
1585  labelListList subCellMap(Pstream::nProcs());
1586  labelListList subFaceMap(Pstream::nProcs());
1587  labelListList subPointMap(Pstream::nProcs());
1588  labelListList subPatchMap(Pstream::nProcs());
1589  // Maps from subsetted mesh to reconstructed mesh
1590  labelListList constructCellMap(Pstream::nProcs());
1591  labelListList constructFaceMap(Pstream::nProcs());
1592  labelListList constructPointMap(Pstream::nProcs());
1593  labelListList constructPatchMap(Pstream::nProcs());
1594 
1595 
1596 
1597 
1598  // Find out schedule
1599  // ~~~~~~~~~~~~~~~~~
1600 
1601  labelListList nSendCells(Pstream::nProcs());
1602  nSendCells[Pstream::myProcNo()] = countCells(distribution);
1603  Pstream::gatherList(nSendCells);
1604  Pstream::scatterList(nSendCells);
1605 
1606 
1607  // Allocate buffers
1609 
1610 
1611  // What to send to neighbouring domains
1612  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1613 
1614  bool oldParRun = UPstream::parRun();
1615  UPstream::parRun() = false;
1616 
1617  forAll(nSendCells[Pstream::myProcNo()], recvProc)
1618  {
1619  if
1620  (
1621  recvProc != Pstream::myProcNo()
1622  && nSendCells[Pstream::myProcNo()][recvProc] > 0
1623  )
1624  {
1625  // Send to recvProc
1626 
1627  if (debug)
1628  {
1629  Pout<< nl
1630  << "SUBSETTING FOR DOMAIN " << recvProc
1631  << " cells to send:"
1632  << nSendCells[Pstream::myProcNo()][recvProc]
1633  << nl << endl;
1634  }
1635 
1636  // Pstream for sending mesh and fields
1637  //OPstream str(Pstream::blocking, recvProc);
1638  UOPstream str(recvProc, pBufs);
1639 
1640  // Mesh subsetting engine
1641  fvMeshSubset subsetter(mesh_);
1642 
1643  // Subset the cells of the current domain.
1644  subsetter.setLargeCellSubset
1645  (
1646  distribution,
1647  recvProc,
1648  oldInternalPatchI, // oldInternalFaces patch
1649  false // no parallel sync
1650  );
1651 
1652  subCellMap[recvProc] = subsetter.cellMap();
1653  subFaceMap[recvProc] = renumber
1654  (
1655  repatchFaceMap,
1656  subsetter.faceMap()
1657  );
1658  subPointMap[recvProc] = subsetter.pointMap();
1659  subPatchMap[recvProc] = subsetter.patchMap();
1660 
1661 
1662  // Subset the boundary fields (owner/neighbour/processor)
1663  labelList procSourceFace;
1664  labelList procSourceProc;
1665  labelList procSourcePatch;
1666  labelList procSourceNewNbrProc;
1667 
1668  subsetBoundaryData
1669  (
1670  subsetter.subMesh(),
1671  subsetter.faceMap(), // from subMesh to mesh
1672  subsetter.cellMap(), // ,, ,,
1673 
1674  distribution, // old mesh distribution
1675  mesh_.faceOwner(), // old owner
1676  mesh_.faceNeighbour(),
1677  mesh_.nInternalFaces(),
1678 
1679  sourceFace,
1680  sourceProc,
1681  sourcePatch,
1682  sourceNewNbrProc,
1683 
1684  procSourceFace,
1685  procSourceProc,
1686  procSourcePatch,
1687  procSourceNewNbrProc
1688  );
1689 
1690 
1691 
1692  // Send to neighbour
1693  sendMesh
1694  (
1695  recvProc,
1696  subsetter.subMesh(),
1697 
1698  pointZoneNames,
1699  faceZoneNames,
1700  cellZoneNames,
1701 
1702  procSourceFace,
1703  procSourceProc,
1704  procSourcePatch,
1705  procSourceNewNbrProc,
1706  str
1707  );
1708  sendFields<volScalarField>(recvProc, volScalars, subsetter, str);
1709  sendFields<volVectorField>(recvProc, volVectors, subsetter, str);
1710  sendFields<volSphericalTensorField>
1711  (
1712  recvProc,
1713  volSphereTensors,
1714  subsetter,
1715  str
1716  );
1717  sendFields<volSymmTensorField>
1718  (
1719  recvProc,
1720  volSymmTensors,
1721  subsetter,
1722  str
1723  );
1724  sendFields<volTensorField>(recvProc, volTensors, subsetter, str);
1725 
1726  sendFields<surfaceScalarField>
1727  (
1728  recvProc,
1729  surfScalars,
1730  subsetter,
1731  str
1732  );
1733  sendFields<surfaceVectorField>
1734  (
1735  recvProc,
1736  surfVectors,
1737  subsetter,
1738  str
1739  );
1740  sendFields<surfaceSphericalTensorField>
1741  (
1742  recvProc,
1743  surfSphereTensors,
1744  subsetter,
1745  str
1746  );
1747  sendFields<surfaceSymmTensorField>
1748  (
1749  recvProc,
1750  surfSymmTensors,
1751  subsetter,
1752  str
1753  );
1754  sendFields<surfaceTensorField>
1755  (
1756  recvProc,
1757  surfTensors,
1758  subsetter,
1759  str
1760  );
1761  }
1762  }
1763 
1764 
1765  UPstream::parRun() = oldParRun;
1766 
1767 
1768  // Start sending&receiving from buffers
1769  pBufs.finishedSends();
1770 
1771 
1772  // Subset the part that stays
1773  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
1774 
1775  {
1776  // Save old mesh maps before changing mesh
1777  const labelList oldFaceOwner(mesh_.faceOwner());
1778  const labelList oldFaceNeighbour(mesh_.faceNeighbour());
1779  const label oldInternalFaces = mesh_.nInternalFaces();
1780 
1781  // Remove cells.
1782  autoPtr<mapPolyMesh> subMap
1783  (
1784  doRemoveCells
1785  (
1786  select(false, distribution, Pstream::myProcNo()),
1787  oldInternalPatchI
1788  )
1789  );
1790 
1791  // Addressing from subsetted mesh
1792  subCellMap[Pstream::myProcNo()] = subMap().cellMap();
1793  subFaceMap[Pstream::myProcNo()] = renumber
1794  (
1795  repatchFaceMap,
1796  subMap().faceMap()
1797  );
1798  subPointMap[Pstream::myProcNo()] = subMap().pointMap();
1799  subPatchMap[Pstream::myProcNo()] = identity(patches.size());
1800 
1801  // Initialize all addressing into current mesh
1802  constructCellMap[Pstream::myProcNo()] = identity(mesh_.nCells());
1803  constructFaceMap[Pstream::myProcNo()] = identity(mesh_.nFaces());
1804  constructPointMap[Pstream::myProcNo()] = identity(mesh_.nPoints());
1805  constructPatchMap[Pstream::myProcNo()] = identity(patches.size());
1806 
1807  // Subset the mesh data: neighbourCell/neighbourProc
1808  // fields
1809  labelList domainSourceFace;
1810  labelList domainSourceProc;
1811  labelList domainSourcePatch;
1812  labelList domainSourceNewNbrProc;
1813 
1814  subsetBoundaryData
1815  (
1816  mesh_, // new mesh
1817  subMap().faceMap(), // from new to original mesh
1818  subMap().cellMap(),
1819 
1820  distribution, // distribution before subsetting
1821  oldFaceOwner, // owner before subsetting
1822  oldFaceNeighbour, // neighbour ,,
1823  oldInternalFaces, // nInternalFaces ,,
1824 
1825  sourceFace,
1826  sourceProc,
1827  sourcePatch,
1828  sourceNewNbrProc,
1829 
1830  domainSourceFace,
1831  domainSourceProc,
1832  domainSourcePatch,
1833  domainSourceNewNbrProc
1834  );
1835 
1836  sourceFace.transfer(domainSourceFace);
1837  sourceProc.transfer(domainSourceProc);
1838  sourcePatch.transfer(domainSourcePatch);
1839  sourceNewNbrProc.transfer(domainSourceNewNbrProc);
1840  }
1841 
1842 
1843  // Print a bit.
1844  if (debug)
1845  {
1846  Pout<< nl << "STARTING MESH:" << endl;
1847  printMeshInfo(mesh_);
1848  printFieldInfo<volScalarField>(mesh_);
1849  printFieldInfo<volVectorField>(mesh_);
1850  printFieldInfo<volSphericalTensorField>(mesh_);
1851  printFieldInfo<volSymmTensorField>(mesh_);
1852  printFieldInfo<volTensorField>(mesh_);
1853  printFieldInfo<surfaceScalarField>(mesh_);
1854  printFieldInfo<surfaceVectorField>(mesh_);
1855  printFieldInfo<surfaceSphericalTensorField>(mesh_);
1856  printFieldInfo<surfaceSymmTensorField>(mesh_);
1857  printFieldInfo<surfaceTensorField>(mesh_);
1858  Pout<< nl << endl;
1859  }
1860 
1861 
1862 
1863  // Receive and add what was sent
1864  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1865 
1866  oldParRun = UPstream::parRun();
1867  UPstream::parRun() = false;
1868 
1869  forAll(nSendCells, sendProc)
1870  {
1871  // Did processor sendProc send anything to me?
1872  if
1873  (
1874  sendProc != Pstream::myProcNo()
1875  && nSendCells[sendProc][Pstream::myProcNo()] > 0
1876  )
1877  {
1878  if (debug)
1879  {
1880  Pout<< nl
1881  << "RECEIVING FROM DOMAIN " << sendProc
1882  << " cells to receive:"
1883  << nSendCells[sendProc][Pstream::myProcNo()]
1884  << nl << endl;
1885  }
1886 
1887 
1888  // Pstream for receiving mesh and fields
1889  UIPstream str(sendProc, pBufs);
1890 
1891 
1892  // Receive from sendProc
1893  labelList domainSourceFace;
1894  labelList domainSourceProc;
1895  labelList domainSourcePatch;
1896  labelList domainSourceNewNbrProc;
1897 
1898  autoPtr<fvMesh> domainMeshPtr;
1909 
1910  // Opposite of sendMesh
1911  {
1912  domainMeshPtr = receiveMesh
1913  (
1914  sendProc,
1915  pointZoneNames,
1916  faceZoneNames,
1917  cellZoneNames,
1918 
1919  const_cast<Time&>(mesh_.time()),
1920  domainSourceFace,
1921  domainSourceProc,
1922  domainSourcePatch,
1923  domainSourceNewNbrProc,
1924  str
1925  );
1926  fvMesh& domainMesh = domainMeshPtr();
1927  // Force construction of various on mesh.
1928  //(void)domainMesh.globalData();
1929 
1930 
1931  // Receive fields. Read as single dictionary because
1932  // of problems reading consecutive fields from single stream.
1933  dictionary fieldDicts(str);
1934 
1935  receiveFields<volScalarField>
1936  (
1937  sendProc,
1938  volScalars,
1939  domainMesh,
1940  vsf,
1941  fieldDicts.subDict(volScalarField::typeName)
1942  );
1943  receiveFields<volVectorField>
1944  (
1945  sendProc,
1946  volVectors,
1947  domainMesh,
1948  vvf,
1949  fieldDicts.subDict(volVectorField::typeName)
1950  );
1951  receiveFields<volSphericalTensorField>
1952  (
1953  sendProc,
1954  volSphereTensors,
1955  domainMesh,
1956  vsptf,
1958  );
1959  receiveFields<volSymmTensorField>
1960  (
1961  sendProc,
1962  volSymmTensors,
1963  domainMesh,
1964  vsytf,
1966  );
1967  receiveFields<volTensorField>
1968  (
1969  sendProc,
1970  volTensors,
1971  domainMesh,
1972  vtf,
1973  fieldDicts.subDict(volTensorField::typeName)
1974  );
1975 
1976  receiveFields<surfaceScalarField>
1977  (
1978  sendProc,
1979  surfScalars,
1980  domainMesh,
1981  ssf,
1983  );
1984  receiveFields<surfaceVectorField>
1985  (
1986  sendProc,
1987  surfVectors,
1988  domainMesh,
1989  svf,
1991  );
1992  receiveFields<surfaceSphericalTensorField>
1993  (
1994  sendProc,
1995  surfSphereTensors,
1996  domainMesh,
1997  ssptf,
1999  );
2000  receiveFields<surfaceSymmTensorField>
2001  (
2002  sendProc,
2003  surfSymmTensors,
2004  domainMesh,
2005  ssytf,
2007  );
2008  receiveFields<surfaceTensorField>
2009  (
2010  sendProc,
2011  surfTensors,
2012  domainMesh,
2013  stf,
2015  );
2016  }
2017  const fvMesh& domainMesh = domainMeshPtr();
2018 
2019 
2020  constructCellMap[sendProc] = identity(domainMesh.nCells());
2021  constructFaceMap[sendProc] = identity(domainMesh.nFaces());
2022  constructPointMap[sendProc] = identity(domainMesh.nPoints());
2023  constructPatchMap[sendProc] =
2024  identity(domainMesh.boundaryMesh().size());
2025 
2026 
2027  // Print a bit.
2028  if (debug)
2029  {
2030  Pout<< nl << "RECEIVED MESH FROM:" << sendProc << endl;
2031  printMeshInfo(domainMesh);
2032  printFieldInfo<volScalarField>(domainMesh);
2033  printFieldInfo<volVectorField>(domainMesh);
2034  printFieldInfo<volSphericalTensorField>(domainMesh);
2035  printFieldInfo<volSymmTensorField>(domainMesh);
2036  printFieldInfo<volTensorField>(domainMesh);
2037  printFieldInfo<surfaceScalarField>(domainMesh);
2038  printFieldInfo<surfaceVectorField>(domainMesh);
2039  printFieldInfo<surfaceSphericalTensorField>(domainMesh);
2040  printFieldInfo<surfaceSymmTensorField>(domainMesh);
2041  printFieldInfo<surfaceTensorField>(domainMesh);
2042  }
2043 
2044 
2045  // Now this mesh we received (from sendProc) needs to be merged
2046  // with the current mesh. On the current mesh we have for all
2047  // boundaryfaces the original face and processor. See if we can
2048  // match these up to the received domainSourceFace and
2049  // domainSourceProc.
2050  labelList masterCoupledFaces;
2051  labelList slaveCoupledFaces;
2052  findCouples
2053  (
2054  mesh_,
2055 
2056  sourceFace,
2057  sourceProc,
2058  sourcePatch,
2059 
2060  sendProc,
2061  domainMesh,
2062  domainSourceFace,
2063  domainSourceProc,
2064  domainSourcePatch,
2065 
2066  masterCoupledFaces,
2067  slaveCoupledFaces
2068  );
2069 
2070  // Generate additional coupling info (points, edges) from
2071  // faces-that-match
2072  faceCoupleInfo couples
2073  (
2074  mesh_,
2075  masterCoupledFaces,
2076  domainMesh,
2077  slaveCoupledFaces,
2078  mergeTol_, // merge tolerance
2079  true, // faces align
2080  true, // couples are ordered already
2081  false
2082  );
2083 
2084 
2085  // Add domainMesh to mesh
2086  // ~~~~~~~~~~~~~~~~~~~~~~
2087 
2089  (
2090  mesh_,
2091  domainMesh,
2092  couples,
2093  false // no parallel comms
2094  );
2095 
2096  // Update mesh data: sourceFace,sourceProc for added
2097  // mesh.
2098 
2099  sourceFace = mapBoundaryData
2100  (
2101  mesh_,
2102  map(),
2103  sourceFace,
2104  domainMesh.nInternalFaces(),
2105  domainSourceFace
2106  );
2107  sourceProc = mapBoundaryData
2108  (
2109  mesh_,
2110  map(),
2111  sourceProc,
2112  domainMesh.nInternalFaces(),
2113  domainSourceProc
2114  );
2115  sourcePatch = mapBoundaryData
2116  (
2117  mesh_,
2118  map(),
2119  sourcePatch,
2120  domainMesh.nInternalFaces(),
2121  domainSourcePatch
2122  );
2123  sourceNewNbrProc = mapBoundaryData
2124  (
2125  mesh_,
2126  map(),
2127  sourceNewNbrProc,
2128  domainMesh.nInternalFaces(),
2129  domainSourceNewNbrProc
2130  );
2131 
2132  // Update all addressing so xxProcAddressing points to correct item
2133  // in masterMesh.
2134  const labelList& oldCellMap = map().oldCellMap();
2135  const labelList& oldFaceMap = map().oldFaceMap();
2136  const labelList& oldPointMap = map().oldPointMap();
2137  const labelList& oldPatchMap = map().oldPatchMap();
2138 
2139  forAll(constructPatchMap, procI)
2140  {
2141  if (procI != sendProc && constructPatchMap[procI].size())
2142  {
2143  // Processor already in mesh (either myProcNo or received)
2144  inplaceRenumber(oldCellMap, constructCellMap[procI]);
2145  inplaceRenumber(oldFaceMap, constructFaceMap[procI]);
2146  inplaceRenumber(oldPointMap, constructPointMap[procI]);
2147  inplaceRenumber(oldPatchMap, constructPatchMap[procI]);
2148  }
2149  }
2150 
2151  // Added processor
2152  inplaceRenumber(map().addedCellMap(), constructCellMap[sendProc]);
2153  inplaceRenumber(map().addedFaceMap(), constructFaceMap[sendProc]);
2154  inplaceRenumber(map().addedPointMap(), constructPointMap[sendProc]);
2155  inplaceRenumber(map().addedPatchMap(), constructPatchMap[sendProc]);
2156 
2157  if (debug)
2158  {
2159  Pout<< nl << "MERGED MESH FROM:" << sendProc << endl;
2160  printMeshInfo(mesh_);
2161  printFieldInfo<volScalarField>(mesh_);
2162  printFieldInfo<volVectorField>(mesh_);
2163  printFieldInfo<volSphericalTensorField>(mesh_);
2164  printFieldInfo<volSymmTensorField>(mesh_);
2165  printFieldInfo<volTensorField>(mesh_);
2166  printFieldInfo<surfaceScalarField>(mesh_);
2167  printFieldInfo<surfaceVectorField>(mesh_);
2168  printFieldInfo<surfaceSphericalTensorField>(mesh_);
2169  printFieldInfo<surfaceSymmTensorField>(mesh_);
2170  printFieldInfo<surfaceTensorField>(mesh_);
2171  Pout<< nl << endl;
2172  }
2173  }
2174  }
2175 
2176  UPstream::parRun() = oldParRun;
2177 
2178  // Print a bit.
2179  if (debug)
2180  {
2181  Pout<< nl << "REDISTRIBUTED MESH:" << endl;
2182  printMeshInfo(mesh_);
2183  printFieldInfo<volScalarField>(mesh_);
2184  printFieldInfo<volVectorField>(mesh_);
2185  printFieldInfo<volSphericalTensorField>(mesh_);
2186  printFieldInfo<volSymmTensorField>(mesh_);
2187  printFieldInfo<volTensorField>(mesh_);
2188  printFieldInfo<surfaceScalarField>(mesh_);
2189  printFieldInfo<surfaceVectorField>(mesh_);
2190  printFieldInfo<surfaceSphericalTensorField>(mesh_);
2191  printFieldInfo<surfaceSymmTensorField>(mesh_);
2192  printFieldInfo<surfaceTensorField>(mesh_);
2193  Pout<< nl << endl;
2194  }
2195 
2196 
2197 
2198  // Add processorPatches
2199  // ~~~~~~~~~~~~~~~~~~~~
2200 
2201  // Per neighbour processor, per originating patch, the patchID
2202  // For faces resulting from internal faces or normal processor patches
2203  // the originating patch is -1. For cyclics this is the cyclic patchID.
2204  List<Map<label> > procPatchID;
2205 
2206  // Add processor and processorCyclic patches.
2207  addProcPatches(sourceNewNbrProc, sourcePatch, procPatchID);
2208 
2209  // Put faces into correct patch. Note that we now have proper
2210  // processorPolyPatches again so repatching will take care of coupled face
2211  // ordering.
2212 
2213  // Get boundary faces to be repatched. Is -1 or new patchID
2214  labelList newPatchID
2215  (
2216  getBoundaryPatch
2217  (
2218  sourceNewNbrProc,
2219  sourcePatch,
2220  procPatchID
2221  )
2222  );
2223 
2224  // Change patches. Since this might change ordering of coupled faces
2225  // we also need to adapt our constructMaps.
2226  // NOTE: there is one very particular problem with this structure.
2227  // We first create the processor patches and use these to merge out
2228  // shared points (see mergeSharedPoints below). So temporarily points
2229  // and edges do not match!
2230  repatch(newPatchID, constructFaceMap);
2231 
2232  // See if any geometrically shared points need to be merged. Note: does
2233  // parallel comms. After this points and edges should again be consistent.
2234  mergeSharedPoints(constructPointMap);
2235 
2236  // Bit of hack: processorFvPatchField does not get reset since created
2237  // from nothing so explicitly reset.
2238  initPatchFields<volScalarField, processorFvPatchField<scalar> >
2239  (
2241  );
2242  initPatchFields<volVectorField, processorFvPatchField<vector> >
2243  (
2245  );
2246  initPatchFields
2247  <
2250  >
2251  (
2253  );
2254  initPatchFields<volSymmTensorField, processorFvPatchField<symmTensor> >
2255  (
2257  );
2258  initPatchFields<volTensorField, processorFvPatchField<tensor> >
2259  (
2261  );
2262 
2263  initPatchFields<surfaceScalarField, processorFvsPatchField<scalar> >
2264  (
2266  );
2267  initPatchFields<surfaceVectorField, processorFvsPatchField<vector> >
2268  (
2270  );
2271  initPatchFields
2272  <
2275  >
2276  (
2278  );
2279  initPatchFields
2280  <
2283  >
2284  (
2286  );
2287  initPatchFields<surfaceTensorField, processorFvsPatchField<tensor> >
2288  (
2290  );
2291 
2292 
2293  mesh_.setInstance(mesh_.time().timeName());
2294 
2295 
2296  // Print a bit
2297  if (debug)
2298  {
2299  Pout<< nl << "FINAL MESH:" << endl;
2300  printMeshInfo(mesh_);
2301  printFieldInfo<volScalarField>(mesh_);
2302  printFieldInfo<volVectorField>(mesh_);
2303  printFieldInfo<volSphericalTensorField>(mesh_);
2304  printFieldInfo<volSymmTensorField>(mesh_);
2305  printFieldInfo<volTensorField>(mesh_);
2306  printFieldInfo<surfaceScalarField>(mesh_);
2307  printFieldInfo<surfaceVectorField>(mesh_);
2308  printFieldInfo<surfaceSphericalTensorField>(mesh_);
2309  printFieldInfo<surfaceSymmTensorField>(mesh_);
2310  printFieldInfo<surfaceTensorField>(mesh_);
2311  Pout<< nl << endl;
2312  }
2313 
2314  // Collect all maps and return
2316  (
2318  (
2319  mesh_,
2320 
2321  nOldPoints,
2322  nOldFaces,
2323  nOldCells,
2324  oldPatchStarts.xfer(),
2325  oldPatchNMeshPoints.xfer(),
2326 
2327  subPointMap.xfer(),
2328  subFaceMap.xfer(),
2329  subCellMap.xfer(),
2330  subPatchMap.xfer(),
2331 
2332  constructPointMap.xfer(),
2333  constructFaceMap.xfer(),
2334  constructCellMap.xfer(),
2335  constructPatchMap.xfer()
2336  )
2337  );
2338 }
2339 
2340 
2341 // ************************************************************************* //
labelList getExposedFaces(const labelList &cellsToRemove) const
Get labels of exposed faces.
Definition: removeCells.C:73
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:469
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:380
void resetMotion() const
Reset motion.
Definition: polyMesh.C:1192
void addFvPatches(const List< polyPatch * > &, const bool validBoundary=true)
Add boundary patches. Constructor helper.
Definition: fvMesh.C:462
const boolList & flipMap() const
Return face flip map.
Definition: faceZone.H:249
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
static void printMeshInfo(const fvMesh &)
Print some info on mesh.
const faceZoneMesh & faceZones() const
Return face zone mesh.
Definition: polyMesh.H:463
static label addPatch(fvMesh &mesh, const polyPatch &patch, const dictionary &patchFieldDict, const word &defaultPatchFieldType, const bool validBoundary)
Add patch. Inserts patch before all processor patches.
Definition: fvMeshTools.C:32
Buffers for inter-processor communications streams (UOPstream, UIPstream).
static const char *const typeName
Definition: Field.H:94
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::surfaceFields.
label size() const
Return the number of elements in the PtrList.
Definition: PtrListI.H:32
label nFaces() const
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
Definition: polyPatch.H:310
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
void setInstance(const fileName &)
Set the instance for mesh files.
Definition: polyMeshIO.C:32
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
This boundary condition enables processor communication across patches.
const word & name() const
Return name.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:390
A subset of mesh cells.
Definition: cellZone.H:61
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:54
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:53
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
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:758
virtual bool owner() const
Does the processor own the patch ?
Xfer< List< T > > xfer()
Transfer contents to the Xfer container.
Definition: ListI.H:90
A subset of mesh points. The labels of points in the zone can be obtained from the addressing() list...
Definition: pointZone.H:62
An STL-conforming hash table.
Definition: HashTable.H:61
bool checkParallelSync(const bool report=false) const
Check whether all procs have all patches and in same order. Return.
void inplaceReorder(const labelUList &oldToNew, ListType &)
Inplace reorder the elements of a list.
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:64
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:316
static void mergePoints(const polyMesh &, const Map< label > &pointToMaster, polyTopoChange &meshMod)
Helper: Merge points.
Foam::fvBoundaryMesh.
A class for handling words, derived from string.
Definition: word.H:59
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
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
static void reorderPatches(fvMesh &, const labelList &oldToNew, const label nPatches, const bool validBoundary)
Definition: fvMeshTools.C:324
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
static void printCoupleInfo(const primitiveMesh &, const labelList &, const labelList &, const labelList &, const labelList &)
Print some info on coupling data.
Given list of cells to remove insert all the topology changes.
Definition: removeCells.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void setRefinement(const labelList &cellsToRemove, const labelList &facesToExpose, const labelList &patchIDs, polyTopoChange &) const
Play commands into polyTopoChange to remove cells.
Definition: removeCells.C:179
label nCells() const
const pointZoneMesh & pointZones() const
Return point zone mesh.
Definition: polyMesh.H:457
const fvMesh & subMesh() const
Return reference to subset mesh.
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:278
const labelList & cellMap() const
Return cell map.
Container for information needed to couple to meshes. When constructed from two meshes and a geometri...
Xfer< T > xferMove(T &)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:638
This boundary condition enables processor communication across cyclic patches.
const labelList & addedFaceMap() const
Namespace for OpenFOAM.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
Definition: faceZone.C:317
patches[0]
void setSize(const label nRows)
Reset size of CompactListList.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
wordList names() const
Return a list of patch names.
const labelList & oldFaceMap() const
dictionary dict
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:74
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:386
static const char nl
Definition: Ostream.H:260
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Neighbour processor patch.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1035
A List with indirect addressing.
Definition: fvMatrix.H:106
Foam::polyBoundaryMesh.
static Map< label > findSharedPoints(const polyMesh &, const scalar mergeTol)
Find topologically and geometrically shared points.
#define forAll(list, i)
Definition: UList.H:421
ListType renumber(const labelUList &oldToNew, const ListType &)
Renumber the values (not the indices) of a list.
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &l)
Swap coupled boundary face values.
Definition: syncTools.H:429
virtual bool owner() const
Does this side own the patch ?
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
Class describing modification of a face.
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
Cyclic plane patch.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const word & name() const
Return name.
Definition: IOobject.H:260
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1073
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:404
const labelList & pointMap() const
Return point map.
const labelList & patchMap() const
Return patch map.
const word & name() const
Return name.
Definition: zone.H:150
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
GeometricField< sphericalTensor, fvsPatchField, surfaceMesh > surfaceSphericalTensorField
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
A packed storage unstructured matrix of objects of type <T> using an offset table for access...
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:55
wordList names() const
Return the list of names of the IOobjects.
bool found
static labelList countCells(const labelList &)
Helper function: count cells per processor in wanted distribution.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
autoPtr< mapDistributePolyMesh > distribute(const labelList &dist)
Send cells to neighbours according to distribution.
label nInternalFaces() const
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:306
#define forAllReverse(list, i)
Definition: UList.H:424
Traits class for primitives.
Definition: pTraits.H:50
List< word > wordList
A List of words.
Definition: fileName.H:54
virtual transformType transform() const
Type of transform.
error FatalError
Xfer< HashTable< T, Key, Hash > > xfer()
Transfer contents to the Xfer container.
Definition: HashTableI.H:102
List< label > labelList
A List of labels.
Definition: labelList.H:56
void clearOut()
Clear all geometry and addressing.
Definition: fvMesh.C:231
void addZones(const List< pointZone * > &pz, const List< faceZone * > &fz, const List< cellZone * > &cz)
Add mesh zones.
Definition: polyMesh.C:971
Pair< label > labelPair
Label pair.
Definition: labelPair.H:48
int neighbProcNo() const
Return neigbour processor number.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
label findZoneID(const word &zoneName) const
Find zone index given a name.
Definition: ZoneMesh.C:354
Direct mesh changes based on v1.3 polyTopoChange syntax.
A List obtained as a section of another List.
Definition: SubList.H:53
GeometricField< sphericalTensor, fvPatchField, volMesh > volSphericalTensorField
Definition: volFieldsFwd.H:57
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1060
wordList names() const
Return a list of zone names.
Definition: ZoneMesh.C:269
void finishedSends(const bool block=true)
Mark all sends as having been done. This will start receives.
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
const labelList & faceMap() const
Return face map.
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1079
Foam::processorFvsPatchField.
label nOldInternalFaces() const
Number of old internal faces.
GeometricField< symmTensor, fvsPatchField, surfaceMesh > surfaceSymmTensorField
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
static autoPtr< mapAddedPolyMesh > add(fvMesh &mesh0, const fvMesh &mesh1, const faceCoupleInfo &coupleInfo, const bool validBoundary=true)
Inplace add mesh to fvMesh. Maps all stored fields. Returns map.
Definition: fvMeshAdder.C:66
void inplaceRenumber(const labelUList &oldToNew, ListType &)
Inplace renumber the values of a list.
autoPtr< mapPolyMesh > changeMesh(polyMesh &mesh, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
label nPoints() const
Combination-Reduction operation for a parallel run. The information from all nodes is collected on th...
defineTypeNameAndDebug(combustionModel, 0)
Post-processing mesh subset tool. Given the original mesh and the list of selected cells...
Definition: fvMeshSubset.H:71
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
const vectorField & faceCentres() const
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:552
Class containing mesh-to-mesh mapping information after a mesh addition where we add a mesh (&#39;added m...