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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2018 OpenFOAM Foundation
9  Copyright (C) 2015-2023 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "fvMeshDistribute.H"
30 #include "fvMeshAdder.H"
31 #include "fvMeshSubset.H"
32 #include "faceCoupleInfo.H"
33 #include "processorFvPatchField.H"
34 #include "processorFvsPatchField.H"
35 #include "processorCyclicPolyPatch.H"
36 #include "processorCyclicFvPatchField.H"
37 #include "polyTopoChange.H"
38 #include "removeCells.H"
39 #include "polyModifyFace.H"
40 #include "polyRemovePoint.H"
41 #include "mapDistributePolyMesh.H"
42 #include "surfaceFields.H"
43 #include "syncTools.H"
44 #include "CompactListList.H"
45 #include "fvMeshTools.H"
46 #include "labelPairHashes.H"
47 #include "ListOps.H"
48 #include "globalIndex.H"
49 #include "cyclicACMIPolyPatch.H"
50 #include "mappedPatchBase.H"
51 
52 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56  defineTypeNameAndDebug(fvMeshDistribute, 0);
57 
58  //- Less function class that can be used for sorting processor patches
59  class lessProcPatches
60  {
61  const labelList& nbrProc_;
62  const labelList& referPatchID_;
63 
64  public:
65 
66  lessProcPatches(const labelList& nbrProc, const labelList& referPatchID)
67  :
68  nbrProc_(nbrProc),
69  referPatchID_(referPatchID)
70  {}
71 
72  bool operator()(const label a, const label b)
73  {
74  if (nbrProc_[a] < nbrProc_[b])
75  {
76  return true;
77  }
78  else if (nbrProc_[a] > nbrProc_[b])
79  {
80  return false;
81  }
82  else
83  {
84  // Equal neighbour processor
85  return referPatchID_[a] < referPatchID_[b];
86  }
87  }
88  };
89 }
90 
91 
92 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
93 
94 void Foam::fvMeshDistribute::inplaceRenumberWithFlip
95 (
96  const labelUList& oldToNew,
97  const bool oldToNewHasFlip,
98  const bool lstHasFlip,
99  labelUList& lst
100 )
101 {
102  if (!lstHasFlip && !oldToNewHasFlip)
103  {
104  Foam::inplaceRenumber(oldToNew, lst);
105  }
106  else
107  {
108  // Either input data or map encodes sign so result encodes sign
109 
110  forAll(lst, elemI)
111  {
112  // Extract old value and sign
113  label val = lst[elemI];
114  label sign = 1;
115  if (lstHasFlip)
116  {
117  if (val > 0)
118  {
119  val = val-1;
120  }
121  else if (val < 0)
122  {
123  val = -val-1;
124  sign = -1;
125  }
126  else
127  {
128  FatalErrorInFunction
129  << "Problem : zero value " << val
130  << " at index " << elemI << " out of " << lst.size()
131  << " list with flip bit" << exit(FatalError);
132  }
133  }
134 
135 
136  // Lookup new value and possibly change sign
137  label newVal = oldToNew[val];
138 
139  if (oldToNewHasFlip)
140  {
141  if (newVal > 0)
142  {
143  newVal = newVal-1;
144  }
145  else if (newVal < 0)
146  {
147  newVal = -newVal-1;
148  sign = -sign;
149  }
150  else
151  {
152  FatalErrorInFunction
153  << "Problem : zero value " << newVal
154  << " at index " << elemI << " out of "
155  << oldToNew.size()
156  << " list with flip bit" << exit(FatalError);
157  }
158  }
159 
160 
161  // Encode new value and sign
162  lst[elemI] = sign*(newVal+1);
163  }
164  }
165 }
166 
167 
168 Foam::labelList Foam::fvMeshDistribute::select
169 (
170  const bool selectEqual,
171  const labelList& values,
172  const label value
173 )
174 {
175  label n = 0;
176 
177  forAll(values, i)
178  {
179  if (selectEqual == (values[i] == value))
180  {
181  n++;
182  }
183  }
184 
185  labelList indices(n);
186  n = 0;
187 
188  forAll(values, i)
189  {
190  if (selectEqual == (values[i] == value))
191  {
192  indices[n++] = i;
193  }
194  }
195  return indices;
196 }
197 
198 
199 Foam::wordList Foam::fvMeshDistribute::mergeWordList(const wordList& procNames)
200 {
201  List<wordList> allNames(Pstream::nProcs());
202  allNames[Pstream::myProcNo()] = procNames;
203  Pstream::gatherList(allNames);
204 
205  // Assume there are few zone names to merge. Use HashSet otherwise (but
206  // maintain ordering)
207  DynamicList<word> mergedNames;
208  if (Pstream::master())
209  {
210  mergedNames = procNames;
211  for (const wordList& names : allNames)
212  {
213  for (const word& name : names)
214  {
215  mergedNames.push_uniq(name);
216  }
217  }
218  }
219  Pstream::broadcast(mergedNames);
220 
221  return mergedNames;
222 }
223 
224 
225 void Foam::fvMeshDistribute::printMeshInfo(const polyMesh& mesh)
226 {
227  Pout<< "Primitives:" << nl
228  << " points :" << mesh.nPoints() << nl
229  << " bb :" << boundBox(mesh.points(), false) << nl
230  << " internalFaces:" << mesh.nInternalFaces() << nl
231  << " faces :" << mesh.nFaces() << nl
232  << " cells :" << mesh.nCells() << nl;
233 
234  Pout<< "Patches:" << endl;
235  for (const polyPatch& pp : mesh.boundaryMesh())
236  {
237  Pout<< " " << pp.index() << " name:" << pp.name()
238  << " size:" << pp.size()
239  << " start:" << pp.start()
240  << " type:" << pp.type()
241  << endl;
242  }
243 
244  if (mesh.pointZones().empty()) Pout<< "PointZones:" << endl;
245  for (const auto& zn : mesh.pointZones())
246  {
247  Pout<< " " << zn.index() << " name:" << zn.name()
248  << " size:" << zn.size() << endl;
249  }
250 
251  if (!mesh.faceZones().empty()) Pout<< "FaceZones:" << endl;
252  for (const auto& zn : mesh.faceZones())
253  {
254  Pout<< " " << zn.index() << " name:" << zn.name()
255  << " size:" << zn.size() << endl;
256  }
257 
258  if (!mesh.cellZones().empty()) Pout<< "CellZones:" << endl;
259  for (const auto& zn : mesh.cellZones())
260  {
261  Pout<< " " << zn.index() << " name:" << zn.name()
262  << " size:" << zn.size() << endl;
263  }
264 }
265 
266 
267 void Foam::fvMeshDistribute::printCoupleInfo
268 (
269  const primitiveMesh& mesh,
270  const labelList& sourceFace,
271  const labelList& sourceProc,
272  const labelList& sourcePatch,
273  const labelList& sourceNewNbrProc
274 )
275 {
276  Pout<< nl
277  << "Current coupling info:"
278  << endl;
279 
280  forAll(sourceFace, bFacei)
281  {
282  label meshFacei = mesh.nInternalFaces() + bFacei;
283 
284  Pout<< " meshFace:" << meshFacei
285  << " fc:" << mesh.faceCentres()[meshFacei]
286  << " connects to proc:" << sourceProc[bFacei]
287  << "/face:" << sourceFace[bFacei]
288  << " which will move to proc:" << sourceNewNbrProc[bFacei]
289  << endl;
290  }
291 }
292 
293 
294 Foam::label Foam::fvMeshDistribute::findNonEmptyPatch() const
295 {
296  // Finds (non-empty) patch that exposed internal and proc faces can be
297  // put into.
298  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
299 
300 
301  // Mark 'special' patches : -coupled, -duplicate faces. These probably
302  // should not be used to (temporarily) store processor faces ...
303 
304  bitSet isCoupledPatch(patches.size());
305  forAll(patches, patchi)
306  {
307  const polyPatch& pp = patches[patchi];
308  const auto* cpp = isA<cyclicACMIPolyPatch>(pp);
309 
310  if (cpp)
311  {
312  isCoupledPatch.set(patchi);
313  const label dupPatchID = cpp->nonOverlapPatchID();
314  if (dupPatchID != -1)
315  {
316  isCoupledPatch.set(dupPatchID);
317  }
318  }
319  else if (pp.coupled())
320  {
321  isCoupledPatch.set(patchi);
322  }
323  }
324 
325  label nonEmptyPatchi = -1;
326 
327  forAllReverse(patches, patchi)
328  {
329  const polyPatch& pp = patches[patchi];
330 
331  if
332  (
333  !isA<emptyPolyPatch>(pp)
334  && !isCoupledPatch(patchi)
335  && !isA<mappedPatchBase>(pp)
336  )
337  {
338  nonEmptyPatchi = patchi;
339  break;
340  }
341  }
342 
343  if (nonEmptyPatchi == -1)
344  {
345  FatalErrorInFunction
346  << "Cannot find a patch which is not of type empty, mapped or"
347  << " coupled in patches " << patches.names() << endl
348  << "There has to be at least one such patch for"
349  << " distribution to work" << abort(FatalError);
350  }
351 
352  if (debug)
353  {
354  Pout<< "findNonEmptyPatch : using patch " << nonEmptyPatchi
355  << " name:" << patches[nonEmptyPatchi].name()
356  << " type:" << patches[nonEmptyPatchi].type()
357  << " to put exposed faces into." << endl;
358  }
359 
360 
361  // Do additional test for processor patches intermingled with non-proc
362  // patches.
363  label procPatchi = -1;
364 
365  forAll(patches, patchi)
366  {
367  if (isA<processorPolyPatch>(patches[patchi]))
368  {
369  procPatchi = patchi;
370  }
371  else if (procPatchi != -1)
372  {
373  FatalErrorInFunction
374  << "Processor patches should be at end of patch list."
375  << endl
376  << "Have processor patch " << procPatchi
377  << " followed by non-processor patch " << patchi
378  << " in patches " << patches.names()
379  << abort(FatalError);
380  }
381  }
382 
383  return nonEmptyPatchi;
384 }
385 
386 
387 Foam::tmp<Foam::surfaceScalarField> Foam::fvMeshDistribute::generateTestField
388 (
389  const fvMesh& mesh
390 )
391 {
392  const vector testNormal = normalised(vector::one);
393 
394  tmp<surfaceScalarField> tfld
395  (
396  new surfaceScalarField
397  (
398  IOobject
399  (
400  "myFlux",
401  mesh.time().timeName(),
402  mesh,
403  IOobject::NO_READ,
404  IOobject::NO_WRITE
405  ),
406  mesh,
407  dimensionedScalar(dimless, Zero)
408  )
409  );
410  surfaceScalarField& fld = tfld.ref();
411 
412  const surfaceVectorField n(mesh.Sf()/mesh.magSf());
413 
414  forAll(fld, facei)
415  {
416  fld[facei] = (n[facei] & testNormal);
417  }
418 
419  surfaceScalarField::Boundary& fluxBf = fld.boundaryFieldRef();
420  const surfaceVectorField::Boundary& nBf = n.boundaryField();
421 
422  forAll(fluxBf, patchi)
423  {
424  fvsPatchScalarField& fvp = fluxBf[patchi];
425 
426  scalarField newPfld(fvp.size());
427  forAll(newPfld, i)
428  {
429  newPfld[i] = (nBf[patchi][i] & testNormal);
430  }
431  fvp == newPfld;
432  }
433 
434  return tfld;
435 }
436 
437 
438 void Foam::fvMeshDistribute::testField(const surfaceScalarField& fld)
439 {
440  const fvMesh& mesh = fld.mesh();
441 
442  const vector testNormal = normalised(vector::one);
443 
444  const surfaceVectorField n(mesh.Sf()/mesh.magSf());
445 
446  forAll(fld, facei)
447  {
448  scalar cos = (n[facei] & testNormal);
449 
450  if (mag(cos - fld[facei]) > 1e-6)
451  {
452  //FatalErrorInFunction
453  WarningInFunction
454  << "On internal face " << facei << " at "
455  << mesh.faceCentres()[facei]
456  << " the field value is " << fld[facei]
457  << " whereas cos angle of " << testNormal
458  << " with mesh normal " << n[facei]
459  << " is " << cos
460  //<< exit(FatalError);
461  << endl;
462  }
463  }
464  forAll(fld.boundaryField(), patchi)
465  {
466  const fvsPatchScalarField& fvp = fld.boundaryField()[patchi];
467  const fvsPatchVectorField& np = n.boundaryField()[patchi];
468 
469  forAll(fvp, i)
470  {
471  scalar cos = (np[i] & testNormal);
472 
473  if (mag(cos - fvp[i]) > 1e-6)
474  {
475  label facei = fvp.patch().start()+i;
476  //FatalErrorInFunction
477  WarningInFunction
478  << "On face " << facei
479  << " on patch " << fvp.patch().name()
480  << " at " << mesh.faceCentres()[facei]
481  << " the field value is " << fvp[i]
482  << " whereas cos angle of " << testNormal
483  << " with mesh normal " << np[i]
484  << " is " << cos
485  //<< exit(FatalError);
486  << endl;
487  }
488  }
489  }
490 }
491 
492 
493 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::deleteProcPatches
494 (
495  const label destinationPatch
496 )
497 {
498  // Delete all processor patches. Move any processor faces into the last
499  // non-processor patch.
500 
501  // New patchID per boundary faces to be repatched. Is -1 (no change)
502  // or new patchID
503  labelList newPatchID(mesh_.nBoundaryFaces(), -1);
504 
505  for (const polyPatch& pp : mesh_.boundaryMesh())
506  {
507  if (isA<processorPolyPatch>(pp))
508  {
509  if (debug)
510  {
511  Pout<< "Moving all faces of patch " << pp.name()
512  << " into patch " << destinationPatch
513  << endl;
514  }
515 
516  SubList<label>
517  (
518  newPatchID,
519  pp.size(),
520  pp.offset()
521  ) = destinationPatch;
522  }
523  }
524 
525  // Note: order of boundary faces been kept the same since the
526  // destinationPatch is at the end and we have visited the patches in
527  // incremental order.
528  labelListList dummyFaceMaps;
529  autoPtr<mapPolyMesh> map = repatch(newPatchID, dummyFaceMaps);
530 
531 
532  // Delete (now empty) processor patches.
533  {
534  labelList oldToNew(identity(mesh_.boundaryMesh().size()));
535  label newi = 0;
536  // Non processor patches first
537  forAll(mesh_.boundaryMesh(), patchi)
538  {
539  if (!isA<processorPolyPatch>(mesh_.boundaryMesh()[patchi]))
540  {
541  oldToNew[patchi] = newi++;
542  }
543  }
544  label nNonProcPatches = newi;
545 
546  // Processor patches as last
547  forAll(mesh_.boundaryMesh(), patchi)
548  {
549  if (isA<processorPolyPatch>(mesh_.boundaryMesh()[patchi]))
550  {
551  oldToNew[patchi] = newi++;
552  }
553  }
554  fvMeshTools::reorderPatches(mesh_, oldToNew, nNonProcPatches, false);
555  }
556 
557  return map;
558 }
559 
560 
561 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::repatch
562 (
563  const labelList& newPatchID, // per boundary face -1 or new patchID
564  labelListList& constructFaceMap
565 )
566 {
567  polyTopoChange meshMod(mesh_);
568 
569  forAll(newPatchID, bFacei)
570  {
571  if (newPatchID[bFacei] != -1)
572  {
573  label facei = mesh_.nInternalFaces() + bFacei;
574 
575  label zoneID = mesh_.faceZones().whichZone(facei);
576  bool zoneFlip = false;
577 
578  if (zoneID >= 0)
579  {
580  const faceZone& fZone = mesh_.faceZones()[zoneID];
581  zoneFlip = fZone.flipMap()[fZone.whichFace(facei)];
582  }
583 
584  meshMod.setAction
585  (
586  polyModifyFace
587  (
588  mesh_.faces()[facei], // modified face
589  facei, // label of face
590  mesh_.faceOwner()[facei], // owner
591  -1, // neighbour
592  false, // face flip
593  newPatchID[bFacei], // patch for face
594  false, // remove from zone
595  zoneID, // zone for face
596  zoneFlip // face flip in zone
597  )
598  );
599  }
600  }
601 
602 
603  // Do mapping of fields from one patchField to the other ourselves since
604  // is currently not supported by updateMesh.
605 
606  // Store boundary fields (we only do this for surfaceFields)
607  PtrList<FieldField<fvsPatchField, scalar>> sFlds;
608  saveBoundaryFields<scalar, surfaceMesh>(sFlds);
609  PtrList<FieldField<fvsPatchField, vector>> vFlds;
610  saveBoundaryFields<vector, surfaceMesh>(vFlds);
611  PtrList<FieldField<fvsPatchField, sphericalTensor>> sptFlds;
612  saveBoundaryFields<sphericalTensor, surfaceMesh>(sptFlds);
613  PtrList<FieldField<fvsPatchField, symmTensor>> sytFlds;
614  saveBoundaryFields<symmTensor, surfaceMesh>(sytFlds);
615  PtrList<FieldField<fvsPatchField, tensor>> tFlds;
616  saveBoundaryFields<tensor, surfaceMesh>(tFlds);
617 
618  // Change the mesh (no inflation). Note: parallel comms allowed.
619  //
620  // NOTE: there is one very particular problem with this ordering.
621  // We first create the processor patches and use these to merge out
622  // shared points (see mergeSharedPoints below). So temporarily points
623  // and edges do not match!
624 
625  // TBD: temporarily unset mesh moving to avoid problems in meshflux
626  // mapping. To be investigated.
627  const bool oldMoving = mesh_.moving(false);
628  autoPtr<mapPolyMesh> mapPtr = meshMod.changeMesh(mesh_, false, true);
629  mesh_.moving(oldMoving);
630  mapPolyMesh& map = *mapPtr;
631 
632 
633  // Update fields. No inflation, parallel sync.
634  mesh_.updateMesh(map);
635 
636  // Map patch fields using stored boundary fields. Note: assumes order
637  // of fields has not changed in object registry!
638  mapBoundaryFields<scalar, surfaceMesh>(map, sFlds);
639  mapBoundaryFields<vector, surfaceMesh>(map, vFlds);
640  mapBoundaryFields<sphericalTensor, surfaceMesh>(map, sptFlds);
641  mapBoundaryFields<symmTensor, surfaceMesh>(map, sytFlds);
642  mapBoundaryFields<tensor, surfaceMesh>(map, tFlds);
643 
644 
645  // Move mesh (since morphing does not do this)
646  if (map.hasMotionPoints())
647  {
648  mesh_.movePoints(map.preMotionPoints());
649  }
650 
651  // Adapt constructMaps.
652 
653  if (debug)
654  {
655  label index = map.reverseFaceMap().find(-1);
656 
657  if (index != -1)
658  {
659  FatalErrorInFunction
660  << "reverseFaceMap contains -1 at index:"
661  << index << endl
662  << "This means that the repatch operation was not just"
663  << " a shuffle?" << abort(FatalError);
664  }
665  }
666 
667  forAll(constructFaceMap, proci)
668  {
669  inplaceRenumberWithFlip
670  (
671  map.reverseFaceMap(),
672  false,
673  true,
674  constructFaceMap[proci]
675  );
676  }
677 
678  return mapPtr;
679 }
680 
681 
682 // Detect shared points. Need processor patches to be present.
683 // Background: when adding bits of mesh one can get points which
684 // share the same position but are only detectable to be topologically
685 // the same point when doing parallel analysis. This routine will
686 // merge those points.
687 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::mergeSharedPoints
688 (
689  const labelList& pointToGlobalMaster,
690  labelListList& constructPointMap
691 )
692 {
693  // Find out which sets of points get merged and create a map from
694  // mesh point to unique point.
695 
696  label nShared = 0;
697  forAll(pointToGlobalMaster, pointi)
698  {
699  if (pointToGlobalMaster[pointi] != -1)
700  {
701  nShared++;
702  }
703  }
704 
705  if (debug)
706  {
707  Pout<< "mergeSharedPoints : found " << nShared
708  << " points on processor boundaries" << nl << endl;
709  }
710 
711  Map<label> globalMasterToLocalMaster(2*nShared);
712  Map<label> pointToMaster(2*nShared);
713  label nMatch = 0;
714 
715  forAll(pointToGlobalMaster, pointi)
716  {
717  label globali = pointToGlobalMaster[pointi];
718  if (globali != -1)
719  {
720  const auto iter = globalMasterToLocalMaster.cfind(globali);
721 
722  if (iter.good())
723  {
724  // Matched to existing master
725  pointToMaster.insert(pointi, *iter);
726  nMatch++;
727  }
728  else
729  {
730  // Found first point. Designate as master
731  globalMasterToLocalMaster.insert(globali, pointi);
732  pointToMaster.insert(pointi, pointi);
733  }
734  }
735  }
736 
737  reduce(nMatch, sumOp<label>());
738 
739  if (debug)
740  {
741  Pout<< "mergeSharedPoints : found "
742  << nMatch << " mergeable points" << nl << endl;
743  }
744 
745 
746  if (nMatch == 0)
747  {
748  return nullptr;
749  }
750 
751 
752  polyTopoChange meshMod(mesh_);
753 
754  fvMeshAdder::mergePoints(mesh_, pointToMaster, meshMod);
755 
756  // Change the mesh (no inflation). Note: parallel comms allowed.
757 
758  // TBD: temporarily unset mesh moving to avoid problems in meshflux
759  // mapping. To be investigated.
760  const bool oldMoving = mesh_.moving(false);
761  autoPtr<mapPolyMesh> mapPtr = meshMod.changeMesh(mesh_, false, true);
762  mesh_.moving(oldMoving);
763  mapPolyMesh& map = *mapPtr;
764 
765  // Update fields. No inflation, parallel sync.
766  mesh_.updateMesh(map);
767 
768  // Adapt constructMaps for merged points.
769  forAll(constructPointMap, proci)
770  {
771  labelList& constructMap = constructPointMap[proci];
772 
773  forAll(constructMap, i)
774  {
775  label oldPointi = constructMap[i];
776 
777  label newPointi = map.reversePointMap()[oldPointi];
778 
779  if (newPointi < -1)
780  {
781  constructMap[i] = -newPointi-2;
782  }
783  else if (newPointi >= 0)
784  {
785  constructMap[i] = newPointi;
786  }
787  else
788  {
789  FatalErrorInFunction
790  << "Problem. oldPointi:" << oldPointi
791  << " newPointi:" << newPointi << abort(FatalError);
792  }
793  }
794  }
795 
796  return mapPtr;
797 }
798 
799 
800 void Foam::fvMeshDistribute::getCouplingData
801 (
802  const labelList& distribution,
803  labelList& sourceFace,
804  labelList& sourceProc,
805  labelList& sourcePatch,
806  labelList& sourceNewNbrProc,
807  labelList& sourcePointMaster
808 ) const
809 {
810  // Construct the coupling information for all (boundary) faces and
811  // points
812 
813  const label nBnd = mesh_.nBoundaryFaces();
814  sourceFace.setSize(nBnd);
815  sourceProc.setSize(nBnd);
816  sourcePatch.setSize(nBnd);
817  sourceNewNbrProc.setSize(nBnd);
818 
819  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
820 
821  // Get neighbouring meshFace labels and new processor of coupled boundaries.
822  labelList nbrFaces(nBnd, -1);
823  labelList nbrNewNbrProc(nBnd, -1);
824 
825  forAll(patches, patchi)
826  {
827  const polyPatch& pp = patches[patchi];
828 
829  if (pp.coupled())
830  {
831  label offset = pp.start() - mesh_.nInternalFaces();
832 
833  // Mesh labels of faces on this side
834  forAll(pp, i)
835  {
836  label bndI = offset + i;
837  nbrFaces[bndI] = pp.start()+i;
838  }
839 
840  // Which processor they will end up on
841  SubList<label>(nbrNewNbrProc, pp.size(), offset) =
842  labelUIndList(distribution, pp.faceCells())();
843  }
844  }
845 
846 
847  // Exchange the boundary data
848  syncTools::swapBoundaryFaceList(mesh_, nbrFaces);
849  syncTools::swapBoundaryFaceList(mesh_, nbrNewNbrProc);
850 
851 
852  forAll(patches, patchi)
853  {
854  const polyPatch& pp = patches[patchi];
855  label offset = pp.start() - mesh_.nInternalFaces();
856 
857  if (isA<processorPolyPatch>(pp))
858  {
859  const processorPolyPatch& procPatch =
860  refCast<const processorPolyPatch>(pp);
861 
862  // Check which of the two faces we store.
863 
864  if (procPatch.owner())
865  {
866  // Use my local face labels
867  forAll(pp, i)
868  {
869  label bndI = offset + i;
870  sourceFace[bndI] = pp.start()+i;
871  sourceProc[bndI] = Pstream::myProcNo();
872  sourceNewNbrProc[bndI] = nbrNewNbrProc[bndI];
873  }
874  }
875  else
876  {
877  // Use my neighbours face labels
878  forAll(pp, i)
879  {
880  label bndI = offset + i;
881  sourceFace[bndI] = nbrFaces[bndI];
882  sourceProc[bndI] = procPatch.neighbProcNo();
883  sourceNewNbrProc[bndI] = nbrNewNbrProc[bndI];
884  }
885  }
886 
887 
888  label patchi = -1;
889  if (isA<processorCyclicPolyPatch>(pp))
890  {
891  patchi = refCast<const processorCyclicPolyPatch>
892  (
893  pp
894  ).referPatchID();
895  }
896 
897  forAll(pp, i)
898  {
899  label bndI = offset + i;
900  sourcePatch[bndI] = patchi;
901  }
902  }
903  else if (isA<cyclicPolyPatch>(pp))
904  {
905  const cyclicPolyPatch& cpp = refCast<const cyclicPolyPatch>(pp);
906 
907  if (cpp.owner())
908  {
909  forAll(pp, i)
910  {
911  label bndI = offset + i;
912  sourceFace[bndI] = pp.start()+i;
913  sourceProc[bndI] = Pstream::myProcNo();
914  sourcePatch[bndI] = patchi;
915  sourceNewNbrProc[bndI] = nbrNewNbrProc[bndI];
916  }
917  }
918  else
919  {
920  forAll(pp, i)
921  {
922  label bndI = offset + i;
923  sourceFace[bndI] = nbrFaces[bndI];
924  sourceProc[bndI] = Pstream::myProcNo();
925  sourcePatch[bndI] = patchi;
926  sourceNewNbrProc[bndI] = nbrNewNbrProc[bndI];
927  }
928  }
929  }
930  else
931  {
932  // Normal (physical) boundary
933  forAll(pp, i)
934  {
935  label bndI = offset + i;
936  sourceFace[bndI] = -1;
937  sourceProc[bndI] = -1;
938  sourcePatch[bndI] = patchi;
939  sourceNewNbrProc[bndI] = -1;
940  }
941  }
942  }
943 
944 
945  // Collect coupled (collocated) points
946  sourcePointMaster.resize_nocopy(mesh_.nPoints());
947  sourcePointMaster = -1;
948  {
949  // Assign global master point
950  const globalIndex globalPoints(mesh_.nPoints());
951 
952  const globalMeshData& gmd = mesh_.globalData();
953  const indirectPrimitivePatch& cpp = gmd.coupledPatch();
954  const labelList& meshPoints = cpp.meshPoints();
955  const mapDistribute& slavesMap = gmd.globalCoPointSlavesMap();
956  const labelListList& slaves = gmd.globalCoPointSlaves();
957 
958  labelList elems(slavesMap.constructSize(), -1);
959  forAll(meshPoints, pointi)
960  {
961  const labelList& slots = slaves[pointi];
962 
963  if (slots.size())
964  {
965  // pointi is a master. Assign a unique label.
966 
967  label globalPointi = globalPoints.toGlobal(meshPoints[pointi]);
968  elems[pointi] = globalPointi;
969  forAll(slots, i)
970  {
971  label sloti = slots[i];
972  if (sloti >= meshPoints.size())
973  {
974  // Filter out local collocated points. We don't want
975  // to merge these
976  elems[slots[i]] = globalPointi;
977  }
978  }
979  }
980  }
981 
982  // Push slave-slot data back to slaves
983  slavesMap.reverseDistribute(elems.size(), elems, false);
984 
985  // Extract back onto mesh
986  forAll(meshPoints, pointi)
987  {
988  sourcePointMaster[meshPoints[pointi]] = elems[pointi];
989  }
990  }
991 }
992 
993 
994 // Subset the neighbourCell/neighbourProc fields
995 void Foam::fvMeshDistribute::subsetCouplingData
996 (
997  const fvMesh& mesh,
998  const labelList& pointMap,
999  const labelList& faceMap,
1000  const labelList& cellMap,
1001 
1002  const labelList& oldDistribution,
1003  const labelList& oldFaceOwner,
1004  const labelList& oldFaceNeighbour,
1005  const label oldInternalFaces,
1006 
1007  const labelList& sourceFace,
1008  const labelList& sourceProc,
1009  const labelList& sourcePatch,
1010  const labelList& sourceNewNbrProc,
1011  const labelList& sourcePointMaster,
1012 
1013  labelList& subFace,
1014  labelList& subProc,
1015  labelList& subPatch,
1016  labelList& subNewNbrProc,
1017  labelList& subPointMaster
1018 )
1019 {
1020  subFace.setSize(mesh.nBoundaryFaces());
1021  subProc.setSize(mesh.nBoundaryFaces());
1022  subPatch.setSize(mesh.nBoundaryFaces());
1023  subNewNbrProc.setSize(mesh.nBoundaryFaces());
1024 
1025  forAll(subFace, newBFacei)
1026  {
1027  label newFacei = newBFacei + mesh.nInternalFaces();
1028 
1029  label oldFacei = faceMap[newFacei];
1030 
1031  // Was oldFacei internal face? If so which side did we get.
1032  if (oldFacei < oldInternalFaces)
1033  {
1034  subFace[newBFacei] = oldFacei;
1035  subProc[newBFacei] = Pstream::myProcNo();
1036  subPatch[newBFacei] = -1;
1037 
1038  label oldOwn = oldFaceOwner[oldFacei];
1039  label oldNei = oldFaceNeighbour[oldFacei];
1040 
1041  if (oldOwn == cellMap[mesh.faceOwner()[newFacei]])
1042  {
1043  // We kept the owner side. Where does the neighbour move to?
1044  subNewNbrProc[newBFacei] = oldDistribution[oldNei];
1045  }
1046  else
1047  {
1048  // We kept the neighbour side.
1049  subNewNbrProc[newBFacei] = oldDistribution[oldOwn];
1050  }
1051  }
1052  else
1053  {
1054  // Was boundary face. Take over boundary information
1055  label oldBFacei = oldFacei - oldInternalFaces;
1056 
1057  subFace[newBFacei] = sourceFace[oldBFacei];
1058  subProc[newBFacei] = sourceProc[oldBFacei];
1059  subPatch[newBFacei] = sourcePatch[oldBFacei];
1060  subNewNbrProc[newBFacei] = sourceNewNbrProc[oldBFacei];
1061  }
1062  }
1063 
1064 
1065  subPointMaster = UIndirectList<label>(sourcePointMaster, pointMap);
1066 }
1067 
1068 
1069 // Find cells on mesh whose faceID/procID match the neighbour cell/proc of
1070 // domainMesh. Store the matching face.
1071 void Foam::fvMeshDistribute::findCouples
1072 (
1073  const primitiveMesh& mesh,
1074  const labelList& sourceFace,
1075  const labelList& sourceProc,
1076  const labelList& sourcePatch,
1077 
1078  const label domain,
1079  const primitiveMesh& domainMesh,
1080  const labelList& domainFace,
1081  const labelList& domainProc,
1082  const labelList& domainPatch,
1083 
1084  labelList& masterCoupledFaces,
1085  labelList& slaveCoupledFaces
1086 )
1087 {
1088  // Store domain neighbour as map so we can easily look for pair
1089  // with same face+proc.
1090  labelPairLookup map(domainFace.size());
1091 
1092  forAll(domainProc, bFacei)
1093  {
1094  if (domainProc[bFacei] != -1 && domainPatch[bFacei] == -1)
1095  {
1096  map.insert
1097  (
1098  labelPair(domainFace[bFacei], domainProc[bFacei]),
1099  bFacei
1100  );
1101  }
1102  }
1103 
1104 
1105  // Try to match mesh data.
1106 
1107  masterCoupledFaces.setSize(domainFace.size());
1108  slaveCoupledFaces.setSize(domainFace.size());
1109  label coupledI = 0;
1110 
1111  forAll(sourceFace, bFacei)
1112  {
1113  if (sourceProc[bFacei] != -1 && sourcePatch[bFacei] == -1)
1114  {
1115  labelPair myData(sourceFace[bFacei], sourceProc[bFacei]);
1116 
1117  const auto iter = map.cfind(myData);
1118 
1119  if (iter.good())
1120  {
1121  label nbrBFacei = *iter;
1122 
1123  masterCoupledFaces[coupledI] = mesh.nInternalFaces() + bFacei;
1124  slaveCoupledFaces[coupledI] =
1125  domainMesh.nInternalFaces()
1126  + nbrBFacei;
1127 
1128  coupledI++;
1129  }
1130  }
1131  }
1132 
1133  masterCoupledFaces.setSize(coupledI);
1134  slaveCoupledFaces.setSize(coupledI);
1135 
1136  if (debug)
1137  {
1138  Pout<< "findCouples : found " << coupledI
1139  << " faces that will be stitched" << nl << endl;
1140  }
1141 }
1142 
1143 
1144 void Foam::fvMeshDistribute::findCouples
1145 (
1146  const UPtrList<polyMesh>& meshes,
1147  const PtrList<labelList>& domainSourceFaces,
1148  const PtrList<labelList>& domainSourceProcs,
1149  const PtrList<labelList>& domainSourcePatchs,
1150 
1151  labelListList& localBoundaryFace,
1152  labelListList& remoteFaceProc,
1153  labelListList& remoteBoundaryFace
1154 )
1155 {
1156  // Collect all matching processor face pairs. These are all the
1157  // faces for which we have both sides
1158 
1159  // Pass 0: count number of inter-processor faces
1160  labelList nProcFaces(meshes.size(), 0);
1161  forAll(meshes, meshi)
1162  {
1163  if (meshes.set(meshi))
1164  {
1165  const labelList& domainProc = domainSourceProcs[meshi];
1166  const labelList& domainPatch = domainSourcePatchs[meshi];
1167 
1168  forAll(domainProc, bFacei)
1169  {
1170  if (domainProc[bFacei] != -1 && domainPatch[bFacei] == -1)
1171  {
1172  nProcFaces[meshi]++;
1173  }
1174  }
1175  }
1176  }
1177 
1178  if (debug)
1179  {
1180  Pout<< "fvMeshDistribute::findCouples : nProcFaces:"
1181  << flatOutput(nProcFaces) << endl;
1182  }
1183 
1184 
1185  // Size
1186  List<DynamicList<label>> dynLocalFace(Pstream::nProcs());
1187  List<DynamicList<label>> dynRemoteProc(Pstream::nProcs());
1188  List<DynamicList<label>> dynRemoteFace(Pstream::nProcs());
1189 
1190  forAll(meshes, meshi)
1191  {
1192  if (meshes.set(meshi))
1193  {
1194  dynLocalFace[meshi].reserve(nProcFaces[meshi]);
1195  dynRemoteProc[meshi].reserve(nProcFaces[meshi]);
1196  dynRemoteFace[meshi].reserve(nProcFaces[meshi]);
1197  }
1198  }
1199 
1200 
1201  // Insert all into big map. Find matches
1202  LabelPairMap<labelPair> map(2*sum(nProcFaces));
1203 
1204  nProcFaces = 0;
1205 
1206  forAll(meshes, meshi)
1207  {
1208  if (meshes.set(meshi))
1209  {
1210  const labelList& domainFace = domainSourceFaces[meshi];
1211  const labelList& domainProc = domainSourceProcs[meshi];
1212  const labelList& domainPatch = domainSourcePatchs[meshi];
1213 
1214  forAll(domainProc, bFacei)
1215  {
1216  if (domainProc[bFacei] != -1 && domainPatch[bFacei] == -1)
1217  {
1218  const labelPair key
1219  (
1220  domainProc[bFacei],
1221  domainFace[bFacei]
1222  );
1223  auto fnd = map.find(key);
1224 
1225  if (!fnd.good())
1226  {
1227  // Insert
1228  map.emplace(key, meshi, bFacei);
1229  }
1230  else
1231  {
1232  // Second occurence. Found match.
1233  const label matchProci = fnd().first();
1234  const label matchFacei = fnd().second();
1235 
1236  dynLocalFace[meshi].append(bFacei);
1237  dynRemoteProc[meshi].append(matchProci);
1238  dynRemoteFace[meshi].append(matchFacei);
1239  nProcFaces[meshi]++;
1240 
1241  dynLocalFace[matchProci].append(matchFacei);
1242  dynRemoteProc[matchProci].append(meshi);
1243  dynRemoteFace[matchProci].append(bFacei);
1244  nProcFaces[matchProci]++;
1245  }
1246  }
1247  }
1248  }
1249  }
1250 
1251  if (debug)
1252  {
1253  Pout<< "fvMeshDistribute::findCouples : stored procFaces:"
1254  << map.size() << endl;
1255  }
1256 
1257  localBoundaryFace.setSize(Pstream::nProcs());
1258  remoteFaceProc.setSize(Pstream::nProcs());
1259  remoteBoundaryFace.setSize(Pstream::nProcs());
1260  forAll(meshes, meshi)
1261  {
1262  if (meshes.set(meshi))
1263  {
1264  localBoundaryFace[meshi] = std::move(dynLocalFace[meshi]);
1265  remoteFaceProc[meshi] = std::move(dynRemoteProc[meshi]);
1266  remoteBoundaryFace[meshi] = std::move(dynRemoteFace[meshi]);
1267  }
1268  }
1269 
1270 
1271  if (debug)
1272  {
1273  Pout<< "fvMeshDistribute::findCouples : found matches:"
1274  << flatOutput(nProcFaces) << endl;
1275  }
1276 }
1277 
1278 
1279 // Map data on boundary faces to new mesh (resulting from adding two meshes)
1280 Foam::labelList Foam::fvMeshDistribute::mapBoundaryData
1281 (
1282  const primitiveMesh& mesh, // mesh after adding
1283  const mapAddedPolyMesh& map,
1284  const labelList& boundaryData0, // on mesh before adding
1285  const label nInternalFaces1,
1286  const labelList& boundaryData1 // on added mesh
1287 )
1288 {
1289  labelList newBoundaryData(mesh.nBoundaryFaces());
1290 
1291  forAll(boundaryData0, oldBFacei)
1292  {
1293  label newFacei = map.oldFaceMap()[oldBFacei + map.nOldInternalFaces()];
1294 
1295  // Face still exists (is necessary?) and still boundary face
1296  if (newFacei >= 0 && newFacei >= mesh.nInternalFaces())
1297  {
1298  newBoundaryData[newFacei - mesh.nInternalFaces()] =
1299  boundaryData0[oldBFacei];
1300  }
1301  }
1302 
1303  forAll(boundaryData1, addedBFacei)
1304  {
1305  label newFacei = map.addedFaceMap()[addedBFacei + nInternalFaces1];
1306 
1307  if (newFacei >= 0 && newFacei >= mesh.nInternalFaces())
1308  {
1309  newBoundaryData[newFacei - mesh.nInternalFaces()] =
1310  boundaryData1[addedBFacei];
1311  }
1312  }
1313 
1314  return newBoundaryData;
1315 }
1316 
1317 
1318 Foam::labelList Foam::fvMeshDistribute::mapPointData
1319 (
1320  const primitiveMesh& mesh, // mesh after adding
1321  const mapAddedPolyMesh& map,
1322  const labelList& boundaryData0, // on mesh before adding
1323  const labelList& boundaryData1 // on added mesh
1324 )
1325 {
1326  labelList newBoundaryData(mesh.nPoints());
1327 
1328  forAll(boundaryData0, oldPointi)
1329  {
1330  label newPointi = map.oldPointMap()[oldPointi];
1331 
1332  // Point still exists (is necessary?)
1333  if (newPointi >= 0)
1334  {
1335  newBoundaryData[newPointi] = boundaryData0[oldPointi];
1336  }
1337  }
1338 
1339  forAll(boundaryData1, addedPointi)
1340  {
1341  label newPointi = map.addedPointMap()[addedPointi];
1342 
1343  if (newPointi >= 0)
1344  {
1345  newBoundaryData[newPointi] = boundaryData1[addedPointi];
1346  }
1347  }
1348 
1349  return newBoundaryData;
1350 }
1351 
1352 
1353 // Remove cells. Add all exposed faces to patch oldInternalPatchi
1354 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::doRemoveCells
1355 (
1356  const labelList& cellsToRemove,
1357  const label oldInternalPatchi
1358 )
1359 {
1360  // Mesh change engine
1361  polyTopoChange meshMod(mesh_);
1362 
1363  // Cell removal topo engine. Do NOT synchronize parallel since
1364  // we are doing a local cell removal.
1365  removeCells cellRemover(mesh_, false);
1366 
1367  // Get all exposed faces
1368  labelList exposedFaces(cellRemover.getExposedFaces(cellsToRemove));
1369 
1370  // Insert the topo changes
1371  cellRemover.setRefinement
1372  (
1373  cellsToRemove,
1374  exposedFaces,
1375  labelList(exposedFaces.size(), oldInternalPatchi), // patch for exposed
1376  // faces.
1377  meshMod
1378  );
1379 
1380 
1382  //tmp<surfaceScalarField> sfld(generateTestField(mesh_));
1383 
1384  // Save internal fields (note: not as DimensionedFields since would
1385  // get mapped)
1386  PtrList<Field<scalar>> sFlds;
1387  saveInternalFields(sFlds);
1388  PtrList<Field<vector>> vFlds;
1389  saveInternalFields(vFlds);
1390  PtrList<Field<sphericalTensor>> sptFlds;
1391  saveInternalFields(sptFlds);
1392  PtrList<Field<symmTensor>> sytFlds;
1393  saveInternalFields(sytFlds);
1394  PtrList<Field<tensor>> tFlds;
1395  saveInternalFields(tFlds);
1396 
1397  // Change the mesh. No inflation. Note: no parallel comms allowed.
1398 
1399  // TBD: temporarily unset mesh moving to avoid problems in meshflux
1400  // mapping. To be investigated.
1401  const bool oldMoving = mesh_.moving(false);
1402  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, false);
1403  mesh_.moving(oldMoving);
1404 
1405  // Update fields
1406  mesh_.updateMesh(map());
1407 
1408 
1409  // Any exposed faces in a surfaceField will not be mapped. Map the value
1410  // of these separately (until there is support in all PatchFields for
1411  // mapping from internal faces ...)
1412 
1413  mapExposedFaces(map(), sFlds);
1414  mapExposedFaces(map(), vFlds);
1415  mapExposedFaces(map(), sptFlds);
1416  mapExposedFaces(map(), sytFlds);
1417  mapExposedFaces(map(), tFlds);
1418 
1419 
1421  //testField(sfld);
1422 
1423 
1424  // Move mesh (since morphing does not do this)
1425  if (map().hasMotionPoints())
1426  {
1427  mesh_.movePoints(map().preMotionPoints());
1428  }
1429 
1430 
1431  return map;
1432 }
1433 
1434 
1435 // Delete and add processor patches. Changes mesh and returns per neighbour proc
1436 // the processor patchID.
1437 void Foam::fvMeshDistribute::addProcPatches
1438 (
1439  const labelList& nbrProc, // processor that neighbour is now on
1440  const labelList& referPatchID, // patchID (or -1) I originated from
1441  List<Map<label>>& procPatchID
1442 )
1443 {
1444  // Now use the neighbourFace/Proc to repatch the mesh. These lists
1445  // contain for all current boundary faces the global patchID (for non-proc
1446  // patch) or the processor.
1447 
1448  // Determine a visit order such that the processor patches get added
1449  // in order of increasing neighbour processor (and for same neighbour
1450  // processor (in case of processor cyclics) in order of increasing
1451  // 'refer' patch)
1452  labelList indices;
1453  sortedOrder(nbrProc, indices, lessProcPatches(nbrProc, referPatchID));
1454 
1455  procPatchID.setSize(Pstream::nProcs());
1456 
1457  forAll(indices, i)
1458  {
1459  label bFacei = indices[i];
1460  label proci = nbrProc[bFacei];
1461 
1462  if (proci != -1 && proci != Pstream::myProcNo())
1463  {
1464  if (!procPatchID[proci].found(referPatchID[bFacei]))
1465  {
1466  // No patch for neighbour yet. Is either a normal processor
1467  // patch or a processorCyclic patch.
1468 
1469  if (referPatchID[bFacei] == -1)
1470  {
1471  // Ordinary processor boundary
1472 
1473  processorPolyPatch pp
1474  (
1475  0, // size
1476  mesh_.nFaces(),
1477  mesh_.boundaryMesh().size(),
1478  mesh_.boundaryMesh(),
1479  Pstream::myProcNo(),
1480  proci
1481  );
1482 
1483  procPatchID[proci].insert
1484  (
1485  referPatchID[bFacei],
1486  fvMeshTools::addPatch
1487  (
1488  mesh_,
1489  pp,
1490  dictionary(), // optional per field patchField
1491  processorFvPatchField<scalar>::typeName,
1492  false // not parallel sync
1493  )
1494  );
1495  }
1496  else
1497  {
1498  const coupledPolyPatch& pcPatch
1499  = refCast<const coupledPolyPatch>
1500  (
1501  mesh_.boundaryMesh()[referPatchID[bFacei]]
1502  );
1503  processorCyclicPolyPatch pp
1504  (
1505  0, // size
1506  mesh_.nFaces(),
1507  mesh_.boundaryMesh().size(),
1508  mesh_.boundaryMesh(),
1509  Pstream::myProcNo(),
1510  proci,
1511  pcPatch.name(),
1512  pcPatch.transform()
1513  );
1514 
1515  procPatchID[proci].insert
1516  (
1517  referPatchID[bFacei],
1518  fvMeshTools::addPatch
1519  (
1520  mesh_,
1521  pp,
1522  dictionary(), // optional per field patchField
1523  processorCyclicFvPatchField<scalar>::typeName,
1524  false // not parallel sync
1525  )
1526  );
1527  }
1528  }
1529  }
1530  }
1531 }
1532 
1533 
1534 // Get boundary faces to be repatched. Is -1 or new patchID
1535 Foam::labelList Foam::fvMeshDistribute::getBoundaryPatch
1536 (
1537  const labelList& nbrProc, // new processor per boundary face
1538  const labelList& referPatchID, // patchID (or -1) I originated from
1539  const List<Map<label>>& procPatchID // per proc the new procPatches
1540 )
1541 {
1542  labelList patchIDs(nbrProc);
1543 
1544  forAll(nbrProc, bFacei)
1545  {
1546  if (nbrProc[bFacei] == Pstream::myProcNo())
1547  {
1548  label origPatchi = referPatchID[bFacei];
1549  patchIDs[bFacei] = origPatchi;
1550  }
1551  else if (nbrProc[bFacei] != -1)
1552  {
1553  label origPatchi = referPatchID[bFacei];
1554  patchIDs[bFacei] = procPatchID[nbrProc[bFacei]][origPatchi];
1555  }
1556  else
1557  {
1558  patchIDs[bFacei] = -1;
1559  }
1560  }
1561  return patchIDs;
1562 }
1563 
1564 
1565 // Send mesh and coupling data.
1566 void Foam::fvMeshDistribute::sendMesh
1567 (
1568  const label domain,
1569  const fvMesh& mesh,
1570 
1571  const wordList& pointZoneNames,
1572  const wordList& faceZoneNames,
1573  const wordList& cellZoneNames,
1574 
1575  const labelList& sourceFace,
1576  const labelList& sourceProc,
1577  const labelList& sourcePatch,
1578  const labelList& sourceNewNbrProc,
1579  const labelList& sourcePointMaster,
1580  Ostream& toDomain
1581 )
1582 {
1583  if (debug)
1584  {
1585  Pout<< "Sending to domain " << domain << nl
1586  << " nPoints:" << mesh.nPoints() << nl
1587  << " nFaces:" << mesh.nFaces() << nl
1588  << " nCells:" << mesh.nCells() << nl
1589  << " nPatches:" << mesh.boundaryMesh().size() << nl
1590  << endl;
1591  }
1592 
1593  // Assume sparse, possibly overlapping point zones. Get contents
1594  // in merged-zone indices.
1595  CompactListList<label> zonePoints;
1596  {
1597  const pointZoneMesh& pointZones = mesh.pointZones();
1598 
1599  labelList rowSizes(pointZoneNames.size(), Zero);
1600 
1601  forAll(pointZoneNames, nameI)
1602  {
1603  label myZoneID = pointZones.findZoneID(pointZoneNames[nameI]);
1604 
1605  if (myZoneID != -1)
1606  {
1607  rowSizes[nameI] = pointZones[myZoneID].size();
1608  }
1609  }
1610  zonePoints.setSize(rowSizes);
1611 
1612  forAll(pointZoneNames, nameI)
1613  {
1614  label myZoneID = pointZones.findZoneID(pointZoneNames[nameI]);
1615 
1616  if (myZoneID != -1)
1617  {
1618  zonePoints[nameI].deepCopy(pointZones[myZoneID]);
1619  }
1620  }
1621  }
1622 
1623  // Assume sparse, possibly overlapping face zones
1624  CompactListList<label> zoneFaces;
1625  CompactListList<bool> zoneFaceFlip;
1626  {
1627  const faceZoneMesh& faceZones = mesh.faceZones();
1628 
1629  labelList rowSizes(faceZoneNames.size(), Zero);
1630 
1631  forAll(faceZoneNames, nameI)
1632  {
1633  label myZoneID = faceZones.findZoneID(faceZoneNames[nameI]);
1634 
1635  if (myZoneID != -1)
1636  {
1637  rowSizes[nameI] = faceZones[myZoneID].size();
1638  }
1639  }
1640 
1641  zoneFaces.setSize(rowSizes);
1642  zoneFaceFlip.setSize(rowSizes);
1643 
1644  forAll(faceZoneNames, nameI)
1645  {
1646  label myZoneID = faceZones.findZoneID(faceZoneNames[nameI]);
1647 
1648  if (myZoneID != -1)
1649  {
1650  zoneFaces[nameI].deepCopy(faceZones[myZoneID]);
1651  zoneFaceFlip[nameI].deepCopy(faceZones[myZoneID].flipMap());
1652  }
1653  }
1654  }
1655 
1656  // Assume sparse, possibly overlapping cell zones
1657  CompactListList<label> zoneCells;
1658  {
1659  const cellZoneMesh& cellZones = mesh.cellZones();
1660 
1661  labelList rowSizes(cellZoneNames.size(), Zero);
1662 
1663  forAll(cellZoneNames, nameI)
1664  {
1665  label myZoneID = cellZones.findZoneID(cellZoneNames[nameI]);
1666 
1667  if (myZoneID != -1)
1668  {
1669  rowSizes[nameI] = cellZones[myZoneID].size();
1670  }
1671  }
1672 
1673  zoneCells.setSize(rowSizes);
1674 
1675  forAll(cellZoneNames, nameI)
1676  {
1677  label myZoneID = cellZones.findZoneID(cellZoneNames[nameI]);
1678 
1679  if (myZoneID != -1)
1680  {
1681  zoneCells[nameI].deepCopy(cellZones[myZoneID]);
1682  }
1683  }
1684  }
1686  //labelList cellZoneID;
1687  //if (hasCellZones)
1688  //{
1689  // cellZoneID.setSize(mesh.nCells());
1690  // cellZoneID = -1;
1691  //
1692  // const cellZoneMesh& cellZones = mesh.cellZones();
1693  //
1694  // forAll(cellZones, zoneI)
1695  // {
1696  // labelUIndList(cellZoneID, cellZones[zoneI]) = zoneI;
1697  // }
1698  //}
1699 
1700  // Send
1701  toDomain
1702  << mesh.points()
1703  // Send as (offsets/values)
1704  << CompactListList<label>::pack<face>(mesh.faces())
1705  << mesh.faceOwner()
1706  << mesh.faceNeighbour()
1707  << mesh.boundaryMesh()
1708 
1709  << zonePoints
1710  << zoneFaces
1711  << zoneFaceFlip
1712  << zoneCells
1713 
1714  << sourceFace
1715  << sourceProc
1716  << sourcePatch
1717  << sourceNewNbrProc
1718  << sourcePointMaster;
1719 
1720 
1721  if (debug)
1722  {
1723  Pout<< "Started sending mesh to domain " << domain
1724  << endl;
1725  }
1726 }
1727 
1728 
1729 // Receive mesh. Opposite of sendMesh
1730 Foam::autoPtr<Foam::fvMesh> Foam::fvMeshDistribute::receiveMesh
1731 (
1732  const label domain,
1733  const wordList& pointZoneNames,
1734  const wordList& faceZoneNames,
1735  const wordList& cellZoneNames,
1736  const Time& runTime,
1737  labelList& domainSourceFace,
1738  labelList& domainSourceProc,
1739  labelList& domainSourcePatch,
1740  labelList& domainSourceNewNbrProc,
1741  labelList& domainSourcePointMaster,
1742  Istream& fromNbr
1743 )
1744 {
1745  pointField domainPoints(fromNbr);
1746  // Receive as (offsets/values), unpack to faceList
1747  faceList domainFaces = CompactListList<label>(fromNbr).unpack<face>();
1748  labelList domainAllOwner(fromNbr);
1749  labelList domainAllNeighbour(fromNbr);
1750  PtrList<entry> patchEntries(fromNbr);
1751 
1752  CompactListList<label> zonePoints(fromNbr);
1753  CompactListList<label> zoneFaces(fromNbr);
1754  CompactListList<bool> zoneFaceFlip(fromNbr);
1755  CompactListList<label> zoneCells(fromNbr);
1756 
1757  fromNbr
1758  >> domainSourceFace
1759  >> domainSourceProc
1760  >> domainSourcePatch
1761  >> domainSourceNewNbrProc
1762  >> domainSourcePointMaster;
1763 
1764  // Construct fvMesh
1765  auto domainMeshPtr = autoPtr<fvMesh>::New
1766  (
1767  IOobject
1768  (
1769  polyMesh::defaultRegion,
1770  runTime.timeName(),
1771  runTime,
1772  IOobject::NO_READ
1773  ),
1774  std::move(domainPoints),
1775  std::move(domainFaces),
1776  std::move(domainAllOwner),
1777  std::move(domainAllNeighbour),
1778  false // no parallel comms
1779  );
1780  fvMesh& domainMesh = *domainMeshPtr;
1781 
1782  polyPatchList patches(patchEntries.size());
1783 
1784  forAll(patchEntries, patchi)
1785  {
1786  patches.set
1787  (
1788  patchi,
1789  polyPatch::New
1790  (
1791  patchEntries[patchi].keyword(),
1792  patchEntries[patchi].dict(),
1793  patchi,
1794  domainMesh.boundaryMesh()
1795  )
1796  );
1797  }
1798  // Add patches; no parallel comms
1799  domainMesh.addFvPatches(patches, false);
1800 
1801  // Construct zones
1802  List<pointZone*> pZonePtrs(pointZoneNames.size());
1803  forAll(pZonePtrs, i)
1804  {
1805  pZonePtrs[i] = new pointZone
1806  (
1807  pointZoneNames[i],
1808  zonePoints[i],
1809  i,
1810  domainMesh.pointZones()
1811  );
1812  }
1813 
1814  List<faceZone*> fZonePtrs(faceZoneNames.size());
1815  forAll(fZonePtrs, i)
1816  {
1817  fZonePtrs[i] = new faceZone
1818  (
1819  faceZoneNames[i],
1820  zoneFaces[i],
1821  zoneFaceFlip[i],
1822  i,
1823  domainMesh.faceZones()
1824  );
1825  }
1826 
1827  List<cellZone*> cZonePtrs(cellZoneNames.size());
1828  forAll(cZonePtrs, i)
1829  {
1830  cZonePtrs[i] = new cellZone
1831  (
1832  cellZoneNames[i],
1833  zoneCells[i],
1834  i,
1835  domainMesh.cellZones()
1836  );
1837  }
1838  domainMesh.addZones(pZonePtrs, fZonePtrs, cZonePtrs);
1839 
1840  return domainMeshPtr;
1841 }
1842 
1843 
1844 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
1845 
1846 Foam::fvMeshDistribute::fvMeshDistribute(fvMesh& mesh)
1847 :
1848  mesh_(mesh)
1849 {}
1850 
1851 
1852 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1853 
1854 Foam::labelList Foam::fvMeshDistribute::countCells
1855 (
1856  const labelList& distribution
1857 )
1858 {
1859  labelList nCells(Pstream::nProcs(), Zero);
1860  forAll(distribution, celli)
1861  {
1862  label newProc = distribution[celli];
1863 
1864  if (newProc < 0 || newProc >= Pstream::nProcs())
1865  {
1866  FatalErrorInFunction
1867  << "Distribution should be in range 0.." << Pstream::nProcs()-1
1868  << endl
1869  << "At index " << celli << " distribution:" << newProc
1870  << abort(FatalError);
1871  }
1872  nCells[newProc]++;
1873  }
1874  return nCells;
1875 }
1876 
1877 
1878 Foam::autoPtr<Foam::mapDistributePolyMesh> Foam::fvMeshDistribute::distribute
1879 (
1880  const labelList& distribution
1881 )
1882 {
1883  // Some checks on distribution
1884  if (distribution.size() != mesh_.nCells())
1885  {
1886  FatalErrorInFunction
1887  << "Size of distribution:"
1888  << distribution.size() << " mesh nCells:" << mesh_.nCells()
1889  << abort(FatalError);
1890  }
1891 
1892 
1893  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
1894 
1895  // Check all processors have same non-proc patches in same order.
1896  if (patches.checkParallelSync(true))
1897  {
1898  FatalErrorInFunction
1899  << "This application requires all non-processor patches"
1900  << " to be present in the same order on all patches" << nl
1901  << "followed by the processor patches (which of course are unique)."
1902  << nl
1903  << "Local patches:" << mesh_.boundaryMesh().names()
1904  << abort(FatalError);
1905  }
1906 
1907  // Save some data for mapping later on
1908  const label nOldPoints(mesh_.nPoints());
1909  const label nOldFaces(mesh_.nFaces());
1910  const label nOldCells(mesh_.nCells());
1911  labelList oldPatchStarts(patches.size());
1912  labelList oldPatchNMeshPoints(patches.size());
1913  forAll(patches, patchi)
1914  {
1915  oldPatchStarts[patchi] = patches[patchi].start();
1916  oldPatchNMeshPoints[patchi] = patches[patchi].nPoints();
1917  }
1918 
1919 
1920 
1921  // Short circuit trivial case.
1922  if (!Pstream::parRun())
1923  {
1924  // Collect all maps and return
1925  return autoPtr<mapDistributePolyMesh>::New
1926  (
1927  mesh_,
1928 
1929  nOldPoints,
1930  nOldFaces,
1931  nOldCells,
1932  std::move(oldPatchStarts),
1933  std::move(oldPatchNMeshPoints),
1934 
1935  labelListList(one(), identity(mesh_.nPoints())), //subPointMap
1936  labelListList(one(), identity(mesh_.nFaces())), //subFaceMap
1937  labelListList(one(), identity(mesh_.nCells())), //subCellMap
1938  labelListList(one(), identity(patches.size())), //subPatchMap
1939 
1940  labelListList(one(), identity(mesh_.nPoints())), //pointMap
1941  labelListList(one(), identity(mesh_.nFaces())), //faceMap
1942  labelListList(one(), identity(mesh_.nCells())), //cellMap
1943  labelListList(one(), identity(patches.size())) //patchMap
1944  );
1945  }
1946 
1947 
1948  // Collect any zone names over all processors and shuffle zones accordingly
1949  // Note that this is not necessary for redistributePar since that already
1950  // checks for it. However other use (e.g. mesh generators) might not.
1951  const wordList pointZoneNames(mergeWordList(mesh_.pointZones().names()));
1952  reorderZones<pointZone>(pointZoneNames, mesh_.pointZones());
1953 
1954  const wordList faceZoneNames(mergeWordList(mesh_.faceZones().names()));
1955  reorderZones<faceZone>(faceZoneNames, mesh_.faceZones());
1956 
1957  const wordList cellZoneNames(mergeWordList(mesh_.cellZones().names()));
1958  reorderZones<cellZone>(cellZoneNames, mesh_.cellZones());
1959 
1960 
1961  // Local environment of all boundary faces
1962  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1963 
1964  // A face is uniquely defined by
1965  // - proc
1966  // - local face no
1967  //
1968  // To glue the parts of meshes which can get sent from anywhere we
1969  // need to know on boundary faces what the above tuple on both sides is.
1970  // So we need to maintain:
1971  // - original face
1972  // - original processor id (= trivial)
1973  // For coupled boundaries (where the faces are 'duplicate') we take the
1974  // lowest numbered processor as the data to store.
1975  //
1976  // Additionally to create the procboundaries we need to know where the owner
1977  // cell on the other side now is: newNeighbourProc.
1978  //
1979 
1980  // physical boundary:
1981  // sourceProc = -1
1982  // sourceNewNbrProc = -1
1983  // sourceFace = -1
1984  // sourcePatch = patchID
1985  // processor boundary:
1986  // sourceProc = proc (on owner side)
1987  // sourceNewNbrProc = distribution of coupled cell
1988  // sourceFace = face (on owner side)
1989  // sourcePatch = -1
1990  // ?cyclic:
1991  // ? sourceProc = proc
1992  // ? sourceNewNbrProc = distribution of coupled cell
1993  // ? sourceFace = face (on owner side)
1994  // ? sourcePatch = patchID
1995  // processor-cyclic boundary:
1996  // sourceProc = proc (on owner side)
1997  // sourceNewNbrProc = distribution of coupled cell
1998  // sourceFace = face (on owner side)
1999  // sourcePatch = patchID
2000 
2001  labelList sourcePatch;
2002  labelList sourceFace;
2003  labelList sourceProc;
2004  labelList sourceNewNbrProc;
2005  labelList sourcePointMaster;
2006  getCouplingData
2007  (
2008  distribution,
2009  sourceFace,
2010  sourceProc,
2011  sourcePatch,
2012  sourceNewNbrProc,
2013  sourcePointMaster
2014  );
2015 
2016 
2017  // Remove meshPhi. Since this would otherwise disappear anyway
2018  // during topo changes and we have to guarantee that all the fields
2019  // can be sent.
2020  mesh_.clearOut();
2021  mesh_.resetMotion();
2022 
2023  // Get data to send. Make sure is synchronised
2024 
2025  HashTable<wordList> allFieldNames;
2026 
2027  getFieldNames<volScalarField>(mesh_, allFieldNames);
2028  getFieldNames<volVectorField>(mesh_, allFieldNames);
2029  getFieldNames<volSphericalTensorField>(mesh_, allFieldNames);
2030  getFieldNames<volSymmTensorField>(mesh_, allFieldNames);
2031  getFieldNames<volTensorField>(mesh_, allFieldNames);
2032 
2033  getFieldNames<surfaceScalarField>(mesh_, allFieldNames);
2034  getFieldNames<surfaceVectorField>(mesh_, allFieldNames);
2035  getFieldNames<surfaceSphericalTensorField>(mesh_, allFieldNames);
2036  getFieldNames<surfaceSymmTensorField>(mesh_, allFieldNames);
2037  getFieldNames<surfaceTensorField>(mesh_, allFieldNames);
2038 
2039  getFieldNames<volScalarField::Internal>
2040  (
2041  mesh_,
2042  allFieldNames,
2043  volScalarField::typeName
2044  );
2045  getFieldNames<volVectorField::Internal>
2046  (
2047  mesh_,
2048  allFieldNames,
2049  volVectorField::typeName
2050  );
2051  getFieldNames<volSphericalTensorField::Internal>
2052  (
2053  mesh_,
2054  allFieldNames,
2055  volSphericalTensorField::typeName
2056  );
2057  getFieldNames<volSymmTensorField::Internal>
2058  (
2059  mesh_,
2060  allFieldNames,
2061  volSymmTensorField::typeName
2062  );
2063  getFieldNames<volTensorField::Internal>
2064  (
2065  mesh_,
2066  allFieldNames,
2067  volTensorField::typeName
2068  );
2069 
2070 
2071  // Find patch to temporarily put exposed and processor faces into.
2072  const label oldInternalPatchi = findNonEmptyPatch();
2073 
2074 
2075  // Delete processor patches, starting from the back. Move all faces into
2076  // oldInternalPatchi.
2077  labelList repatchFaceMap;
2078  {
2079  autoPtr<mapPolyMesh> repatchMap = deleteProcPatches(oldInternalPatchi);
2080 
2081  // Store face map (only face ordering that changed)
2082  repatchFaceMap = repatchMap().faceMap();
2083 
2084  // Reorder all boundary face data (sourceProc, sourceFace etc.)
2085  labelList bFaceMap
2086  (
2087  SubList<label>
2088  (
2089  repatchMap().reverseFaceMap(),
2090  mesh_.nBoundaryFaces(),
2091  mesh_.nInternalFaces()
2092  )
2093  - mesh_.nInternalFaces()
2094  );
2095 
2096  inplaceReorder(bFaceMap, sourceFace);
2097  inplaceReorder(bFaceMap, sourceProc);
2098  inplaceReorder(bFaceMap, sourcePatch);
2099  inplaceReorder(bFaceMap, sourceNewNbrProc);
2100  }
2101 
2102 
2103 
2104  // Print a bit.
2105  if (debug)
2106  {
2107  Pout<< nl << "MESH WITH PROC PATCHES DELETED:" << endl;
2108  printMeshInfo(mesh_);
2109  printFieldInfo<volScalarField>(mesh_);
2110  printFieldInfo<volVectorField>(mesh_);
2111  printFieldInfo<volSphericalTensorField>(mesh_);
2112  printFieldInfo<volSymmTensorField>(mesh_);
2113  printFieldInfo<volTensorField>(mesh_);
2114  printFieldInfo<surfaceScalarField>(mesh_);
2115  printFieldInfo<surfaceVectorField>(mesh_);
2116  printFieldInfo<surfaceSphericalTensorField>(mesh_);
2117  printFieldInfo<surfaceSymmTensorField>(mesh_);
2118  printFieldInfo<surfaceTensorField>(mesh_);
2119  printIntFieldInfo<volScalarField::Internal>(mesh_);
2120  printIntFieldInfo<volVectorField::Internal>(mesh_);
2121  printIntFieldInfo<volSphericalTensorField::Internal>(mesh_);
2122  printIntFieldInfo<volSymmTensorField::Internal>(mesh_);
2123  printIntFieldInfo<volTensorField::Internal>(mesh_);
2124  Pout<< nl << endl;
2125  }
2126 
2127 
2128 
2129  // Maps from subsetted mesh (that is sent) back to original maps
2130  labelListList subCellMap(Pstream::nProcs());
2131  labelListList subFaceMap(Pstream::nProcs());
2132  labelListList subPointMap(Pstream::nProcs());
2133  labelListList subPatchMap(Pstream::nProcs());
2134  // Maps from subsetted mesh to reconstructed mesh
2135  labelListList constructCellMap(Pstream::nProcs());
2136  labelListList constructFaceMap(Pstream::nProcs());
2137  labelListList constructPointMap(Pstream::nProcs());
2138  labelListList constructPatchMap(Pstream::nProcs());
2139 
2140 
2141  // Find out schedule
2142  // ~~~~~~~~~~~~~~~~~
2143 
2144  labelList nSendCells(countCells(distribution));
2145  labelList nRecvCells(Pstream::nProcs());
2146  UPstream::allToAll(nSendCells, nRecvCells);
2147 
2148  // Allocate buffers
2149  PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
2150 
2151 
2152  // What to send to neighbouring domains
2153  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2154 
2155  // Disable parallel.
2156  const bool oldParRun = UPstream::parRun(false);
2157 
2158  forAll(nSendCells, recvProc)
2159  {
2160  if (recvProc != Pstream::myProcNo() && nSendCells[recvProc] > 0)
2161  {
2162  // Send to recvProc
2163 
2164  if (debug)
2165  {
2166  Pout<< nl
2167  << "SUBSETTING FOR DOMAIN " << recvProc
2168  << " cells to send:"
2169  << nSendCells[recvProc]
2170  << nl << endl;
2171  }
2172 
2173  // Pstream for sending mesh and fields
2174  //OPstream str(Pstream::commsTypes::blocking, recvProc);
2175  UOPstream str(recvProc, pBufs);
2176 
2177  // Mesh subsetting engine - subset the cells of the current domain.
2178  fvMeshSubset subsetter
2179  (
2180  mesh_,
2181  recvProc,
2182  distribution,
2183  oldInternalPatchi, // oldInternalFaces patch
2184  false // no parallel sync
2185  );
2186 
2187  subCellMap[recvProc] = subsetter.cellMap();
2188  subFaceMap[recvProc] = subsetter.faceFlipMap();
2189  inplaceRenumberWithFlip
2190  (
2191  repatchFaceMap,
2192  false, // oldToNew has flip
2193  true, // subFaceMap has flip
2194  subFaceMap[recvProc]
2195  );
2196  subPointMap[recvProc] = subsetter.pointMap();
2197  subPatchMap[recvProc] = subsetter.patchMap();
2198 
2199 
2200  // Subset the boundary fields (owner/neighbour/processor)
2201  labelList procSourceFace;
2202  labelList procSourceProc;
2203  labelList procSourcePatch;
2204  labelList procSourceNewNbrProc;
2205  labelList procSourcePointMaster;
2206 
2207  subsetCouplingData
2208  (
2209  subsetter.subMesh(),
2210  subsetter.pointMap(), // from subMesh to mesh
2211  subsetter.faceMap(), // ,, ,,
2212  subsetter.cellMap(), // ,, ,,
2213 
2214  distribution, // old mesh distribution
2215  mesh_.faceOwner(), // old owner
2216  mesh_.faceNeighbour(),
2217  mesh_.nInternalFaces(),
2218 
2219  sourceFace,
2220  sourceProc,
2221  sourcePatch,
2222  sourceNewNbrProc,
2223  sourcePointMaster,
2224 
2225  procSourceFace,
2226  procSourceProc,
2227  procSourcePatch,
2228  procSourceNewNbrProc,
2229  procSourcePointMaster
2230  );
2231 
2232 
2233  // Send to neighbour
2234  sendMesh
2235  (
2236  recvProc,
2237  subsetter.subMesh(),
2238 
2239  pointZoneNames,
2240  faceZoneNames,
2241  cellZoneNames,
2242 
2243  procSourceFace,
2244  procSourceProc,
2245  procSourcePatch,
2246  procSourceNewNbrProc,
2247  procSourcePointMaster,
2248 
2249  str
2250  );
2251 
2252  // volFields
2253  sendFields<volScalarField>
2254  (
2255  recvProc,
2256  allFieldNames,
2257  subsetter,
2258  str
2259  );
2260  sendFields<volVectorField>
2261  (
2262  recvProc,
2263  allFieldNames,
2264  subsetter,
2265  str
2266  );
2267  sendFields<volSphericalTensorField>
2268  (
2269  recvProc,
2270  allFieldNames,
2271  subsetter,
2272  str
2273  );
2274  sendFields<volSymmTensorField>
2275  (
2276  recvProc,
2277  allFieldNames,
2278  subsetter,
2279  str
2280  );
2281  sendFields<volTensorField>
2282  (
2283  recvProc,
2284  allFieldNames,
2285  subsetter,
2286  str
2287  );
2288 
2289  // surfaceFields
2290  sendFields<surfaceScalarField>
2291  (
2292  recvProc,
2293  allFieldNames,
2294  subsetter,
2295  str
2296  );
2297  sendFields<surfaceVectorField>
2298  (
2299  recvProc,
2300  allFieldNames,
2301  subsetter,
2302  str
2303  );
2304  sendFields<surfaceSphericalTensorField>
2305  (
2306  recvProc,
2307  allFieldNames,
2308  subsetter,
2309  str
2310  );
2311  sendFields<surfaceSymmTensorField>
2312  (
2313  recvProc,
2314  allFieldNames,
2315  subsetter,
2316  str
2317  );
2318  sendFields<surfaceTensorField>
2319  (
2320  recvProc,
2321  allFieldNames,
2322  subsetter,
2323  str
2324  );
2325 
2326  // Dimensioned fields
2327  sendFields<volScalarField::Internal>
2328  (
2329  recvProc,
2330  allFieldNames,
2331  subsetter,
2332  str
2333  );
2334  sendFields<volVectorField::Internal>
2335  (
2336  recvProc,
2337  allFieldNames,
2338  subsetter,
2339  str
2340  );
2341  sendFields<volSphericalTensorField::Internal>
2342  (
2343  recvProc,
2344  allFieldNames,
2345  subsetter,
2346  str
2347  );
2348  sendFields<volSymmTensorField::Internal>
2349  (
2350  recvProc,
2351  allFieldNames,
2352  subsetter,
2353  str
2354  );
2355  sendFields<volTensorField::Internal>
2356  (
2357  recvProc,
2358  allFieldNames,
2359  subsetter,
2360  str
2361  );
2362  }
2363  }
2364 
2365 
2366  UPstream::parRun(oldParRun); // Restore parallel state
2367 
2368  if (debug)
2369  {
2370  Pout<< "Starting sending" << endl;
2371  }
2372 
2373  pBufs.finishedSends();
2374 
2375  if (debug)
2376  {
2377  Pout<< "Finished sending and receiving : "
2378  << flatOutput(pBufs.recvDataCounts()) << endl;
2379  }
2380 
2381 
2382  // Subset the part that stays
2383  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
2384 
2385  {
2386  // Save old mesh maps before changing mesh
2387  const labelList oldFaceOwner(mesh_.faceOwner());
2388  const labelList oldFaceNeighbour(mesh_.faceNeighbour());
2389  const label oldInternalFaces = mesh_.nInternalFaces();
2390 
2391  // Remove cells.
2392  autoPtr<mapPolyMesh> subMap
2393  (
2394  doRemoveCells
2395  (
2396  select(false, distribution, Pstream::myProcNo()),
2397  oldInternalPatchi
2398  )
2399  );
2400 
2401  // Addressing from subsetted mesh
2402  subCellMap[Pstream::myProcNo()] = subMap().cellMap();
2403  subFaceMap[Pstream::myProcNo()] = renumber
2404  (
2405  repatchFaceMap,
2406  subMap().faceMap()
2407  );
2408  // Insert the sign bit from face flipping
2409  labelList& faceMap = subFaceMap[Pstream::myProcNo()];
2410  forAll(faceMap, faceI)
2411  {
2412  faceMap[faceI] += 1;
2413  }
2414  const labelHashSet& flip = subMap().flipFaceFlux();
2415  for (const label facei : flip)
2416  {
2417  faceMap[facei] = -faceMap[facei];
2418  }
2419  subPointMap[Pstream::myProcNo()] = subMap().pointMap();
2420  subPatchMap[Pstream::myProcNo()] = identity(patches.size());
2421 
2422  // Subset the mesh data: neighbourCell/neighbourProc fields
2423  labelList domainSourceFace;
2424  labelList domainSourceProc;
2425  labelList domainSourcePatch;
2426  labelList domainSourceNewNbrProc;
2427  labelList domainSourcePointMaster;
2428 
2429  subsetCouplingData
2430  (
2431  mesh_, // new mesh
2432  subMap().pointMap(), // from new to original mesh
2433  subMap().faceMap(), // from new to original mesh
2434  subMap().cellMap(),
2435 
2436  distribution, // distribution before subsetting
2437  oldFaceOwner, // owner before subsetting
2438  oldFaceNeighbour, // neighbour ,,
2439  oldInternalFaces, // nInternalFaces ,,
2440 
2441  sourceFace,
2442  sourceProc,
2443  sourcePatch,
2444  sourceNewNbrProc,
2445  sourcePointMaster,
2446 
2447  domainSourceFace,
2448  domainSourceProc,
2449  domainSourcePatch,
2450  domainSourceNewNbrProc,
2451  domainSourcePointMaster
2452  );
2453 
2454  sourceFace.transfer(domainSourceFace);
2455  sourceProc.transfer(domainSourceProc);
2456  sourcePatch.transfer(domainSourcePatch);
2457  sourceNewNbrProc.transfer(domainSourceNewNbrProc);
2458  sourcePointMaster.transfer(domainSourcePointMaster);
2459  }
2460 
2461 
2462  // Print a bit.
2463  if (debug)
2464  {
2465  Pout<< nl << "STARTING MESH:" << endl;
2466  printMeshInfo(mesh_);
2467  printFieldInfo<volScalarField>(mesh_);
2468  printFieldInfo<volVectorField>(mesh_);
2469  printFieldInfo<volSphericalTensorField>(mesh_);
2470  printFieldInfo<volSymmTensorField>(mesh_);
2471  printFieldInfo<volTensorField>(mesh_);
2472  printFieldInfo<surfaceScalarField>(mesh_);
2473  printFieldInfo<surfaceVectorField>(mesh_);
2474  printFieldInfo<surfaceSphericalTensorField>(mesh_);
2475  printFieldInfo<surfaceSymmTensorField>(mesh_);
2476  printFieldInfo<surfaceTensorField>(mesh_);
2477  printIntFieldInfo<volScalarField::Internal>(mesh_);
2478  printIntFieldInfo<volVectorField::Internal>(mesh_);
2479  printIntFieldInfo<volSphericalTensorField::Internal>(mesh_);
2480  printIntFieldInfo<volSymmTensorField::Internal>(mesh_);
2481  printIntFieldInfo<volTensorField::Internal>(mesh_);
2482  Pout<< nl << endl;
2483  }
2484 
2485 
2486 
2487  // Receive and add what was sent
2488  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2489 
2490  // Disable parallel. Original state already known.
2491  UPstream::parRun(false);
2492 
2493  PtrList<labelList> domainSourceFaces(Pstream::nProcs());
2494  PtrList<labelList> domainSourceProcs(Pstream::nProcs());
2495  PtrList<labelList> domainSourcePatchs(Pstream::nProcs());
2496  PtrList<labelList> domainSourceNewNbrProcs(Pstream::nProcs());
2497  PtrList<labelList> domainSourcePointMasters(Pstream::nProcs());
2498 
2499  PtrList<fvMesh> domainMeshPtrs(Pstream::nProcs());
2500 
2501  // Define field storage lists
2502  #undef doLocalCode
2503  #define doLocalCode(FieldType, Variable) \
2504  PtrList<PtrList<FieldType>> Variable(UPstream::nProcs());
2505 
2506  doLocalCode(volScalarField, vsfs);
2507  doLocalCode(volVectorField, vvfs);
2508  doLocalCode(volSphericalTensorField, vsptfs);
2509  doLocalCode(volSymmTensorField, vsytfs);
2510  doLocalCode(volTensorField, vtfs);
2511 
2512  doLocalCode(surfaceScalarField, ssfs);
2513  doLocalCode(surfaceVectorField, svfs);
2514  doLocalCode(surfaceSphericalTensorField, ssptfs);
2515  doLocalCode(surfaceSymmTensorField, ssytfs);
2516  doLocalCode(surfaceTensorField, stfs);
2517 
2518  doLocalCode(volScalarField::Internal, dsfs);
2519  doLocalCode(volVectorField::Internal, dvfs);
2520  doLocalCode(volSphericalTensorField::Internal, dsptfs)
2521  doLocalCode(volSymmTensorField::Internal, dsytfs);
2522  doLocalCode(volTensorField::Internal, dtfs);
2523 
2524  #undef doLocalCode
2525 
2526 
2527  forAll(nRecvCells, sendProc)
2528  {
2529  // Did processor sendProc send anything to me?
2530  if (sendProc != UPstream::myProcNo() && nRecvCells[sendProc] > 0)
2531  {
2532  if (debug)
2533  {
2534  Pout<< nl
2535  << "RECEIVING FROM DOMAIN " << sendProc
2536  << " cells to receive:"
2537  << nRecvCells[sendProc]
2538  << nl << endl;
2539  }
2540 
2541 
2542  // Pstream for receiving mesh and fields
2543  UIPstream str(sendProc, pBufs);
2544 
2545 
2546  // Receive from sendProc - opposite of sendMesh
2547  {
2548  autoPtr<fvMesh> domainMeshPtr = receiveMesh
2549  (
2550  sendProc,
2551  pointZoneNames,
2552  faceZoneNames,
2553  cellZoneNames,
2554 
2555  const_cast<Time&>(mesh_.time()),
2556 
2557  domainSourceFaces.emplace_set(sendProc),
2558  domainSourceProcs.emplace_set(sendProc),
2559  domainSourcePatchs.emplace_set(sendProc),
2560  domainSourceNewNbrProcs.emplace_set(sendProc),
2561  domainSourcePointMasters.emplace_set(sendProc),
2562  str
2563  );
2564  domainMeshPtrs.set(sendProc, std::move(domainMeshPtr));
2565  fvMesh& domainMesh = domainMeshPtrs[sendProc];
2566  // Force construction of various on mesh.
2567  //(void)domainMesh.globalData();
2568 
2569 
2570  // Receive fields. Read as single dictionary because
2571  // of problems reading consecutive fields from single stream.
2572  dictionary fieldDicts(str);
2573 
2574  #undef doLocalCode
2575  #define doLocalCode(FieldType, Variable) \
2576  receiveFields<FieldType> \
2577  ( \
2578  sendProc, \
2579  allFieldNames, \
2580  domainMesh, \
2581  Variable.emplace_set(sendProc), \
2582  fieldDicts \
2583  )
2584 
2585  // Volume Fields
2586  doLocalCode(volScalarField, vsfs);
2587  doLocalCode(volVectorField, vvfs);
2588  doLocalCode(volSphericalTensorField, vsptfs);
2589  doLocalCode(volSymmTensorField, vsytfs);
2590  doLocalCode(volTensorField, vtfs);
2591 
2592  // Surface Fields
2593  doLocalCode(surfaceScalarField, ssfs);
2594  doLocalCode(surfaceVectorField, svfs);
2595  doLocalCode(surfaceSphericalTensorField, ssptfs);
2596  doLocalCode(surfaceSymmTensorField, ssytfs);
2597  doLocalCode(surfaceTensorField, stfs);
2598 
2599  // Dimensioned Fields
2600  doLocalCode(volScalarField::Internal, dsfs);
2601  doLocalCode(volVectorField::Internal, dvfs);
2602  doLocalCode(volSphericalTensorField::Internal, dsptfs);
2603  doLocalCode(volSymmTensorField::Internal, dsytfs);
2604  doLocalCode(volTensorField::Internal, dtfs);
2605 
2606  #undef doLocalCode
2607  }
2608  }
2609  }
2610 
2611  // Clear out storage
2612  pBufs.clear();
2613 
2614 
2615  // Set up pointers to meshes so we can include our mesh_
2616  UPtrList<polyMesh> meshes(domainMeshPtrs.size());
2617  UPtrList<fvMesh> fvMeshes(domainMeshPtrs.size());
2618  forAll(domainMeshPtrs, proci)
2619  {
2620  if (domainMeshPtrs.set(proci))
2621  {
2622  meshes.set(proci, &domainMeshPtrs[proci]);
2623  fvMeshes.set(proci, &domainMeshPtrs[proci]);
2624  }
2625  }
2626 
2627  // 'Receive' from myself
2628  {
2629  meshes.set(UPstream::myProcNo(), &mesh_);
2630  fvMeshes.set(UPstream::myProcNo(), &mesh_);
2631 
2632  domainSourceFaces.emplace_set(UPstream::myProcNo()) = sourceFace;
2633  domainSourceProcs.emplace_set(UPstream::myProcNo()) = sourceProc;
2634  domainSourcePatchs.emplace_set(UPstream::myProcNo()) = sourcePatch;
2635 
2636  domainSourceNewNbrProcs.emplace_set(UPstream::myProcNo()) =
2637  sourceNewNbrProc;
2638 
2639  domainSourcePointMasters.emplace_set(UPstream::myProcNo()) =
2640  sourcePointMaster;
2641  }
2642 
2643 
2644  // Find matching faces that need to be stitched
2645  labelListList localBoundaryFace(UPstream::nProcs());
2646  labelListList remoteFaceProc(UPstream::nProcs());
2647  labelListList remoteBoundaryFace(UPstream::nProcs());
2648  findCouples
2649  (
2650  meshes,
2651  domainSourceFaces,
2652  domainSourceProcs,
2653  domainSourcePatchs,
2654 
2655  localBoundaryFace,
2656  remoteFaceProc,
2657  remoteBoundaryFace
2658  );
2659 
2660 
2661  const label nOldInternalFaces = mesh_.nInternalFaces();
2662  const labelList oldFaceOwner(mesh_.faceOwner());
2663 
2664  // TBD: temporarily unset mesh moving to avoid problems in meshflux
2665  // mapping. To be investigated.
2666  const bool oldMoving = mesh_.moving(false);
2667 
2668  fvMeshAdder::add
2669  (
2670  Pstream::myProcNo(), // index of mesh to modify (== mesh_)
2671  fvMeshes,
2672  oldFaceOwner,
2673 
2674  // Coupling info
2675  localBoundaryFace,
2676  remoteFaceProc,
2677  remoteBoundaryFace,
2678 
2679  constructPatchMap,
2680  constructCellMap,
2681  constructFaceMap,
2682  constructPointMap
2683  );
2684 
2685  mesh_.moving(oldMoving);
2686 
2687 
2688  if (debug)
2689  {
2690  Pout<< nl << "ADDED REMOTE MESHES:" << endl;
2691  printMeshInfo(mesh_);
2692  printFieldInfo<volScalarField>(mesh_);
2693  printFieldInfo<volVectorField>(mesh_);
2694  printFieldInfo<volSphericalTensorField>(mesh_);
2695  printFieldInfo<volSymmTensorField>(mesh_);
2696  printFieldInfo<volTensorField>(mesh_);
2697  printFieldInfo<surfaceScalarField>(mesh_);
2698  printFieldInfo<surfaceVectorField>(mesh_);
2699  printFieldInfo<surfaceSphericalTensorField>(mesh_);
2700  printFieldInfo<surfaceSymmTensorField>(mesh_);
2701  printFieldInfo<surfaceTensorField>(mesh_);
2702  printIntFieldInfo<volScalarField::Internal>(mesh_);
2703  printIntFieldInfo<volVectorField::Internal>(mesh_);
2704  printIntFieldInfo<volSphericalTensorField::Internal>(mesh_);
2705  printIntFieldInfo<volSymmTensorField::Internal>(mesh_);
2706  printIntFieldInfo<volTensorField::Internal>(mesh_);
2707  Pout<< nl << endl;
2708  }
2709 
2710  {
2711  //- Combine sourceProc, sourcePatch, sourceFace
2712  sourceProc.resize_nocopy(mesh_.nBoundaryFaces());
2713  sourceProc = -1;
2714  sourcePatch.resize_nocopy(mesh_.nBoundaryFaces());
2715  sourcePatch = -1;
2716  sourceFace.resize_nocopy(mesh_.nBoundaryFaces());
2717  sourceFace = -1;
2718  sourceNewNbrProc.resize_nocopy(mesh_.nBoundaryFaces());
2719  sourceNewNbrProc = -1;
2720  sourcePointMaster.resize_nocopy(mesh_.nPoints());
2721  sourcePointMaster = -1;
2722 
2723  if (mesh_.nPoints() > 0)
2724  {
2725  forAll(meshes, meshi)
2726  {
2727  if (domainSourceFaces.set(meshi))
2728  {
2729  const label nIntFaces =
2730  (
2731  meshi == Pstream::myProcNo()
2732  ? nOldInternalFaces
2733  : meshes[meshi].nInternalFaces()
2734  );
2735  const labelList& faceOwner
2736  (
2737  meshi == Pstream::myProcNo()
2738  ? oldFaceOwner
2739  : meshes[meshi].faceOwner()
2740  );
2741 
2742  labelList& faceMap = constructFaceMap[meshi];
2743  const labelList& cellMap = constructCellMap[meshi];
2744 
2745  const labelList& domainSourceFace =
2746  domainSourceFaces[meshi];
2747  const labelList& domainSourceProc =
2748  domainSourceProcs[meshi];
2749  const labelList& domainSourcePatch =
2750  domainSourcePatchs[meshi];
2751  const labelList& domainSourceNewNbr =
2752  domainSourceNewNbrProcs[meshi];
2753  UIndirectList<label>
2754  (
2755  sourcePointMaster,
2756  constructPointMap[meshi]
2757  ) = domainSourcePointMasters[meshi];
2758 
2759 
2760  forAll(domainSourceFace, bFacei)
2761  {
2762  const label oldFacei = bFacei+nIntFaces;
2763  const label allFacei = faceMap[oldFacei];
2764  const label allbFacei = allFacei-mesh_.nInternalFaces();
2765 
2766  if (allbFacei >= 0)
2767  {
2768  sourceProc[allbFacei] = domainSourceProc[bFacei];
2769  sourcePatch[allbFacei] = domainSourcePatch[bFacei];
2770  sourceFace[allbFacei] = domainSourceFace[bFacei];
2771  sourceNewNbrProc[allbFacei] =
2772  domainSourceNewNbr[bFacei];
2773  }
2774  }
2775 
2776 
2777  // Add flip to constructFaceMap
2778  forAll(faceMap, oldFacei)
2779  {
2780  const label allFacei = faceMap[oldFacei];
2781  const label allOwn = mesh_.faceOwner()[allFacei];
2782 
2783  if (cellMap[faceOwner[oldFacei]] == allOwn)
2784  {
2785  // Master face
2786  faceMap[oldFacei] += 1;
2787  }
2788  else
2789  {
2790  // Slave face. Flip.
2791  faceMap[oldFacei] = -faceMap[oldFacei] - 1;
2792  }
2793  }
2794  }
2795  }
2796  }
2797  }
2798 
2799 
2800  UPstream::parRun(oldParRun); // Restore parallel state
2801 
2802 
2803  // Print a bit.
2804  if (debug)
2805  {
2806  Pout<< nl << "REDISTRIBUTED MESH:" << endl;
2807  printMeshInfo(mesh_);
2808  printFieldInfo<volScalarField>(mesh_);
2809  printFieldInfo<volVectorField>(mesh_);
2810  printFieldInfo<volSphericalTensorField>(mesh_);
2811  printFieldInfo<volSymmTensorField>(mesh_);
2812  printFieldInfo<volTensorField>(mesh_);
2813  printFieldInfo<surfaceScalarField>(mesh_);
2814  printFieldInfo<surfaceVectorField>(mesh_);
2815  printFieldInfo<surfaceSphericalTensorField>(mesh_);
2816  printFieldInfo<surfaceSymmTensorField>(mesh_);
2817  printFieldInfo<surfaceTensorField>(mesh_);
2818  printIntFieldInfo<volScalarField::Internal>(mesh_);
2819  printIntFieldInfo<volVectorField::Internal>(mesh_);
2820  printIntFieldInfo<volSphericalTensorField::Internal>(mesh_);
2821  printIntFieldInfo<volSymmTensorField::Internal>(mesh_);
2822  printIntFieldInfo<volTensorField::Internal>(mesh_);
2823  Pout<< nl << endl;
2824  }
2825 
2826 
2827  // See if any originally shared points need to be merged. Note: does
2828  // parallel comms. After this points and edges should again be consistent.
2829  mergeSharedPoints(sourcePointMaster, constructPointMap);
2830 
2831 
2832  // Add processorPatches
2833  // ~~~~~~~~~~~~~~~~~~~~
2834 
2835  // Per neighbour processor, per originating patch, the patchID
2836  // For faces resulting from internal faces or normal processor patches
2837  // the originating patch is -1. For cyclics this is the cyclic patchID.
2838  List<Map<label>> procPatchID;
2839 
2840  // Add processor and processorCyclic patches.
2841  addProcPatches(sourceNewNbrProc, sourcePatch, procPatchID);
2842 
2843  // Put faces into correct patch. Note that we now have proper
2844  // processorPolyPatches again so repatching will take care of coupled face
2845  // ordering.
2846 
2847  // Get boundary faces to be repatched. Is -1 or new patchID
2848  labelList newPatchID
2849  (
2850  getBoundaryPatch
2851  (
2852  sourceNewNbrProc,
2853  sourcePatch,
2854  procPatchID
2855  )
2856  );
2857 
2858  // Change patches. Since this might change ordering of coupled faces
2859  // we also need to adapt our constructMaps.
2860  repatch(newPatchID, constructFaceMap);
2861 
2862  // Bit of hack: processorFvPatchField does not get reset since created
2863  // from nothing so explicitly reset.
2864  initPatchFields<volScalarField, processorFvPatchField<scalar>>
2865  (
2866  Zero
2867  );
2868  initPatchFields<volVectorField, processorFvPatchField<vector>>
2869  (
2870  Zero
2871  );
2872  initPatchFields
2873  <
2874  volSphericalTensorField,
2875  processorFvPatchField<sphericalTensor>
2876  >
2877  (
2878  Zero
2879  );
2880  initPatchFields<volSymmTensorField, processorFvPatchField<symmTensor>>
2881  (
2882  Zero
2883  );
2884  initPatchFields<volTensorField, processorFvPatchField<tensor>>
2885  (
2886  Zero
2887  );
2888 
2889 
2890  mesh_.setInstance(mesh_.time().timeName());
2891 
2892 
2893  // Print a bit
2894  if (debug)
2895  {
2896  Pout<< nl << "FINAL MESH:" << endl;
2897  printMeshInfo(mesh_);
2898  printFieldInfo<volScalarField>(mesh_);
2899  printFieldInfo<volVectorField>(mesh_);
2900  printFieldInfo<volSphericalTensorField>(mesh_);
2901  printFieldInfo<volSymmTensorField>(mesh_);
2902  printFieldInfo<volTensorField>(mesh_);
2903  printFieldInfo<surfaceScalarField>(mesh_);
2904  printFieldInfo<surfaceVectorField>(mesh_);
2905  printFieldInfo<surfaceSphericalTensorField>(mesh_);
2906  printFieldInfo<surfaceSymmTensorField>(mesh_);
2907  printFieldInfo<surfaceTensorField>(mesh_);
2908  printIntFieldInfo<volScalarField::Internal>(mesh_);
2909  printIntFieldInfo<volVectorField::Internal>(mesh_);
2910  printIntFieldInfo<volSphericalTensorField::Internal>(mesh_);
2911  printIntFieldInfo<volSymmTensorField::Internal>(mesh_);
2912  printIntFieldInfo<volTensorField::Internal>(mesh_);
2913  Pout<< nl << endl;
2914  }
2915 
2916  // Collect all maps and return
2917  return autoPtr<mapDistributePolyMesh>::New
2918  (
2919  mesh_,
2920 
2921  nOldPoints,
2922  nOldFaces,
2923  nOldCells,
2924  std::move(oldPatchStarts),
2925  std::move(oldPatchNMeshPoints),
2926 
2927  std::move(subPointMap),
2928  std::move(subFaceMap),
2929  std::move(subCellMap),
2930  std::move(subPatchMap),
2931 
2932  std::move(constructPointMap),
2933  std::move(constructFaceMap),
2934  std::move(constructCellMap),
2935  std::move(constructPatchMap),
2936 
2937  true, // subFaceMap has flip
2938  true // constructFaceMap has flip
2939  );
2940 }
2941 
2942 
2943 // ************************************************************************* //
surfaceFields.H
Foam::surfaceFields.
Foam::sign
dimensionedScalar sign(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:159
Foam::fvsPatchVectorField
fvsPatchField< vector > fvsPatchVectorField
Definition: fvsPatchFieldsFwd.H:42
patchIDs
const labelList patchIDs(pbm.indices(polyPatchNames, true))
Foam::fvsPatchScalarField
fvsPatchField< scalar > fvsPatchScalarField
Definition: fvsPatchFieldsFwd.H:39
Foam::volSymmTensorField
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:85
dict
dictionary dict
Definition: searchingEngine.H:11
Foam::fvMeshTools::reorderPatches
static void reorderPatches(fvMesh &, const labelList &oldToNew, const label nPatches, const bool validBoundary)
Reorder and remove trailing patches.
Definition: fvMeshTools.C:317
Foam::faceZoneMesh
ZoneMesh< faceZone, polyMesh > faceZoneMesh
A ZoneMesh with the type faceZone.
Definition: faceZoneMeshFwd.H:38
Foam::fvMesh::Sf
const surfaceVectorField & Sf() const
Return cell face area vectors.
Definition: fvMeshGeometry.C:300
Foam::sum
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &f1)
Definition: DimensionedFieldFunctions.C:295
Foam::PtrListOps::names
List< word > names(const UPtrList< T > &list, const UnaryMatchPredicate &matcher)
List of names generated by calling name() for each list item and filtered for matches.
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:125
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::primitiveMesh::nPoints
label nPoints() const noexcept
Number of mesh points.
Definition: primitiveMeshI.H:30
Foam::FatalError
error FatalError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL ERROR&#39; header text and sta...
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:598
Foam::polyMesh::faceNeighbour
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1122
Foam::sortedOrder
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
Foam::IOobject::name
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:195
Foam::tmp::ref
T & ref() const
Return non-const reference to the contents of a non-null managed pointer.
Definition: tmpI.H:235
Foam::volTensorField
GeometricField< tensor, fvPatchField, volMesh > volTensorField
Definition: volFieldsFwd.H:86
Foam::renumber
IntListType renumber(const labelUList &oldToNew, const IntListType &input)
Renumber the values (not the indices) of a list.
Definition: ListOpsTemplates.C:26
Foam::BitOps::select
List< bool > select(const label n, const labelUList &locations)
Construct a selection list of bools (all false) with the given pre-size, subsequently add specified l...
Definition: BitOps.C:134
Foam::nl
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:50
Foam::inplaceReorder
void inplaceReorder(const labelUList &oldToNew, ListType &input, const bool prune=false)
Inplace reorder the elements of a list.
Definition: ListOpsTemplates.C:113
Foam::polyBoundaryMesh::start
label start() const noexcept
The start label of boundary faces in the polyMesh face list.
Definition: polyBoundaryMesh.C:708
Foam::labelUIndList
UIndirectList< label > labelUIndList
UIndirectList of labels.
Definition: IndirectList.H:65
Foam::surfaceVectorField
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
Definition: surfaceFieldsFwd.H:65
Foam::surfaceSphericalTensorField
GeometricField< sphericalTensor, fvsPatchField, surfaceMesh > surfaceSphericalTensorField
Definition: surfaceFieldsFwd.H:67
runTime
engineTime & runTime
Definition: createEngineTime.H:13
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:531
Foam::UPstream::parRun
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:1049
Foam::lessProcPatches
Less function class that can be used for sorting processor patches.
Definition: fvMeshDistribute.C:54
Foam::boundBox
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:63
Foam::CompactListList::pack
static CompactListList< T > pack(const UList< SubListType > &lists, const bool checkOverflow=false)
Construct by packing together the list of lists.
Foam::GeometricField< scalar, fvsPatchField, surfaceMesh >::Boundary
GeometricBoundaryField< scalar, fvsPatchField, surfaceMesh > Boundary
Type of boundary fields.
Definition: GeometricField.H:99
Foam::IOobjectOption::NO_WRITE
Ignore writing from objectRegistry::writeObject()
Definition: IOobjectOption.H:86
Foam::normalised
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Definition: quaternionI.H:674
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:182
Foam::UPstream::myProcNo
static int myProcNo(const label communicator=worldComm)
Rank of this process in the communicator (starting from masterNo()). Can be negative if the process i...
Definition: UPstream.H:1074
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:360
Foam::primitiveMesh::nFaces
label nFaces() const noexcept
Number of mesh faces.
Definition: primitiveMeshI.H:83
Foam::indirectPrimitivePatch
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
A PrimitivePatch with an IndirectList for the faces, const reference for the point field...
Definition: indirectPrimitivePatch.H:51
Foam::volVectorField
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:82
Foam::OSstream::name
virtual const fileName & name() const override
Get the name of the output serial stream. (eg, the name of the Fstream file name) ...
Definition: OSstream.H:128
Foam::labelListList
List< labelList > labelListList
List of labelList.
Definition: labelList.H:38
Foam::GeometricField< scalar, fvPatchField, volMesh >::Internal
DimensionedField< scalar, volMesh > Internal
The internal field type from which this GeometricField is derived.
Definition: GeometricField.H:94
doLocalCode
#define doLocalCode(FieldType, Variable)
Foam::volSphericalTensorField
GeometricField< sphericalTensor, fvPatchField, volMesh > volSphericalTensorField
Definition: volFieldsFwd.H:84
Foam::labelUList
UList< label > labelUList
A UList of labels.
Definition: UList.H:78
Foam::faceMap
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
Definition: blockMeshMergeTopological.C:87
Foam::Pstream::broadcast
static void broadcast(Type &value, const label comm=UPstream::worldComm)
Broadcast content (contiguous or non-contiguous) to all communicator ranks. Does nothing in non-paral...
Definition: PstreamBroadcast.C:28
ListOps.H
Various functions to operate on Lists.
Foam::polyMesh::points
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1078
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:421
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:81
Foam::labelHashSet
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
Definition: HashSet.H:85
Foam::FieldBase::typeName
static const char *const typeName
Typename for Field.
Definition: Field.H:86
Foam::faceList
List< face > faceList
List of faces.
Definition: faceListFwd.H:39
Foam::HashTableOps::values
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
Definition: HashOps.H:164
Foam::polyPatchList
PtrList< polyPatch > polyPatchList
Store lists of polyPatch as a PtrList.
Definition: polyPatch.H:56
Foam::HashTable::size
label size() const noexcept
The number of elements in table.
Definition: HashTable.H:342
Foam::UPstream::nProcs
static label nProcs(const label communicator=worldComm)
Number of ranks in parallel run (for given communicator). It is 1 for serial run. ...
Definition: UPstream.H:1065
Foam::fvMeshDistribute
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
Definition: fvMeshDistribute.H:66
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:38
Foam::Pstream::gatherList
static void gatherList(const List< commsStruct > &comms, List< T > &values, const int tag, const label comm)
Gather data, but keep individual values separate. Uses the specified communication schedule...
Definition: PstreamGatherList.C:41
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::cos
dimensionedScalar cos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:258
Foam::surfaceSymmTensorField
GeometricField< symmTensor, fvsPatchField, surfaceMesh > surfaceSymmTensorField
Definition: surfaceFieldsFwd.H:69
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
Definition: exprTraits.C:127
Foam::lessProcPatches::operator()
bool operator()(const label a, const label b)
Definition: fvMeshDistribute.C:67
Foam::identity
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
Definition: labelLists.C:44
Foam::constant::physicoChemical::b
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
Foam::polyMeshAdder::mergePoints
static void mergePoints(const polyMesh &, const Map< label > &pointToMaster, polyTopoChange &meshMod)
Helper: Merge points.
Definition: polyMeshAdder.C:2147
Foam::polyMesh::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
Definition: polyMesh.H:608
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:46
Foam::fvMeshTools::addPatch
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:38
Foam::polyMesh::defaultRegion
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:405
Foam::polyBoundaryMesh::names
wordList names() const
Return a list of patch names.
Definition: polyBoundaryMesh.C:646
Foam::UPtrList::size
label size() const noexcept
The number of entries in the list.
Definition: UPtrListI.H:106
Foam::surfaceTensorField
GeometricField< tensor, fvsPatchField, surfaceMesh > surfaceTensorField
Definition: surfaceFieldsFwd.H:70
meshes
Foam::PtrList< Foam::fvMesh > meshes(regionNames.size())
Foam::fvMeshAdder::add
static autoPtr< mapAddedPolyMesh > add(fvMesh &mesh0, const fvMesh &mesh1, const faceCoupleInfo &coupleInfo, const bool validBoundary=true, const bool fullyMapped=false)
Inplace add mesh to fvMesh. Maps all stored fields. Returns map.
Definition: fvMeshAdder.C:67
Foam::polyMesh::faceOwner
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1116
Foam::primitiveMesh::nInternalFaces
label nInternalFaces() const noexcept
Number of internal faces.
Definition: primitiveMeshI.H:71
Foam::pointZoneMesh
ZoneMesh< pointZone, polyMesh > pointZoneMesh
A ZoneMesh with the type pointZone.
Definition: pointZoneMeshFwd.H:38
Foam::polyMesh::faces
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1103
Foam::distribution
Accumulating histogram of values. Specified bin resolution automatic generation of bins...
Definition: distribution.H:57
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:139
Foam::PtrList::set
const T * set(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
Definition: PtrList.H:159
labelPairHashes.H
A HashTable to objects of type <T> with a labelPair key. The hashing is based on labelPair (FixedList...
Foam::fvMesh::magSf
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
Definition: fvMeshGeometry.C:311
Foam::polyBoundaryMesh::checkParallelSync
bool checkParallelSync(const bool report=false) const
Check whether all procs have all patches and in same order.
Definition: polyBoundaryMesh.C:1148
Foam::UPstream::allToAll
static void allToAll(const UList< int32_t > &sendData, UList< int32_t > &recvData, const label communicator=worldComm)
Exchange int32_t data with all ranks in communicator.
Definition: UPstreamAllToAll.C:42
Foam::vector
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Return a time name for the given scalar time value formatted with the given precision.
Definition: Time.C:714
Foam::ensightOutput::debug
int debug
Static debugging option.
Foam::labelPair
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:51
Foam::polyMesh::faceZones
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
Definition: polyMesh.H:670
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::fvMeshDistribute::distribute
autoPtr< mapDistributePolyMesh > distribute(const labelList &dist)
Send cells to neighbours according to distribution.
Definition: fvMeshDistribute.C:1874
Foam::DimensionedField::typeName
const word typeName("volScalarField::Internal")
Foam::fvMeshDistribute::printMeshInfo
static void printMeshInfo(const polyMesh &)
Print some info on mesh.
Definition: fvMeshDistribute.C:220
fld
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputLagrangian.H:23
Foam::primitiveMesh::faceCentres
const vectorField & faceCentres() const
Definition: primitiveMeshFaceCentresAndAreas.C:71
Foam::UPtrList::empty
bool empty() const noexcept
True if the list is empty (ie, size() is zero)
Definition: UPtrListI.H:99
Foam::polyMesh::pointZones
const pointZoneMesh & pointZones() const noexcept
Return point zone mesh.
Definition: polyMesh.H:662
Foam::wordList
List< word > wordList
List of word.
Definition: fileName.H:59
Foam::inplaceRenumber
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
Definition: ListOpsTemplates.C:50
Foam::fvMeshDistribute::generateTestField
static tmp< surfaceScalarField > generateTestField(const fvMesh &)
Generate a test field on faces.
Definition: fvMeshDistribute.C:383
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:414
newPointi
label newPointi
Definition: readKivaGrid.H:496
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:33
Foam::primitiveMesh::nCells
label nCells() const noexcept
Number of mesh cells.
Definition: primitiveMeshI.H:89
Foam::glTF::key
auto key(const Type &t) -> typename std::enable_if< std::is_enum< Type >::value, typename std::underlying_type< Type >::type >::type
Definition: foamGltfBase.H:103
Foam::fvMeshDistribute::printCoupleInfo
static void printCoupleInfo(const primitiveMesh &, const labelList &, const labelList &, const labelList &, const labelList &)
Print some info on coupling data.
Definition: fvMeshDistribute.C:263
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::UPstream::master
static bool master(const label communicator=worldComm)
True if process corresponds to the master rank in the communicator.
Definition: UPstream.H:1082
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:54
Foam::IOobjectOption::NO_READ
Nothing to be read.
Definition: IOobjectOption.H:60
Foam::lessProcPatches::lessProcPatches
lessProcPatches(const labelList &nbrProc, const labelList &referPatchID)
Definition: fvMeshDistribute.C:61
Foam::labelPairLookup
HashTable< label, labelPair, Foam::Hash< labelPair > > labelPairLookup
This is a Map of a labelPair to a label. Used for e.g. for face1, face2 to shared edge...
Definition: labelPairHashes.H:63
Foam::polyMesh::cellZones
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
Definition: polyMesh.H:678
Foam::reduce
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
Reduce inplace (cf. MPI Allreduce) using specified communication schedule.
Definition: PstreamReduceOps.H:49
Foam::UPstream::commsTypes::nonBlocking
"nonBlocking" : (MPI_Isend, MPI_Irecv)
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::polyPatch::New
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return pointer to a new patch created on freestore from components.
Definition: polyPatchNew.C:28
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
Definition: HashPtrTable.H:48
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
forAllReverse
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Definition: stdFoam.H:437
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:62
Foam::tmp
A class for managing temporary objects.
Definition: HashPtrTable.H:50
Foam::surfaceScalarField
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Definition: surfaceFieldsFwd.H:64
Foam::autoPtr::New
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
Definition: autoPtr.H:178
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:69
Foam::primitiveMesh::nBoundaryFaces
label nBoundaryFaces() const noexcept
Number of boundary faces (== nFaces - nInternalFaces)
Definition: primitiveMeshI.H:77
Foam::fvMeshDistribute::countCells
static labelList countCells(const labelList &)
Helper function: count cells per processor in wanted distribution.
Definition: fvMeshDistribute.C:1850
Foam::fvMeshDistribute::testField
static void testField(const surfaceScalarField &)
Check whether field consistent with face orientation.
Definition: fvMeshDistribute.C:433
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported...
Definition: IOobject.H:172
found
bool found
Definition: TABSMDCalcMethod2.H:32
Foam::Pout
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
Foam::cellZoneMesh
ZoneMesh< cellZone, polyMesh > cellZoneMesh
A ZoneMesh with the type cellZone.
Definition: cellZoneMeshFwd.H:38
Foam::HashTable::set
bool set(const Key &key, const T &obj)
Copy assign a new entry, overwriting existing entries.
Definition: HashTableI.H:174
pp
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:26
Foam::syncTools::swapBoundaryFaceList
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &faceValues)
Swap coupled boundary face values. Uses eqOp.
Definition: syncTools.H:485
Foam::flatOutput
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
Definition: FlatOutput.H:225
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:127