redistributePar.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-2017 OpenFOAM Foundation
9  Copyright (C) 2015-2020 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 Application
28  redistributePar
29 
30 Group
31  grpParallelUtilities
32 
33 Description
34  Redistributes existing decomposed mesh and fields according to the current
35  settings in the decomposeParDict file.
36 
37  Must be run on maximum number of source and destination processors.
38  Balances mesh and writes new mesh to new time directory.
39 
40  Can optionally run in decompose/reconstruct mode to decompose/reconstruct
41  mesh and fields.
42 
43 Usage
44  \b redistributePar [OPTION]
45 
46  Options:
47  - \par -decompose
48  Remove any existing \a processor subdirectories and decomposes the
49  mesh. Equivalent to running without processor subdirectories.
50 
51  - \par -reconstruct
52  Reconstruct mesh and fields (like reconstructParMesh+reconstructPar).
53 
54  - \par -newTimes
55  (in combination with -reconstruct) reconstruct only new times.
56 
57  - \par -dry-run
58  (not in combination with -reconstruct) Test without actually
59  decomposing.
60 
61  - \par -cellDist
62  not in combination with -reconstruct) Write the cell distribution
63  as a labelList, for use with 'manual'
64  decomposition method and as a volScalarField for visualization.
65 
66  - \par -region <regionName>
67  Distribute named region.
68 
69  - \par -allRegions
70  Distribute all regions in regionProperties. Does not check for
71  existence of processor*.
72 
73 \*---------------------------------------------------------------------------*/
74 
75 #include "argList.H"
76 #include "sigFpe.H"
77 #include "Time.H"
78 #include "fvMesh.H"
79 #include "fvMeshTools.H"
80 #include "fvMeshDistribute.H"
81 #include "decompositionMethod.H"
82 #include "decompositionModel.H"
83 #include "timeSelector.H"
84 #include "PstreamReduceOps.H"
85 #include "volFields.H"
86 #include "surfaceFields.H"
88 #include "IOobjectList.H"
89 #include "globalIndex.H"
90 #include "loadOrCreateMesh.H"
91 #include "processorFvPatchField.H"
93 #include "topoSet.H"
94 #include "regionProperties.H"
95 #include "basicFvGeometryScheme.H"
96 
100 #include "hexRef8Data.H"
101 #include "meshRefinement.H"
102 #include "pointFields.H"
103 
104 #include "cyclicACMIFvPatch.H"
105 
106 using namespace Foam;
107 
108 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
109 
110 // Tolerance (as fraction of the bounding box). Needs to be fairly lax since
111 // usually meshes get written with limited precision (6 digits)
112 static const scalar defaultMergeTol = 1e-6;
113 
114 
115 // Get merging distance when matching face centres
116 scalar getMergeDistance
117 (
118  const argList& args,
119  const Time& runTime,
120  const boundBox& bb
121 )
122 {
123  const scalar mergeTol =
124  args.getOrDefault<scalar>("mergeTol", defaultMergeTol);
125 
126  const scalar writeTol =
127  Foam::pow(scalar(10), -scalar(IOstream::defaultPrecision()));
128 
129  Info<< "Merge tolerance : " << mergeTol << nl
130  << "Write tolerance : " << writeTol << endl;
131 
132  if (runTime.writeFormat() == IOstream::ASCII && mergeTol < writeTol)
133  {
135  << "Your current settings specify ASCII writing with "
136  << IOstream::defaultPrecision() << " digits precision." << nl
137  << "Your merging tolerance (" << mergeTol << ") is finer than this."
138  << nl
139  << "Please change your writeFormat to binary"
140  << " or increase the writePrecision" << endl
141  << "or adjust the merge tolerance (-mergeTol)."
142  << exit(FatalError);
143  }
144 
145  const scalar mergeDist = mergeTol * bb.mag();
146 
147  Info<< "Overall meshes bounding box : " << bb << nl
148  << "Relative tolerance : " << mergeTol << nl
149  << "Absolute matching distance : " << mergeDist << nl
150  << endl;
151 
152  return mergeDist;
153 }
154 
155 
156 void setBasicGeometry(fvMesh& mesh)
157 {
158  // Set the fvGeometryScheme to basic since it does not require
159  // any parallel communication to construct the geometry. During
160  // redistributePar one might temporarily end up with processors
161  // with zero procBoundaries. Normally procBoundaries trigger geometry
162  // calculation (e.g. send over cellCentres) so on the processors without
163  // procBoundaries this will not happen. The call to the geometry calculation
164  // is not synchronised and this might lead to a hang for geometry schemes
165  // that do require synchronisation
166 
167  tmp<fvGeometryScheme> basicGeometry
168  (
170  (
171  mesh,
172  dictionary(),
173  basicFvGeometryScheme::typeName
174  )
175  );
176  mesh.geometry(basicGeometry);
177 }
178 
179 
180 void printMeshData(const polyMesh& mesh)
181 {
182  // Collect all data on master
183 
184  globalIndex globalCells(mesh.nCells());
185  labelListList patchNeiProcNo(Pstream::nProcs());
186  labelListList patchSize(Pstream::nProcs());
187  const labelList& pPatches = mesh.globalData().processorPatches();
188  patchNeiProcNo[Pstream::myProcNo()].setSize(pPatches.size());
189  patchSize[Pstream::myProcNo()].setSize(pPatches.size());
190  forAll(pPatches, i)
191  {
192  const processorPolyPatch& ppp = refCast<const processorPolyPatch>
193  (
194  mesh.boundaryMesh()[pPatches[i]]
195  );
196  patchNeiProcNo[Pstream::myProcNo()][i] = ppp.neighbProcNo();
197  patchSize[Pstream::myProcNo()][i] = ppp.size();
198  }
199  Pstream::gatherList(patchNeiProcNo);
200  Pstream::gatherList(patchSize);
201 
202 
203  // Print stats
204 
205  globalIndex globalBoundaryFaces(mesh.nBoundaryFaces());
206 
207  label maxProcCells = 0;
208  label totProcFaces = 0;
209  label maxProcPatches = 0;
210  label totProcPatches = 0;
211  label maxProcFaces = 0;
212 
213  for (const int procI : Pstream::allProcs())
214  {
215  Info<< nl
216  << "Processor " << procI << nl
217  << " Number of cells = " << globalCells.localSize(procI)
218  << endl;
219 
220  label nProcFaces = 0;
221 
222  const labelList& nei = patchNeiProcNo[procI];
223 
224  forAll(patchNeiProcNo[procI], i)
225  {
226  Info<< " Number of faces shared with processor "
227  << patchNeiProcNo[procI][i] << " = " << patchSize[procI][i]
228  << endl;
229 
230  nProcFaces += patchSize[procI][i];
231  }
232 
233  Info<< " Number of processor patches = " << nei.size() << nl
234  << " Number of processor faces = " << nProcFaces << nl
235  << " Number of boundary faces = "
236  << globalBoundaryFaces.localSize(procI)-nProcFaces << endl;
237 
238  maxProcCells = max(maxProcCells, globalCells.localSize(procI));
239  totProcFaces += nProcFaces;
240  totProcPatches += nei.size();
241  maxProcPatches = max(maxProcPatches, nei.size());
242  maxProcFaces = max(maxProcFaces, nProcFaces);
243  }
244 
245  // Stats
246 
247  scalar avgProcCells = scalar(globalCells.size())/Pstream::nProcs();
248  scalar avgProcPatches = scalar(totProcPatches)/Pstream::nProcs();
249  scalar avgProcFaces = scalar(totProcFaces)/Pstream::nProcs();
250 
251  // In case of all faces on one processor. Just to avoid division by 0.
252  if (totProcPatches == 0)
253  {
254  avgProcPatches = 1;
255  }
256  if (totProcFaces == 0)
257  {
258  avgProcFaces = 1;
259  }
260 
261  Info<< nl
262  << "Number of processor faces = " << totProcFaces/2 << nl
263  << "Max number of cells = " << maxProcCells
264  << " (" << 100.0*(maxProcCells-avgProcCells)/avgProcCells
265  << "% above average " << avgProcCells << ")" << nl
266  << "Max number of processor patches = " << maxProcPatches
267  << " (" << 100.0*(maxProcPatches-avgProcPatches)/avgProcPatches
268  << "% above average " << avgProcPatches << ")" << nl
269  << "Max number of faces between processors = " << maxProcFaces
270  << " (" << 100.0*(maxProcFaces-avgProcFaces)/avgProcFaces
271  << "% above average " << avgProcFaces << ")" << nl
272  << endl;
273 }
274 
275 
276 // Debugging: write volScalarField with decomposition for post processing.
277 void writeDecomposition
278 (
279  const word& name,
280  const fvMesh& mesh,
281  const labelUList& decomp
282 )
283 {
284  // Write the decomposition as labelList for use with 'manual'
285  // decomposition method.
286  labelIOList cellDecomposition
287  (
288  IOobject
289  (
290  "cellDecomposition",
291  mesh.facesInstance(), // mesh read from facesInstance
292  mesh,
295  false
296  ),
297  decomp
298  );
299  cellDecomposition.write();
300 
301  Info<< "Writing wanted cell distribution to volScalarField " << name
302  << " for postprocessing purposes." << nl << endl;
303 
304  volScalarField procCells
305  (
306  IOobject
307  (
308  name,
309  mesh.time().timeName(),
310  mesh,
313  false // do not register
314  ),
315  mesh,
317  zeroGradientFvPatchScalarField::typeName
318  );
319 
320  forAll(procCells, cI)
321  {
322  procCells[cI] = decomp[cI];
323  }
324 
325  procCells.correctBoundaryConditions();
326  procCells.write();
327 }
328 
329 
330 void determineDecomposition
331 (
332  const Time& baseRunTime,
333  const fileName& decompDictFile, // optional location for decomposeParDict
334  const bool decompose, // decompose, i.e. read from undecomposed case
335  const fileName& proc0CaseName,
336  const fvMesh& mesh,
337  const bool writeCellDist,
338 
339  label& nDestProcs,
340  labelList& decomp
341 )
342 {
343  // Read decomposeParDict (on all processors)
345  (
346  mesh,
347  decompDictFile
348  );
349 
350  decompositionMethod& decomposer = method.decomposer();
351 
352  if (!decomposer.parallelAware())
353  {
355  << "You have selected decomposition method "
356  << decomposer.typeName
357  << " which does" << nl
358  << "not synchronise the decomposition across"
359  << " processor patches." << nl
360  << " You might want to select a decomposition method"
361  << " which is aware of this. Continuing."
362  << endl;
363  }
364 
365  if (Pstream::master() && decompose)
366  {
367  Info<< "Setting caseName to " << baseRunTime.caseName()
368  << " to read decomposeParDict" << endl;
369  const_cast<Time&>(mesh.time()).caseName() = baseRunTime.caseName();
370  }
371 
372  scalarField cellWeights;
373  if (method.found("weightField"))
374  {
375  word weightName = method.get<word>("weightField");
376 
377  volScalarField weights
378  (
379  IOobject
380  (
381  weightName,
382  mesh.time().timeName(),
383  mesh,
386  ),
387  mesh
388  );
389  cellWeights = weights.internalField();
390  }
391 
392  nDestProcs = decomposer.nDomains();
393  decomp = decomposer.decompose(mesh, cellWeights);
394 
395  if (Pstream::master() && decompose)
396  {
397  Info<< "Restoring caseName to " << proc0CaseName << endl;
398  const_cast<Time&>(mesh.time()).caseName() = proc0CaseName;
399  }
400 
401  // Dump decomposition to volScalarField
402  if (writeCellDist)
403  {
404  // Note: on master make sure to write to processor0
405  if (decompose)
406  {
407  if (Pstream::master())
408  {
409  Info<< "Setting caseName to " << baseRunTime.caseName()
410  << " to write undecomposed cellDist" << endl;
411 
412  Time& tm = const_cast<Time&>(mesh.time());
413 
414  tm.caseName() = baseRunTime.caseName();
415  writeDecomposition("cellDist", mesh, decomp);
416  Info<< "Restoring caseName to " << proc0CaseName << endl;
417  tm.caseName() = proc0CaseName;
418  }
419  }
420  else
421  {
422  writeDecomposition("cellDist", mesh, decomp);
423  }
424  }
425 }
426 
427 
428 // Write addressing if decomposing (1 to many) or reconstructing (many to 1)
429 void writeProcAddressing
430 (
431  const fvMesh& mesh,
432  const mapDistributePolyMesh& map,
433  const bool decompose
434 )
435 {
436  Info<< "Writing procAddressing files to " << mesh.facesInstance()
437  << endl;
438 
439  labelIOList cellMap
440  (
441  IOobject
442  (
443  "cellProcAddressing",
446  mesh,
448  ),
449  0
450  );
451 
453  (
454  IOobject
455  (
456  "faceProcAddressing",
459  mesh,
461  ),
462  0
463  );
464 
465  labelIOList pointMap
466  (
467  IOobject
468  (
469  "pointProcAddressing",
472  mesh,
474  ),
475  0
476  );
477 
478  labelIOList patchMap
479  (
480  IOobject
481  (
482  "boundaryProcAddressing",
485  mesh,
487  ),
488  0
489  );
490 
491  // Decomposing: see how cells moved from undecomposed case
492  if (decompose)
493  {
494  cellMap = identity(map.nOldCells());
495  map.distributeCellData(cellMap);
496 
497  faceMap = identity(map.nOldFaces());
498  {
499  const mapDistribute& faceDistMap = map.faceMap();
500 
501  if (faceDistMap.subHasFlip() || faceDistMap.constructHasFlip())
502  {
503  // Offset by 1
504  faceMap = faceMap + 1;
505  }
506  // Apply face flips
508  (
510  List<labelPair>(),
511  faceDistMap.constructSize(),
512  faceDistMap.subMap(),
513  faceDistMap.subHasFlip(),
514  faceDistMap.constructMap(),
515  faceDistMap.constructHasFlip(),
516  faceMap,
517  flipLabelOp()
518  );
519  }
520 
521  pointMap = identity(map.nOldPoints());
522  map.distributePointData(pointMap);
523 
524  patchMap = identity(map.oldPatchSizes().size());
525  const mapDistribute& patchDistMap = map.patchMap();
526  // Use explicit distribute since we need to provide a null value
527  // (for new patches) and this is the only call that allow us to
528  // provide one ...
530  (
532  List<labelPair>(),
533  patchDistMap.constructSize(),
534  patchDistMap.subMap(),
535  patchDistMap.subHasFlip(),
536  patchDistMap.constructMap(),
537  patchDistMap.constructHasFlip(),
538  patchMap,
539  label(-1),
540  eqOp<label>(),
541  flipOp(),
543  );
544  }
545  else // reconstruct
546  {
547  cellMap = identity(mesh.nCells());
548  map.cellMap().reverseDistribute(map.nOldCells(), cellMap);
549 
551  {
552  const mapDistribute& faceDistMap = map.faceMap();
553 
554  if (faceDistMap.subHasFlip() || faceDistMap.constructHasFlip())
555  {
556  // Offset by 1
557  faceMap = faceMap + 1;
558  }
559 
561  (
563  List<labelPair>(),
564  map.nOldFaces(),
565  faceDistMap.constructMap(),
566  faceDistMap.constructHasFlip(),
567  faceDistMap.subMap(),
568  faceDistMap.subHasFlip(),
569  faceMap,
570  flipLabelOp()
571  );
572  }
573 
574  pointMap = identity(mesh.nPoints());
575  map.pointMap().reverseDistribute(map.nOldPoints(), pointMap);
576 
577  const mapDistribute& patchDistMap = map.patchMap();
578  patchMap = identity(mesh.boundaryMesh().size());
579  patchDistMap.reverseDistribute
580  (
581  map.oldPatchSizes().size(),
582  label(-1),
583  patchMap
584  );
585  }
586 
587  const bool cellOk = cellMap.write();
588  const bool faceOk = faceMap.write();
589  const bool pointOk = pointMap.write();
590  const bool patchOk = patchMap.write();
591 
592  if (!cellOk || !faceOk || !pointOk || !patchOk)
593  {
595  << "Failed to write " << cellMap.objectPath()
596  << ", " << faceMap.objectPath()
597  << ", " << pointMap.objectPath()
598  << ", " << patchMap.objectPath()
599  << endl;
600  }
601 }
602 
603 
604 // Remove addressing
605 void removeProcAddressing(const polyMesh& mesh)
606 {
607  for (const auto prefix : {"boundary", "cell", "face", "point"})
608  {
609  IOobject io
610  (
611  prefix + word("ProcAddressing"),
614  mesh
615  );
616 
617  const fileName procFile(io.objectPath());
618  rm(procFile);
619  }
620 }
621 
622 
623 // Generic mesh-based field reading
624 template<class GeoField>
625 void readField
626 (
627  const IOobject& io,
628  const fvMesh& mesh,
629  const label i,
631 )
632 {
633  fields.set(i, new GeoField(io, mesh));
634 }
635 
636 
637 // Definition of readField for GeometricFields only
638 template<class Type, template<class> class PatchField, class GeoMesh>
639 void readField
640 (
641  const IOobject& io,
642  const fvMesh& mesh,
643  const label i,
645 )
646 {
647  fields.set
648  (
649  i,
651  );
652 }
653 
654 
655 // Read vol or surface fields
656 template<class GeoField>
657 void readFields
658 (
659  const boolList& haveMesh,
660  const fvMesh& mesh,
661  const autoPtr<fvMeshSubset>& subsetterPtr,
662  IOobjectList& allObjects,
664 )
665 {
666  // Get my objects of type
667  IOobjectList objects(allObjects.lookupClass(GeoField::typeName));
668 
669  // Check that we all have all objects
670  wordList objectNames = objects.sortedNames();
671 
672  // Get master names
673  wordList masterNames(objectNames);
674  Pstream::scatter(masterNames);
675 
676  if (haveMesh[Pstream::myProcNo()] && objectNames != masterNames)
677  {
679  << "Objects not synchronised across processors." << nl
680  << "Master has " << flatOutput(masterNames) << nl
681  << "Processor " << Pstream::myProcNo()
682  << " has " << flatOutput(objectNames)
683  << exit(FatalError);
684  }
685 
686  fields.setSize(masterNames.size());
687 
688  // Have master send all fields to processors that don't have a mesh
689  if (Pstream::master())
690  {
691  forAll(masterNames, i)
692  {
693  const word& name = masterNames[i];
694  IOobject& io = *objects[name];
696 
697  // Load field (but not oldTime)
698  readField(io, mesh, i, fields);
699  // Create zero sized field and send
700  if (subsetterPtr)
701  {
702  tmp<GeoField> tsubfld = subsetterPtr->interpolate(fields[i]);
703 
704  // Send to all processors that don't have a mesh
705  for (const int procI : Pstream::subProcs())
706  {
707  if (!haveMesh[procI])
708  {
710  toProc<< tsubfld();
711  }
712  }
713  }
714  }
715  }
716  else if (!haveMesh[Pstream::myProcNo()])
717  {
718  // Don't have mesh (nor fields). Receive empty field from master.
719 
720  forAll(masterNames, i)
721  {
722  const word& name = masterNames[i];
723 
724  // Receive field
725  IPstream fromMaster
726  (
729  );
730  dictionary fieldDict(fromMaster);
731 
732  fields.set
733  (
734  i,
735  new GeoField
736  (
737  IOobject
738  (
739  name,
740  mesh.time().timeName(),
741  mesh,
744  ),
745  mesh,
746  fieldDict
747  )
748  );
749 
751  //fields[i].write();
752  }
753  }
754  else
755  {
756  // Have mesh so just try to load
757  forAll(masterNames, i)
758  {
759  const word& name = masterNames[i];
760  IOobject& io = *objects[name];
762 
763  // Load field (but not oldtime)
764  readField(io, mesh, i, fields);
765  }
766  }
767 }
768 
769 
770 // Variant of GeometricField::correctBoundaryConditions that only
771 // evaluates selected patch fields
772 template<class GeoField, class CoupledPatchType>
773 void correctCoupledBoundaryConditions(fvMesh& mesh)
774 {
776  (
777  mesh.objectRegistry::lookupClass<GeoField>()
778  );
779 
780  forAllIters(flds, iter)
781  {
782  GeoField& fld = *iter();
783 
784  typename GeoField::Boundary& bfld = fld.boundaryFieldRef();
785  if
786  (
789  )
790  {
791  const label nReq = Pstream::nRequests();
792 
793  forAll(bfld, patchi)
794  {
795  auto& pfld = bfld[patchi];
796  const auto& fvp = mesh.boundary()[patchi];
797 
798  if (fvp.coupled() && !isA<cyclicACMIFvPatch>(fvp))
799  {
800  pfld.initEvaluate(Pstream::defaultCommsType);
801  }
802  }
803 
804  // Block for any outstanding requests
805  if
806  (
809  )
810  {
811  Pstream::waitRequests(nReq);
812  }
813 
814  for (auto& pfld : bfld)
815  {
816  if (pfld.patch().coupled())
817  {
818  pfld.evaluate(Pstream::defaultCommsType);
819  }
820  }
821  }
823  {
824  const lduSchedule& patchSchedule =
825  fld.mesh().globalData().patchSchedule();
826 
827  forAll(patchSchedule, patchEvali)
828  {
829  const label patchi = patchSchedule[patchEvali].patch;
830  const auto& fvp = mesh.boundary()[patchi];
831  auto& pfld = bfld[patchi];
832 
833  if (fvp.coupled() && !isA<cyclicACMIFvPatch>(fvp))
834  {
835  if (patchSchedule[patchEvali].init)
836  {
837  pfld.initEvaluate(Pstream::commsTypes::scheduled);
838  }
839  else
840  {
841  pfld.evaluate(Pstream::commsTypes::scheduled);
842  }
843  }
844  }
845  }
846  else
847  {
849  << "Unsupported communications type "
851  << exit(FatalError);
852  }
853  }
854 }
855 
856 
857 // Inplace redistribute mesh and any fields
858 autoPtr<mapDistributePolyMesh> redistributeAndWrite
859 (
860  const Time& baseRunTime,
861  const scalar tolDim,
862  const boolList& haveMesh,
863  const fileName& meshSubDir,
864  const bool doReadFields,
865  const bool decompose, // decompose, i.e. read from undecomposed case
866  const bool reconstruct,
867  const bool overwrite,
868  const fileName& proc0CaseName,
869  const label nDestProcs,
870  const labelList& decomp,
871  const fileName& masterInstDir,
872  fvMesh& mesh
873 )
874 {
875  Time& runTime = const_cast<Time&>(mesh.time());
876 
878  //Info<< "Before distribution:" << endl;
879  //printMeshData(mesh);
880 
881 
882  PtrList<volScalarField> volScalarFields;
883  PtrList<volVectorField> volVectorFields;
884  PtrList<volSphericalTensorField> volSphereTensorFields;
885  PtrList<volSymmTensorField> volSymmTensorFields;
886  PtrList<volTensorField> volTensorFields;
887 
888  PtrList<surfaceScalarField> surfScalarFields;
889  PtrList<surfaceVectorField> surfVectorFields;
890  PtrList<surfaceSphericalTensorField> surfSphereTensorFields;
891  PtrList<surfaceSymmTensorField> surfSymmTensorFields;
892  PtrList<surfaceTensorField> surfTensorFields;
893 
899 
900  DynamicList<word> pointFieldNames;
901 
902 
903  if (doReadFields)
904  {
905  // Create 0 sized mesh to do all the generation of zero sized
906  // fields on processors that have zero sized meshes. Note that this is
907  // only necessary on master but since polyMesh construction with
908  // Pstream::parRun does parallel comms we have to do it on all
909  // processors
910  autoPtr<fvMeshSubset> subsetterPtr;
911 
912  const bool allHaveMesh = !haveMesh.found(false);
913  if (!allHaveMesh)
914  {
915  // Find last non-processor patch.
917 
918  const label nonProcI = (patches.nNonProcessor() - 1);
919 
920  if (nonProcI < 0)
921  {
923  << "Cannot find non-processor patch on processor "
924  << Pstream::myProcNo() << nl
925  << " Current patches:" << patches.names()
926  << abort(FatalError);
927  }
928 
929  // Subset 0 cells, no parallel comms.
930  // This is used to create zero-sized fields.
931  subsetterPtr.reset
932  (
933  new fvMeshSubset(mesh, bitSet(), nonProcI, false)
934  );
935  }
936 
937 
938  // Get original objects (before incrementing time!)
939  if (Pstream::master() && decompose)
940  {
941  runTime.caseName() = baseRunTime.caseName();
942  }
943  IOobjectList objects(mesh, runTime.timeName());
944  if (Pstream::master() && decompose)
945  {
946  runTime.caseName() = proc0CaseName;
947  }
948 
949  Info<< "From time " << runTime.timeName()
950  << " have objects:" << objects.names() << endl;
951 
952  // We don't want to map the decomposition (mapping already tested when
953  // mapping the cell centre field)
954  auto iter = objects.find("cellDist");
955  if (iter.found())
956  {
957  objects.erase(iter);
958  }
959 
960 
961  // volFields
962 
963  if (Pstream::master() && decompose)
964  {
965  runTime.caseName() = baseRunTime.caseName();
966  }
967  readFields
968  (
969  haveMesh,
970  mesh,
971  subsetterPtr,
972  objects,
973  volScalarFields
974  );
975 
976  readFields
977  (
978  haveMesh,
979  mesh,
980  subsetterPtr,
981  objects,
982  volVectorFields
983  );
984 
985  readFields
986  (
987  haveMesh,
988  mesh,
989  subsetterPtr,
990  objects,
991  volSphereTensorFields
992  );
993 
994  readFields
995  (
996  haveMesh,
997  mesh,
998  subsetterPtr,
999  objects,
1000  volSymmTensorFields
1001  );
1002 
1003  readFields
1004  (
1005  haveMesh,
1006  mesh,
1007  subsetterPtr,
1008  objects,
1009  volTensorFields
1010  );
1011 
1012 
1013  // surfaceFields
1014 
1015  readFields
1016  (
1017  haveMesh,
1018  mesh,
1019  subsetterPtr,
1020  objects,
1021  surfScalarFields
1022  );
1023 
1024  readFields
1025  (
1026  haveMesh,
1027  mesh,
1028  subsetterPtr,
1029  objects,
1030  surfVectorFields
1031  );
1032 
1033  readFields
1034  (
1035  haveMesh,
1036  mesh,
1037  subsetterPtr,
1038  objects,
1039  surfSphereTensorFields
1040  );
1041 
1042  readFields
1043  (
1044  haveMesh,
1045  mesh,
1046  subsetterPtr,
1047  objects,
1048  surfSymmTensorFields
1049  );
1050 
1051  readFields
1052  (
1053  haveMesh,
1054  mesh,
1055  subsetterPtr,
1056  objects,
1057  surfTensorFields
1058  );
1059 
1060 
1061  // Dimensioned internal fields
1062  readFields
1063  (
1064  haveMesh,
1065  mesh,
1066  subsetterPtr,
1067  objects,
1068  dimScalarFields
1069  );
1070 
1071  readFields
1072  (
1073  haveMesh,
1074  mesh,
1075  subsetterPtr,
1076  objects,
1077  dimVectorFields
1078  );
1079 
1080  readFields
1081  (
1082  haveMesh,
1083  mesh,
1084  subsetterPtr,
1085  objects,
1086  dimSphereTensorFields
1087  );
1088 
1089  readFields
1090  (
1091  haveMesh,
1092  mesh,
1093  subsetterPtr,
1094  objects,
1095  dimSymmTensorFields
1096  );
1097 
1098  readFields
1099  (
1100  haveMesh,
1101  mesh,
1102  subsetterPtr,
1103  objects,
1104  dimTensorFields
1105  );
1106 
1107 
1108  // pointFields currently not supported. Read their names so we
1109  // can delete them.
1110  {
1111  // Get my objects of type
1112  pointFieldNames.append
1113  (
1114  objects.lookupClass(pointScalarField::typeName).sortedNames()
1115  );
1116  pointFieldNames.append
1117  (
1118  objects.lookupClass(pointVectorField::typeName).sortedNames()
1119  );
1120  pointFieldNames.append
1121  (
1122  objects.lookupClass
1123  (
1124  pointSphericalTensorField::typeName
1125  ).sortedNames()
1126  );
1127  pointFieldNames.append
1128  (
1129  objects.lookupClass
1130  (
1131  pointSymmTensorField::typeName
1132  ).sortedNames()
1133  );
1134  pointFieldNames.append
1135  (
1136  objects.lookupClass(pointTensorField::typeName).sortedNames()
1137  );
1138 
1139  // Make sure all processors have the same set
1140  Pstream::scatter(pointFieldNames);
1141  }
1142 
1143  if (Pstream::master() && decompose)
1144  {
1145  runTime.caseName() = proc0CaseName;
1146  }
1147  }
1148 
1149 
1150  // Mesh distribution engine
1151  fvMeshDistribute distributor(mesh, tolDim);
1152 
1153  // Do all the distribution of mesh and fields
1154  autoPtr<mapDistributePolyMesh> rawMap = distributor.distribute(decomp);
1155 
1156  // Print some statistics
1157  Info<< "After distribution:" << endl;
1158  printMeshData(mesh);
1159 
1160  // Get other side of processor boundaries
1161  correctCoupledBoundaryConditions
1162  <
1165  >(mesh);
1166  correctCoupledBoundaryConditions
1167  <
1170  >(mesh);
1171  correctCoupledBoundaryConditions
1172  <
1175  >(mesh);
1176  correctCoupledBoundaryConditions
1177  <
1180  >(mesh);
1181  correctCoupledBoundaryConditions
1182  <
1185  >(mesh);
1186  // No update surface fields
1187 
1188 
1189  // Set the minimum write precision
1191 
1192 
1193  if (!overwrite)
1194  {
1195  ++runTime;
1197  }
1198  else
1199  {
1200  mesh.setInstance(masterInstDir);
1201  }
1202 
1203 
1205  (
1206  IOobject
1207  (
1208  "procAddressing",
1209  mesh.facesInstance(),
1211  mesh,
1214  )
1215  );
1216  map.transfer(rawMap());
1217 
1218 
1219  if (reconstruct)
1220  {
1221  if (Pstream::master())
1222  {
1223  Info<< "Setting caseName to " << baseRunTime.caseName()
1224  << " to write reconstructed mesh and fields." << endl;
1225  runTime.caseName() = baseRunTime.caseName();
1226 
1227  mesh.write();
1229  for (const word& fieldName : pointFieldNames)
1230  {
1231  IOobject io
1232  (
1233  fieldName,
1234  runTime.timeName(),
1235  mesh
1236  );
1237 
1238  const fileName fieldFile(io.objectPath());
1239  if (topoSet::debug) DebugVar(fieldFile);
1240  rm(fieldFile);
1241  }
1242 
1243  // Now we've written all. Reset caseName on master
1244  Info<< "Restoring caseName to " << proc0CaseName << endl;
1245  runTime.caseName() = proc0CaseName;
1246  }
1247  }
1248  else
1249  {
1250  mesh.write();
1252  for (const word& fieldName : pointFieldNames)
1253  {
1254  IOobject io
1255  (
1256  fieldName,
1257  runTime.timeName(),
1258  mesh
1259  );
1260 
1261  const fileName fieldFile(io.objectPath());
1262  if (topoSet::debug) DebugVar(fieldFile);
1263  rm(fieldFile);
1264  }
1265  }
1266  Info<< "Written redistributed mesh to " << mesh.facesInstance() << nl
1267  << endl;
1268 
1269 
1270  if (decompose || reconstruct)
1271  {
1272  // Decompose (1 -> N) or reconstruct (N -> 1)
1273  // so {boundary,cell,face,point}ProcAddressing have meaning
1274  writeProcAddressing(mesh, map, decompose);
1275  }
1276  else
1277  {
1278  // Redistribute (N -> M)
1279  // {boundary,cell,face,point}ProcAddressing would be incorrect
1280  // - can either remove or redistribute previous
1281  removeProcAddressing(mesh);
1282  }
1283 
1284 
1285  // Refinement data
1286  {
1287 
1288  // Read refinement data
1289  if (Pstream::master() && decompose)
1290  {
1291  runTime.caseName() = baseRunTime.caseName();
1292  }
1293  IOobject io
1294  (
1295  "dummy",
1296  mesh.facesInstance(),
1298  mesh,
1301  false
1302  );
1303 
1304  hexRef8Data refData(io);
1305  if (Pstream::master() && decompose)
1306  {
1307  runTime.caseName() = proc0CaseName;
1308  }
1309 
1310  // Make sure all processors have valid data (since only some will
1311  // read)
1312  refData.sync(io);
1313 
1314  // Distribute
1315  refData.distribute(map);
1316 
1317 
1318  // Now we've read refinement data we can remove it
1320 
1321  if (reconstruct)
1322  {
1323  if (Pstream::master())
1324  {
1325  Info<< "Setting caseName to " << baseRunTime.caseName()
1326  << " to write reconstructed refinement data." << endl;
1327  runTime.caseName() = baseRunTime.caseName();
1328 
1329  refData.write();
1330 
1331  // Now we've written all. Reset caseName on master
1332  Info<< "Restoring caseName to " << proc0CaseName << endl;
1333  runTime.caseName() = proc0CaseName;
1334  }
1335  }
1336  else
1337  {
1338  refData.write();
1339  }
1340  }
1341 
1343  //{
1344  // // Read sets
1345  // if (Pstream::master() && decompose)
1346  // {
1347  // runTime.caseName() = baseRunTime.caseName();
1348  // }
1349  // IOobjectList objects(mesh, mesh.facesInstance(), "polyMesh/sets");
1350  //
1351  // PtrList<cellSet> cellSets;
1352  // ReadFields(objects, cellSets);
1353  //
1354  // if (Pstream::master() && decompose)
1355  // {
1356  // runTime.caseName() = proc0CaseName;
1357  // }
1358  //
1359  // forAll(cellSets, i)
1360  // {
1361  // cellSets[i].distribute(map);
1362  // }
1363  //
1364  // if (reconstruct)
1365  // {
1366  // if (Pstream::master())
1367  // {
1368  // Info<< "Setting caseName to " << baseRunTime.caseName()
1369  // << " to write reconstructed refinement data." << endl;
1370  // runTime.caseName() = baseRunTime.caseName();
1371  //
1372  // forAll(cellSets, i)
1373  // {
1374  // cellSets[i].distribute(map);
1375  // }
1376  //
1377  // // Now we've written all. Reset caseName on master
1378  // Info<< "Restoring caseName to " << proc0CaseName << endl;
1379  // runTime.caseName() = proc0CaseName;
1380  // }
1381  // }
1382  // else
1383  // {
1384  // forAll(cellSets, i)
1385  // {
1386  // cellSets[i].distribute(map);
1387  // }
1388  // }
1389  //}
1390 
1391 
1392  return autoPtr<mapDistributePolyMesh>::New(std::move(map));
1393 }
1394 
1395 
1396 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1397 //
1398 // Field Mapping
1399 //
1400 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1401 
1402 autoPtr<mapDistributePolyMesh> createReconstructMap
1403 (
1404  const autoPtr<fvMesh>& baseMeshPtr,
1405  const fvMesh& mesh,
1406  const labelList& cellProcAddressing,
1408  const labelList& pointProcAddressing,
1409  const labelList& boundaryProcAddressing
1410 )
1411 {
1412  // Send addressing to master
1413  labelListList cellAddressing(Pstream::nProcs());
1414  cellAddressing[Pstream::myProcNo()] = cellProcAddressing;
1415  Pstream::gatherList(cellAddressing);
1416 
1417  labelListList faceAddressing(Pstream::nProcs());
1418  faceAddressing[Pstream::myProcNo()] = faceProcAddressing;
1419  Pstream::gatherList(faceAddressing);
1420 
1421  labelListList pointAddressing(Pstream::nProcs());
1422  pointAddressing[Pstream::myProcNo()] = pointProcAddressing;
1423  Pstream::gatherList(pointAddressing);
1424 
1425  labelListList boundaryAddressing(Pstream::nProcs());
1426  {
1427  // Remove -1 entries
1428  DynamicList<label> patchProcAddressing(boundaryProcAddressing.size());
1429  forAll(boundaryProcAddressing, i)
1430  {
1431  if (boundaryProcAddressing[i] != -1)
1432  {
1433  patchProcAddressing.append(boundaryProcAddressing[i]);
1434  }
1435  }
1436  boundaryAddressing[Pstream::myProcNo()] = patchProcAddressing;
1437  Pstream::gatherList(boundaryAddressing);
1438  }
1439 
1440 
1442 
1443  if (baseMeshPtr && baseMeshPtr->nCells())
1444  {
1445  const fvMesh& baseMesh = *baseMeshPtr;
1446 
1447  labelListList cellSubMap(Pstream::nProcs());
1448  cellSubMap[Pstream::masterNo()] = identity(mesh.nCells());
1449 
1450  mapDistribute cellMap
1451  (
1452  baseMesh.nCells(),
1453  std::move(cellSubMap),
1454  std::move(cellAddressing)
1455  );
1456 
1457  labelListList faceSubMap(Pstream::nProcs());
1458  faceSubMap[Pstream::masterNo()] = identity(mesh.nFaces());
1459 
1461  (
1462  baseMesh.nFaces(),
1463  std::move(faceSubMap),
1464  std::move(faceAddressing),
1465  false, //subHasFlip
1466  true //constructHasFlip
1467  );
1468 
1469  labelListList pointSubMap(Pstream::nProcs());
1470  pointSubMap[Pstream::masterNo()] = identity(mesh.nPoints());
1471 
1472  mapDistribute pointMap
1473  (
1474  baseMesh.nPoints(),
1475  std::move(pointSubMap),
1476  std::move(pointAddressing)
1477  );
1478 
1479  labelListList patchSubMap(Pstream::nProcs());
1480  // Send (filtered) patches to master
1481  patchSubMap[Pstream::masterNo()] =
1482  boundaryAddressing[Pstream::myProcNo()];
1483 
1484  mapDistribute patchMap
1485  (
1486  baseMesh.boundaryMesh().size(),
1487  std::move(patchSubMap),
1488  std::move(boundaryAddressing)
1489  );
1490 
1491  const label nOldPoints = mesh.nPoints();
1492  const label nOldFaces = mesh.nFaces();
1493  const label nOldCells = mesh.nCells();
1494 
1495  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
1496  labelList oldPatchStarts(pbm.size());
1497  labelList oldPatchNMeshPoints(pbm.size());
1498  forAll(pbm, patchI)
1499  {
1500  oldPatchStarts[patchI] = pbm[patchI].start();
1501  oldPatchNMeshPoints[patchI] = pbm[patchI].nPoints();
1502  }
1503 
1504  mapPtr.reset
1505  (
1507  (
1508  nOldPoints,
1509  nOldFaces,
1510  nOldCells,
1511  std::move(oldPatchStarts),
1512  std::move(oldPatchNMeshPoints),
1513  std::move(pointMap),
1514  std::move(faceMap),
1515  std::move(cellMap),
1516  std::move(patchMap)
1517  )
1518  );
1519  }
1520  else
1521  {
1522  labelListList cellSubMap(Pstream::nProcs());
1523  cellSubMap[Pstream::masterNo()] = identity(mesh.nCells());
1524  labelListList cellConstructMap(Pstream::nProcs());
1525 
1526  mapDistribute cellMap
1527  (
1528  0,
1529  std::move(cellSubMap),
1530  std::move(cellConstructMap)
1531  );
1532 
1533  labelListList faceSubMap(Pstream::nProcs());
1534  faceSubMap[Pstream::masterNo()] = identity(mesh.nFaces());
1535  labelListList faceConstructMap(Pstream::nProcs());
1536 
1538  (
1539  0,
1540  std::move(faceSubMap),
1541  std::move(faceConstructMap),
1542  false, //subHasFlip
1543  true //constructHasFlip
1544  );
1545 
1546  labelListList pointSubMap(Pstream::nProcs());
1547  pointSubMap[Pstream::masterNo()] = identity(mesh.nPoints());
1548  labelListList pointConstructMap(Pstream::nProcs());
1549 
1550  mapDistribute pointMap
1551  (
1552  0,
1553  std::move(pointSubMap),
1554  std::move(pointConstructMap)
1555  );
1556 
1557  labelListList patchSubMap(Pstream::nProcs());
1558  // Send (filtered) patches to master
1559  patchSubMap[Pstream::masterNo()] =
1560  boundaryAddressing[Pstream::myProcNo()];
1561  labelListList patchConstructMap(Pstream::nProcs());
1562 
1563  mapDistribute patchMap
1564  (
1565  0,
1566  std::move(patchSubMap),
1567  std::move(patchConstructMap)
1568  );
1569 
1570  const label nOldPoints = mesh.nPoints();
1571  const label nOldFaces = mesh.nFaces();
1572  const label nOldCells = mesh.nCells();
1573 
1574  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
1575  labelList oldPatchStarts(pbm.size());
1576  labelList oldPatchNMeshPoints(pbm.size());
1577  forAll(pbm, patchI)
1578  {
1579  oldPatchStarts[patchI] = pbm[patchI].start();
1580  oldPatchNMeshPoints[patchI] = pbm[patchI].nPoints();
1581  }
1582 
1583  mapPtr.reset
1584  (
1586  (
1587  nOldPoints,
1588  nOldFaces,
1589  nOldCells,
1590  std::move(oldPatchStarts),
1591  std::move(oldPatchNMeshPoints),
1592  std::move(pointMap),
1593  std::move(faceMap),
1594  std::move(cellMap),
1595  std::move(patchMap)
1596  )
1597  );
1598  }
1599 
1600  return mapPtr;
1601 }
1602 
1603 
1604 void readProcAddressing
1605 (
1606  const fvMesh& mesh,
1607  const autoPtr<fvMesh>& baseMeshPtr,
1609 )
1610 {
1611  //IOobject io
1612  //(
1613  // "procAddressing",
1614  // mesh.facesInstance(),
1615  // polyMesh::meshSubDir,
1616  // mesh,
1617  // IOobject::MUST_READ
1618  //);
1619  //if (io.typeHeaderOk<labelIOList>(true))
1620  //{
1621  // Pout<< "Reading addressing from " << io.name() << " at "
1622  // << mesh.facesInstance() << nl << endl;
1623  // distMap.reset(new IOmapDistributePolyMesh(io));
1624  //}
1625  //else
1626  {
1627  Info<< "Reading addressing from procXXXAddressing at "
1628  << mesh.facesInstance() << nl << endl;
1629  labelIOList cellProcAddressing
1630  (
1631  IOobject
1632  (
1633  "cellProcAddressing",
1634  mesh.facesInstance(),
1636  mesh,
1638  ),
1639  labelList()
1640  );
1642  (
1643  IOobject
1644  (
1645  "faceProcAddressing",
1646  mesh.facesInstance(),
1648  mesh,
1650  ),
1651  labelList()
1652  );
1653  labelIOList pointProcAddressing
1654  (
1655  IOobject
1656  (
1657  "pointProcAddressing",
1658  mesh.facesInstance(),
1660  mesh,
1662  ),
1663  labelList()
1664  );
1665  labelIOList boundaryProcAddressing
1666  (
1667  IOobject
1668  (
1669  "boundaryProcAddressing",
1670  mesh.facesInstance(),
1672  mesh,
1674  ),
1675  labelList()
1676  );
1677 
1678 
1679  if
1680  (
1681  mesh.nCells() != cellProcAddressing.size()
1682  || mesh.nPoints() != pointProcAddressing.size()
1683  || mesh.nFaces() != faceProcAddressing.size()
1684  || mesh.boundaryMesh().size() != boundaryProcAddressing.size()
1685  )
1686  {
1688  << "Read addressing inconsistent with mesh sizes" << nl
1689  << "cells:" << mesh.nCells()
1690  << " addressing:" << cellProcAddressing.objectPath()
1691  << " size:" << cellProcAddressing.size() << nl
1692  << "faces:" << mesh.nFaces()
1693  << " addressing:" << faceProcAddressing.objectPath()
1694  << " size:" << faceProcAddressing.size() << nl
1695  << "points:" << mesh.nPoints()
1696  << " addressing:" << pointProcAddressing.objectPath()
1697  << " size:" << pointProcAddressing.size()
1698  << "patches:" << mesh.boundaryMesh().size()
1699  << " addressing:" << boundaryProcAddressing.objectPath()
1700  << " size:" << boundaryProcAddressing.size()
1701  << exit(FatalError);
1702  }
1703 
1704  distMap.clear();
1705  distMap = createReconstructMap
1706  (
1707  baseMeshPtr,
1708  mesh,
1709  cellProcAddressing,
1711  pointProcAddressing,
1712  boundaryProcAddressing
1713  );
1714  }
1715 }
1716 
1717 
1718 void reconstructMeshFields
1719 (
1720  const parFvFieldReconstructor& fvReconstructor,
1721  const IOobjectList& objects,
1722  const wordRes& selectedFields
1723 )
1724 {
1725  // Dimensioned fields
1726 
1727  fvReconstructor.reconstructFvVolumeInternalFields<scalar>
1728  (
1729  objects,
1730  selectedFields
1731  );
1732  fvReconstructor.reconstructFvVolumeInternalFields<vector>
1733  (
1734  objects,
1735  selectedFields
1736  );
1738  (
1739  objects,
1740  selectedFields
1741  );
1743  (
1744  objects,
1745  selectedFields
1746  );
1747  fvReconstructor.reconstructFvVolumeInternalFields<tensor>
1748  (
1749  objects,
1750  selectedFields
1751  );
1752 
1753 
1754  // volFields
1755 
1756  fvReconstructor.reconstructFvVolumeFields<scalar>
1757  (
1758  objects,
1759  selectedFields
1760  );
1761  fvReconstructor.reconstructFvVolumeFields<vector>
1762  (
1763  objects,
1764  selectedFields
1765  );
1767  (
1768  objects,
1769  selectedFields
1770  );
1771  fvReconstructor.reconstructFvVolumeFields<symmTensor>
1772  (
1773  objects,
1774  selectedFields
1775  );
1776  fvReconstructor.reconstructFvVolumeFields<tensor>
1777  (
1778  objects,
1779  selectedFields
1780  );
1781 
1782 
1783  // surfaceFields
1784 
1785  fvReconstructor.reconstructFvSurfaceFields<scalar>
1786  (
1787  objects,
1788  selectedFields
1789  );
1790  fvReconstructor.reconstructFvSurfaceFields<vector>
1791  (
1792  objects,
1793  selectedFields
1794  );
1796  (
1797  objects,
1798  selectedFields
1799  );
1800  fvReconstructor.reconstructFvSurfaceFields<symmTensor>
1801  (
1802  objects,
1803  selectedFields
1804  );
1805  fvReconstructor.reconstructFvSurfaceFields<tensor>
1806  (
1807  objects,
1808  selectedFields
1809  );
1810 }
1811 
1812 
1813 void reconstructLagrangian
1814 (
1815  autoPtr<parLagrangianRedistributor>& lagrangianReconstructorPtr,
1816  const fvMesh& baseMesh,
1817  const fvMesh& mesh,
1818  const mapDistributePolyMesh& distMap,
1819  const wordRes& selectedLagrangianFields
1820 )
1821 {
1822  // Clouds (note: might not be present on all processors)
1823 
1825  List<wordList> fieldNames;
1826  // Find all cloudNames on all processors
1828 
1829  if (cloudNames.size())
1830  {
1831  if (!lagrangianReconstructorPtr)
1832  {
1833  lagrangianReconstructorPtr.reset
1834  (
1836  (
1837  mesh,
1838  baseMesh,
1839  mesh.nCells(), // range of cell indices in clouds
1840  distMap
1841  )
1842  );
1843  }
1845  *lagrangianReconstructorPtr;
1846 
1847  for (const word& cloudName : cloudNames)
1848  {
1849  Info<< "Reconstructing lagrangian fields for cloud "
1850  << cloudName << nl << endl;
1851 
1852  autoPtr<mapDistributeBase> lagrangianMapPtr =
1853  lagrangianReconstructor.redistributeLagrangianPositions
1854  (
1855  cloudName
1856  );
1857 
1858  const mapDistributeBase& lagrangianMap = *lagrangianMapPtr;
1859 
1860  IOobjectList cloudObjs
1861  (
1862  mesh,
1863  mesh.time().timeName(),
1865  );
1866 
1867  lagrangianReconstructor.redistributeFields<label>
1868  (
1869  lagrangianMap,
1870  cloudName,
1871  cloudObjs,
1872  selectedLagrangianFields
1873  );
1874  lagrangianReconstructor.redistributeFieldFields<label>
1875  (
1876  lagrangianMap,
1877  cloudName,
1878  cloudObjs,
1879  selectedLagrangianFields
1880  );
1881  lagrangianReconstructor.redistributeFields<scalar>
1882  (
1883  lagrangianMap,
1884  cloudName,
1885  cloudObjs,
1886  selectedLagrangianFields
1887  );
1888  lagrangianReconstructor.redistributeFieldFields<scalar>
1889  (
1890  lagrangianMap,
1891  cloudName,
1892  cloudObjs,
1893  selectedLagrangianFields
1894  );
1895  lagrangianReconstructor.redistributeFields<vector>
1896  (
1897  lagrangianMap,
1898  cloudName,
1899  cloudObjs,
1900  selectedLagrangianFields
1901  );
1902  lagrangianReconstructor.redistributeFieldFields<vector>
1903  (
1904  lagrangianMap,
1905  cloudName,
1906  cloudObjs,
1907  selectedLagrangianFields
1908  );
1909  lagrangianReconstructor.redistributeFields
1910  <sphericalTensor>
1911  (
1912  lagrangianMap,
1913  cloudName,
1914  cloudObjs,
1915  selectedLagrangianFields
1916  );
1917  lagrangianReconstructor.redistributeFieldFields
1918  <sphericalTensor>
1919  (
1920  lagrangianMap,
1921  cloudName,
1922  cloudObjs,
1923  selectedLagrangianFields
1924  );
1925  lagrangianReconstructor.redistributeFields<symmTensor>
1926  (
1927  lagrangianMap,
1928  cloudName,
1929  cloudObjs,
1930  selectedLagrangianFields
1931  );
1932  lagrangianReconstructor.redistributeFieldFields
1933  <symmTensor>
1934  (
1935  lagrangianMap,
1936  cloudName,
1937  cloudObjs,
1938  selectedLagrangianFields
1939  );
1940  lagrangianReconstructor.redistributeFields<tensor>
1941  (
1942  lagrangianMap,
1943  cloudName,
1944  cloudObjs,
1945  selectedLagrangianFields
1946  );
1947  lagrangianReconstructor.redistributeFieldFields<tensor>
1948  (
1949  lagrangianMap,
1950  cloudName,
1951  cloudObjs,
1952  selectedLagrangianFields
1953  );
1954  }
1955  }
1956 }
1957 
1958 
1959 // Read clouds (note: might not be present on all processors)
1960 void readLagrangian
1961 (
1962  const fvMesh& mesh,
1963  const wordList& cloudNames,
1964  const wordRes& selectedLagrangianFields,
1966 )
1967 {
1968  (void)mesh.tetBasePtIs();
1969 
1970  forAll(cloudNames, i)
1971  {
1972  //Pout<< "Loading cloud " << cloudNames[i] << endl;
1973  clouds.set
1974  (
1975  i,
1977  );
1978 
1979 
1980  //for (passivePositionParticle& p : clouds[i]))
1981  //{
1982  // Pout<< "Particle position:" << p.position()
1983  // << " cell:" << p.cell()
1984  // << " with cc:" << mesh.cellCentres()[p.cell()]
1985  // << endl;
1986  //}
1987 
1988 
1989  IOobjectList cloudObjs(clouds[i], clouds[i].time().timeName());
1990 
1991  //Pout<< "Found clould objects:" << cloudObjs.names() << endl;
1992 
1994  <IOField<label>>
1995  (
1996  clouds[i],
1997  cloudObjs,
1998  selectedLagrangianFields
1999  );
2002  (
2003  clouds[i],
2004  cloudObjs,
2005  selectedLagrangianFields
2006  );
2008  <CompactIOField<Field<label>, label>>
2009  (
2010  clouds[i],
2011  cloudObjs,
2012  selectedLagrangianFields
2013  );
2014 
2015 
2017  <IOField<scalar>>
2018  (
2019  clouds[i],
2020  cloudObjs,
2021  selectedLagrangianFields
2022  );
2025  (
2026  clouds[i],
2027  cloudObjs,
2028  selectedLagrangianFields
2029  );
2031  <CompactIOField<Field<scalar>, scalar>>
2032  (
2033  clouds[i],
2034  cloudObjs,
2035  selectedLagrangianFields
2036  );
2037 
2038 
2040  <IOField<vector>>
2041  (
2042  clouds[i],
2043  cloudObjs,
2044  selectedLagrangianFields
2045  );
2048  (
2049  clouds[i],
2050  cloudObjs,
2051  selectedLagrangianFields
2052  );
2055  (
2056  clouds[i],
2057  cloudObjs,
2058  selectedLagrangianFields
2059  );
2060 
2061 
2064  (
2065  clouds[i],
2066  cloudObjs,
2067  selectedLagrangianFields
2068  );
2071  (
2072  clouds[i],
2073  cloudObjs,
2074  selectedLagrangianFields
2075  );
2078  (
2079  clouds[i],
2080  cloudObjs,
2081  selectedLagrangianFields
2082  );
2083 
2084 
2087  (
2088  clouds[i],
2089  cloudObjs,
2090  selectedLagrangianFields
2091  );
2094  (
2095  clouds[i],
2096  cloudObjs,
2097  selectedLagrangianFields
2098  );
2101  (
2102  clouds[i],
2103  cloudObjs,
2104  selectedLagrangianFields
2105  );
2106 
2107 
2109  <IOField<tensor>>
2110  (
2111  clouds[i],
2112  cloudObjs,
2113  selectedLagrangianFields
2114  );
2117  (
2118  clouds[i],
2119  cloudObjs,
2120  selectedLagrangianFields
2121  );
2124  (
2125  clouds[i],
2126  cloudObjs,
2127  selectedLagrangianFields
2128  );
2129  }
2130 }
2131 
2132 
2133 void redistributeLagrangian
2134 (
2135  autoPtr<parLagrangianRedistributor>& lagrangianReconstructorPtr,
2136  const fvMesh& mesh,
2137  const label nOldCells,
2138  const mapDistributePolyMesh& distMap,
2140 )
2141 {
2142  if (clouds.size())
2143  {
2144  if (!lagrangianReconstructorPtr)
2145  {
2146  lagrangianReconstructorPtr.reset
2147  (
2149  (
2150  mesh,
2151  mesh,
2152  nOldCells, // range of cell indices in clouds
2153  distMap
2154  )
2155  );
2156  }
2157  const parLagrangianRedistributor& distributor =
2158  lagrangianReconstructorPtr();
2159 
2160  forAll(clouds, i)
2161  {
2162  autoPtr<mapDistributeBase> lagrangianMapPtr =
2163  distributor.redistributeLagrangianPositions(clouds[i]);
2164  const mapDistributeBase& lagrangianMap = *lagrangianMapPtr;
2165 
2166  distributor.redistributeStoredFields
2167  <IOField<label>>
2168  (
2169  lagrangianMap,
2170  clouds[i]
2171  );
2172  distributor.redistributeStoredFields
2174  (
2175  lagrangianMap,
2176  clouds[i]
2177  );
2178  distributor.redistributeStoredFields
2179  <CompactIOField<Field<label>, label>>
2180  (
2181  lagrangianMap,
2182  clouds[i]
2183  );
2184 
2185 
2186  distributor.redistributeStoredFields
2187  <IOField<scalar>>
2188  (
2189  lagrangianMap,
2190  clouds[i]
2191  );
2192  distributor.redistributeStoredFields
2194  (
2195  lagrangianMap,
2196  clouds[i]
2197  );
2198  distributor.redistributeStoredFields
2199  <CompactIOField<Field<scalar>, scalar>>
2200  (
2201  lagrangianMap,
2202  clouds[i]
2203  );
2204 
2205 
2206  distributor.redistributeStoredFields
2207  <IOField<vector>>
2208  (
2209  lagrangianMap,
2210  clouds[i]
2211  );
2212  distributor.redistributeStoredFields
2214  (
2215  lagrangianMap,
2216  clouds[i]
2217  );
2218  distributor.redistributeStoredFields
2220  (
2221  lagrangianMap,
2222  clouds[i]
2223  );
2224 
2225 
2226  distributor.redistributeStoredFields
2228  (
2229  lagrangianMap,
2230  clouds[i]
2231  );
2232  distributor.redistributeStoredFields
2234  (
2235  lagrangianMap,
2236  clouds[i]
2237  );
2238  distributor.redistributeStoredFields
2240  (
2241  lagrangianMap,
2242  clouds[i]
2243  );
2244 
2245 
2246  distributor.redistributeStoredFields
2248  (
2249  lagrangianMap,
2250  clouds[i]
2251  );
2252  distributor.redistributeStoredFields
2254  (
2255  lagrangianMap,
2256  clouds[i]
2257  );
2258  distributor.redistributeStoredFields
2260  (
2261  lagrangianMap,
2262  clouds[i]
2263  );
2264 
2265 
2266  distributor.redistributeStoredFields
2267  <IOField<tensor>>
2268  (
2269  lagrangianMap,
2270  clouds[i]
2271  );
2272  distributor.redistributeStoredFields
2274  (
2275  lagrangianMap,
2276  clouds[i]
2277  );
2278  distributor.redistributeStoredFields
2280  (
2281  lagrangianMap,
2282  clouds[i]
2283  );
2284  }
2285  }
2286 }
2287 
2288 
2289 int main(int argc, char *argv[])
2290 {
2292  (
2293  "Redistribute decomposed mesh and fields according"
2294  " to the decomposeParDict settings.\n"
2295  "Optionally run in decompose/reconstruct mode"
2296  );
2297 
2298  argList::noFunctionObjects(); // Never use function objects
2299 
2300  // enable -constant ... if someone really wants it
2301  // enable -zeroTime to prevent accidentally trashing the initial fields
2302  timeSelector::addOptions(true, true);
2303  #include "addRegionOption.H"
2305  (
2306  "allRegions",
2307  "operate on all regions in regionProperties"
2308  );
2309  #include "addOverwriteOption.H"
2310  argList::addBoolOption("decompose", "Decompose case");
2311  argList::addBoolOption("reconstruct", "Reconstruct case");
2313  (
2314  "dry-run",
2315  "Test without writing the decomposition. "
2316  "Changes -cellDist to only write volScalarField."
2317  );
2319  (
2320  "mergeTol",
2321  "scalar",
2322  "The merge distance relative to the bounding box size (default 1e-6)"
2323  );
2325  (
2326  "cellDist",
2327  "Write cell distribution as a labelList - for use with 'manual' "
2328  "decomposition method or as a volScalarField for post-processing."
2329  );
2331  (
2332  "newTimes",
2333  "Only reconstruct new times (i.e. that do not exist already)"
2334  );
2335 
2336 
2337  // Handle arguments
2338  // ~~~~~~~~~~~~~~~~
2339  // (replacement for setRootCase that does not abort)
2340 
2341  argList args(argc, argv);
2342 
2343 
2344  // As much as possible avoid synchronised operation. To be looked at more
2345  // closely for the three scenarios:
2346  // - decompose - reads on master (and from parent directory) and sends
2347  // dictionary to slaves
2348  // - distribute - reads on potentially a different number of processors
2349  // than it writes to
2350  // - reconstruct - reads parallel, write on master only and to parent
2351  // directory
2352  const_cast<fileOperation&>(fileHandler()).distributed(true);
2353 
2354 
2355  #include "foamDlOpenLibs.H"
2356 
2357  const bool reconstruct = args.found("reconstruct");
2358  const bool writeCellDist = args.found("cellDist");
2359  const bool dryrun = args.found("dry-run");
2360  const bool newTimes = args.found("newTimes");
2361 
2362  bool decompose = args.found("decompose");
2363  bool overwrite = args.found("overwrite");
2364 
2365  // Disable NaN setting and floating point error trapping. This is to avoid
2366  // any issues inside the field redistribution inside fvMeshDistribute
2367  // which temporarily moves processor faces into existing patches. These
2368  // will now not have correct values. After all bits have been assembled
2369  // the processor fields will be restored but until then there might
2370  // be uninitialised values which might upset some patch field constructors.
2371  // Workaround by disabling floating point error trapping. TBD: have
2372  // actual field redistribution instead of subsetting inside
2373  // fvMeshDistribute.
2374  Foam::sigFpe::unset(true);
2375 
2376  const wordRes selectedFields;
2377  const wordRes selectedLagrangianFields;
2378 
2379 
2380  if (decompose)
2381  {
2382  Info<< "Decomposing case (like decomposePar)" << nl << endl;
2383  if (reconstruct)
2384  {
2386  << "Cannot specify both -decompose and -reconstruct"
2387  << exit(FatalError);
2388  }
2389  }
2390  else if (reconstruct)
2391  {
2392  Info<< "Reconstructing case (like reconstructParMesh)" << nl << endl;
2393  }
2394 
2395 
2396  if (decompose || reconstruct)
2397  {
2398  if (!overwrite)
2399  {
2401  << "Working in decompose or reconstruction mode automatically"
2402  << " implies -overwrite" << nl << endl;
2403  overwrite = true;
2404  }
2405  }
2406 
2407 
2408  if (!Pstream::parRun())
2409  {
2411  << ": This utility can only be run parallel"
2412  << exit(FatalError);
2413  }
2414 
2415 
2416  if (!isDir(args.rootPath()))
2417  {
2419  << ": cannot open root directory " << args.rootPath()
2420  << exit(FatalError);
2421  }
2422 
2423  // Detect if running data-distributed (multiple roots)
2424  bool nfs = true;
2425  {
2426  List<fileName> roots(1, args.rootPath());
2428  nfs = (roots.size() == 1);
2429  }
2430 
2431  if (!nfs)
2432  {
2433  Info<< "Detected multiple roots i.e. non-nfs running"
2434  << nl << endl;
2435  }
2436 
2437  if (isDir(args.path()))
2438  {
2439  if (decompose)
2440  {
2441  Info<< "Removing existing processor directories" << endl;
2442  rmDir(args.path());
2443  }
2444  }
2445  else
2446  {
2447  // Directory does not exist. If this happens on master -> decompose mode
2448  decompose = true;
2449  Info<< "No processor directories; switching on decompose mode"
2450  << nl << endl;
2451  }
2452  // If master changed to decompose mode make sure all nodes know about it
2453  Pstream::scatter(decompose);
2454 
2455 
2456  // If running distributed we have problem of new processors not finding
2457  // a system/controlDict. However if we switch on the master-only reading
2458  // the problem becomes that the time directories are differing sizes and
2459  // e.g. latestTime will pick up a different time (which causes createTime.H
2460  // to abort). So for now make sure to have master times on all
2461  // processors
2462  {
2463  Info<< "Creating time directories on all processors" << nl << endl;
2464  instantList timeDirs;
2465  if (Pstream::master())
2466  {
2467  const bool oldParRun = Pstream::parRun(false);
2468  timeDirs = Time::findTimes(args.path(), "constant");
2469  Pstream::parRun(oldParRun); // Restore parallel state
2470  }
2471  Pstream::scatter(timeDirs);
2472  for (const instant& t : timeDirs)
2473  {
2474  mkDir(args.path()/t.name());
2475  }
2476  }
2477 
2478 
2479  // Construct time
2480  // ~~~~~~~~~~~~~~
2481 
2482  #include "createTime.H"
2483  runTime.functionObjects().off(); // Extra safety?
2484 
2485 
2486  // Save local processor0 casename
2487  const fileName proc0CaseName = runTime.caseName();
2488 
2489 
2490  // Construct undecomposed Time
2491  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
2492  // This will read the same controlDict but might have a different
2493  // set of times so enforce same times
2494 
2495  if (!nfs)
2496  {
2497  Info<< "Creating time directories for undecomposed Time"
2498  << " on all processors" << nl << endl;
2499  instantList timeDirs;
2500 
2501  const fileName basePath(args.globalPath());
2502 
2503  if (Pstream::master())
2504  {
2505  const bool oldParRun = Pstream::parRun(false);
2506  timeDirs = Time::findTimes(basePath, "constant");
2507  Pstream::parRun(oldParRun); // Restore parallel state
2508  }
2509  Pstream::scatter(timeDirs);
2510  for (const instant& t : timeDirs)
2511  {
2512  mkDir(basePath/t.name());
2513  }
2514  }
2515 
2516 
2517  Info<< "Create undecomposed database"<< nl << endl;
2518  Time baseRunTime
2519  (
2520  runTime.controlDict(),
2521  runTime.rootPath(),
2523  runTime.system(),
2524  runTime.constant(),
2525  false // enableFunctionObjects
2526  );
2527 
2528 
2529  wordHashSet masterTimeDirSet;
2530  if (newTimes)
2531  {
2532  instantList baseTimeDirs(baseRunTime.times());
2533  for (const instant& t : baseTimeDirs)
2534  {
2535  masterTimeDirSet.insert(t.name());
2536  }
2537  }
2538 
2539 
2540  // Allow override of decomposeParDict location
2541  const fileName decompDictFile =
2542  args.getOrDefault<fileName>("decomposeParDict", "");
2543 
2544  // Get all region names
2545  wordList regionNames;
2546  if (args.found("allRegions"))
2547  {
2548  regionNames = regionProperties(runTime).names();
2549 
2550  Info<< "Decomposing all regions in regionProperties" << nl
2551  << " " << flatOutput(regionNames) << nl << endl;
2552  }
2553  else
2554  {
2555  regionNames.resize(1);
2556  regionNames.first() =
2558  }
2559 
2560 
2561  // Demand driven lagrangian mapper
2562  autoPtr<parLagrangianRedistributor> lagrangianReconstructorPtr;
2563 
2564 
2565  if (reconstruct)
2566  {
2567  // use the times list from the master processor
2568  // and select a subset based on the command-line options
2570  Pstream::scatter(timeDirs);
2571 
2572  if (timeDirs.empty())
2573  {
2575  << "No times selected"
2576  << exit(FatalError);
2577  }
2578 
2579 
2580  // Pass1 : reconstruct mesh and addressing
2581  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2582 
2583 
2584  Info<< nl
2585  << "Pass1 : reconstructing mesh and addressing" << nl << endl;
2586 
2587 
2588  forAll(regionNames, regioni)
2589  {
2590  const word& regionName = regionNames[regioni];
2591  const fileName meshSubDir
2592  (
2595  : regionNames[regioni]/polyMesh::meshSubDir
2596  );
2597 
2598  Info<< "\n\nReconstructing mesh " << regionName << nl << endl;
2599 
2600  // Loop over all times
2601  forAll(timeDirs, timeI)
2602  {
2603  // Set time for global database
2604  runTime.setTime(timeDirs[timeI], timeI);
2605  baseRunTime.setTime(timeDirs[timeI], timeI);
2606 
2607  Info<< "Time = " << runTime.timeName() << endl << endl;
2608 
2609 
2610  // See where the mesh is
2611  //const bool oldParRun = Pstream::parRun(false);
2612  fileName facesInstance = runTime.findInstance
2613  (
2614  meshSubDir,
2615  "faces",
2617  );
2618  //Pstream::parRun(oldParRun);
2619  //Pout<< "facesInstance:" << facesInstance << endl;
2620 
2621  Pstream::scatter(facesInstance);
2622 
2623  // Check who has a mesh
2624  const fileName meshPath =
2625  runTime.path()/facesInstance/meshSubDir/"faces";
2626 
2627  Info<< "Checking for mesh in " << meshPath << nl << endl;
2628 
2629  boolList haveMesh(Pstream::nProcs(), false);
2630  haveMesh[Pstream::myProcNo()] = isFile(meshPath);
2631  Pstream::gatherList(haveMesh);
2632  Pstream::scatterList(haveMesh);
2633  Info<< "Per processor mesh availability:" << nl
2634  << " " << flatOutput(haveMesh) << nl << endl;
2635 
2636 
2637  // Addressing back to reconstructed mesh as xxxProcAddressing.
2638  // - all processors have consistent faceProcAddressing
2639  // - processors with no mesh don't need faceProcAddressing
2640 
2641 
2642  // Note: filePath searches up on processors that don't have
2643  // processor if instance = constant so explicitly check
2644  // found filename.
2645  bool haveAddressing = false;
2646  if (haveMesh[Pstream::myProcNo()])
2647  {
2648  // Read faces (just to know their size)
2649  faceCompactIOList faces
2650  (
2651  IOobject
2652  (
2653  "faces",
2654  facesInstance,
2655  meshSubDir,
2656  runTime,
2658  )
2659  );
2660 
2661  // Check faceProcAddressing
2663  (
2664  IOobject
2665  (
2666  "faceProcAddressing",
2667  facesInstance,
2668  meshSubDir,
2669  runTime,
2671  ),
2672  labelList()
2673  );
2674  if
2675  (
2676  faceProcAddressing.headerOk()
2677  && faceProcAddressing.size() == faces.size()
2678  )
2679  {
2680  haveAddressing = true;
2681  }
2682  }
2683  else
2684  {
2685  // Have no mesh. Don't need addressing
2686  haveAddressing = true;
2687  }
2688 
2689  if (!returnReduce(haveAddressing, andOp<bool>()))
2690  {
2691  Info<< "loading mesh from " << facesInstance << endl;
2693  (
2694  IOobject
2695  (
2696  regionName,
2697  facesInstance,
2698  runTime,
2700  )
2701  );
2702  fvMesh& mesh = meshPtr();
2703 
2704  // Use basic geometry calculation to avoid synchronisation
2705  // problems. See comment in routine
2706  setBasicGeometry(mesh);
2707 
2708  // Global matching tolerance
2709  const scalar tolDim = getMergeDistance
2710  (
2711  args,
2712  runTime,
2713  mesh.bounds()
2714  );
2715 
2716 
2717  // Determine decomposition
2718  // ~~~~~~~~~~~~~~~~~~~~~~~
2719 
2720  Info<< "Reconstructing mesh for time " << facesInstance
2721  << endl;
2722 
2723  label nDestProcs = 1;
2724  labelList finalDecomp = labelList(mesh.nCells(), Zero);
2725 
2726  redistributeAndWrite
2727  (
2728  baseRunTime,
2729  tolDim,
2730  haveMesh,
2731  meshSubDir,
2732  false, // do not read fields
2733  false, // do not read undecomposed case on proc0
2734  true, // write redistributed files to proc0
2735  overwrite,
2736  proc0CaseName,
2737  nDestProcs,
2738  finalDecomp,
2739  facesInstance,
2740  mesh
2741  );
2742  }
2743  }
2744 
2745 
2746  // Pass2 : read mesh and addressing and reconstruct fields
2747  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2748 
2749  Info<< nl
2750  << "Pass2 : reconstructing fields" << nl << endl;
2751 
2752  runTime.setTime(timeDirs[0], 0);
2753  baseRunTime.setTime(timeDirs[0], 0);
2754  Info<< "Time = " << runTime.timeName() << endl << endl;
2755 
2756 
2757  // Read undecomposed mesh on master and 'empty' mesh
2758  // (zero faces, point, cells but valid patches and zones) on slaves.
2759  // This is a bit of tricky code and hidden inside fvMeshTools for
2760  // now.
2761  Info<< "Reading undecomposed mesh (on master)" << endl;
2763  (
2764  IOobject
2765  (
2766  regionName,
2767  baseRunTime.timeName(),
2768  baseRunTime,
2770  ),
2771  true // read on master only
2772  );
2773  setBasicGeometry(baseMeshPtr());
2774 
2775 
2776  Info<< "Reading local, decomposed mesh" << endl;
2778  (
2779  IOobject
2780  (
2781  regionName,
2782  baseMeshPtr().facesInstance(),
2783  runTime,
2785  )
2786  );
2787  fvMesh& mesh = meshPtr();
2788 
2789 
2790  // Read addressing back to base mesh
2792  readProcAddressing(mesh, baseMeshPtr, distMap);
2793 
2794  // Construct field mapper
2795  autoPtr<parFvFieldReconstructor> fvReconstructorPtr
2796  (
2798  (
2799  baseMeshPtr(),
2800  mesh,
2801  distMap(),
2802  Pstream::master() // do I need to write?
2803  )
2804  );
2805 
2806 
2807 
2808  // Since we start from Times[0] and not runTime.timeName() we
2809  // might overlook point motion in the first timestep
2810  // (since mesh.readUpdate() below will not be triggered). Instead
2811  // detect points by hand
2813  {
2814  Info<< " Detected initial mesh motion;"
2815  << " reconstructing points" << nl
2816  << endl;
2817  fvReconstructorPtr().reconstructPoints();
2818  }
2819 
2820 
2821  // Loop over all times
2822  forAll(timeDirs, timeI)
2823  {
2824  if (newTimes && masterTimeDirSet.found(timeDirs[timeI].name()))
2825  {
2826  Info<< "Skipping time " << timeDirs[timeI].name()
2827  << endl << endl;
2828  continue;
2829  }
2830 
2831  // Set time for global database
2832  runTime.setTime(timeDirs[timeI], timeI);
2833  baseRunTime.setTime(timeDirs[timeI], timeI);
2834 
2835  Info<< "Time = " << runTime.timeName() << endl << endl;
2836 
2837 
2838  // Check if any new meshes need to be read.
2840 
2841  if (procStat == fvMesh::POINTS_MOVED)
2842  {
2843  Info<< " Dected mesh motion; reconstructing points" << nl
2844  << endl;
2845  fvReconstructorPtr().reconstructPoints();
2846  }
2847  else if
2848  (
2849  procStat == fvMesh::TOPO_CHANGE
2850  || procStat == fvMesh::TOPO_PATCH_CHANGE
2851  )
2852  {
2853  Info<< " Detected topology change;"
2854  << " reconstructing addressing" << nl << endl;
2855 
2856  if (baseMeshPtr)
2857  {
2858  // Cannot do a baseMesh::readUpdate() since not all
2859  // processors will have mesh files. So instead just
2860  // recreate baseMesh
2861  baseMeshPtr.clear();
2862  baseMeshPtr = fvMeshTools::newMesh
2863  (
2864  IOobject
2865  (
2866  regionName,
2867  baseRunTime.timeName(),
2868  baseRunTime,
2870  ),
2871  true // read on master only
2872  );
2873  }
2874 
2875  // Re-read procXXXaddressing
2876  readProcAddressing(mesh, baseMeshPtr, distMap);
2877 
2878  // Reset field mapper
2879  fvReconstructorPtr.reset
2880  (
2882  (
2883  baseMeshPtr(),
2884  mesh,
2885  distMap(),
2886  Pstream::master()
2887  )
2888  );
2889  lagrangianReconstructorPtr.clear();
2890  }
2891 
2892 
2893  // Get list of objects
2894  IOobjectList objects(mesh, runTime.timeName());
2895 
2896 
2897  // Mesh fields (vol, surface, volInternal)
2898  reconstructMeshFields
2899  (
2900  fvReconstructorPtr(),
2901  objects,
2902  selectedFields
2903  );
2904 
2905  // Clouds (note: might not be present on all processors)
2906  reconstructLagrangian
2907  (
2908  lagrangianReconstructorPtr,
2909  baseMeshPtr(),
2910  mesh,
2911  distMap(),
2912  selectedLagrangianFields
2913  );
2914 
2915  // If there are any "uniform" directories copy them from
2916  // the master processor
2917  if (Pstream::master())
2918  {
2919  fileName uniformDir0 = runTime.timePath()/"uniform";
2920  if (isDir(uniformDir0))
2921  {
2922  Info<< "Detected additional non-decomposed files in "
2923  << uniformDir0 << endl;
2924  cp(uniformDir0, baseRunTime.timePath());
2925  }
2926  }
2927  }
2928  }
2929  }
2930  else
2931  {
2932  // Time coming from processor0 (or undecomposed if no processor0)
2933  scalar masterTime;
2934  if (decompose)
2935  {
2936  // Use base time. This is to handle e.g. startTime = latestTime
2937  // which will not do anything if there are no processor directories
2938  masterTime = timeSelector::selectIfPresent
2939  (
2940  baseRunTime,
2941  args
2942  )[0].value();
2943  }
2944  else
2945  {
2946  masterTime = timeSelector::selectIfPresent
2947  (
2948  runTime,
2949  args
2950  )[0].value();
2951  }
2952  Pstream::scatter(masterTime);
2953  Info<< "Setting time to that of master or undecomposed case : "
2954  << masterTime << endl;
2955  runTime.setTime(masterTime, 0);
2956 
2957 
2958 
2959  forAll(regionNames, regioni)
2960  {
2961  const word& regionName = regionNames[regioni];
2962  const fileName meshSubDir
2963  (
2966  : regionNames[regioni]/polyMesh::meshSubDir
2967  );
2968 
2969  if (decompose)
2970  {
2971  Info<< "\n\nDecomposing mesh " << regionName << nl << endl;
2972  }
2973  else
2974  {
2975  Info<< "\n\nRedistributing mesh " << regionName << nl << endl;
2976  }
2977 
2978 
2979  // Get time instance directory
2980  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
2981  // At this point we should be able to read at least a mesh on
2982  // processor0. Note the changing of the processor0 casename to
2983  // enforce it to read/write from the undecomposed case
2984 
2985  fileName masterInstDir;
2986  if (Pstream::master())
2987  {
2988  if (decompose)
2989  {
2990  Info<< "Setting caseName to " << baseRunTime.caseName()
2991  << " to find undecomposed mesh" << endl;
2992  runTime.caseName() = baseRunTime.caseName();
2993  }
2994 
2995  const bool oldParRun = Pstream::parRun(false);
2996  masterInstDir = runTime.findInstance
2997  (
2998  meshSubDir,
2999  "faces",
3001  );
3002  Pstream::parRun(oldParRun);
3003 
3004  if (decompose)
3005  {
3006  Info<< "Restoring caseName to " << proc0CaseName << endl;
3007  runTime.caseName() = proc0CaseName;
3008  }
3009  }
3010  Pstream::scatter(masterInstDir);
3011 
3012  // Check who has a mesh
3013  const fileName meshPath =
3014  runTime.path()/masterInstDir/meshSubDir/"faces";
3015 
3016  Info<< "Checking for mesh in " << meshPath << nl << endl;
3017 
3018 
3019  boolList haveMesh(Pstream::nProcs(), false);
3020  haveMesh[Pstream::myProcNo()] = isFile(meshPath);
3021  Pstream::gatherList(haveMesh);
3022  Pstream::scatterList(haveMesh);
3023  Info<< "Per processor mesh availability:" << nl
3024  << " " << flatOutput(haveMesh) << nl << endl;
3025 
3026  // Load mesh (or create dummy one)
3027  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
3028 
3029  if (Pstream::master() && decompose)
3030  {
3031  Info<< "Setting caseName to " << baseRunTime.caseName()
3032  << " to read undecomposed mesh" << endl;
3033  runTime.caseName() = baseRunTime.caseName();
3034  }
3035 
3037  (
3038  IOobject
3039  (
3040  regionName,
3041  masterInstDir,
3042  runTime,
3044  )
3045  );
3046 
3047  if (Pstream::master() && decompose)
3048  {
3049  Info<< "Restoring caseName to " << proc0CaseName << endl;
3050  runTime.caseName() = proc0CaseName;
3051  }
3052 
3053  fvMesh& mesh = meshPtr();
3054 
3055 
3056  const label nOldCells = mesh.nCells();
3057  //Pout<< "Loaded mesh : nCells:" << nOldCells
3058  // << " nPatches:" << mesh.boundaryMesh().size() << endl;
3059 
3060 
3061  // Global matching tolerance
3062  const scalar tolDim = getMergeDistance
3063  (
3064  args,
3065  runTime,
3066  mesh.bounds()
3067  );
3068 
3069  // Determine decomposition
3070  // ~~~~~~~~~~~~~~~~~~~~~~~
3071 
3072  label nDestProcs;
3073  labelList finalDecomp;
3074  determineDecomposition
3075  (
3076  baseRunTime,
3077  decompDictFile,
3078  decompose,
3079  proc0CaseName,
3080  mesh,
3081  writeCellDist,
3082 
3083  nDestProcs,
3084  finalDecomp
3085  );
3086 
3087  if (dryrun)
3088  {
3089  if (!Pstream::master() && !haveMesh[Pstream::myProcNo()])
3090  {
3091  // Remove dummy mesh created by loadOrCreateMesh
3092  const bool oldParRun = Pstream::parRun(false);
3093  mesh.removeFiles();
3094  rmDir(mesh.objectRegistry::objectPath());
3095  Pstream::parRun(oldParRun); // Restore parallel state
3096  }
3097  continue;
3098  }
3099 
3100 
3101 
3103  List<wordList> fieldNames;
3104 
3105  // Detect lagrangian fields
3106  if (Pstream::master() && decompose)
3107  {
3108  runTime.caseName() = baseRunTime.caseName();
3109  }
3111  (
3112  mesh,
3113  cloudNames,
3114  fieldNames
3115  );
3116 
3117  // Read lagrangian fields and store on cloud (objectRegistry)
3119  (
3120  cloudNames.size()
3121  );
3122  readLagrangian
3123  (
3124  mesh,
3125  cloudNames,
3126  selectedLagrangianFields,
3127  clouds
3128  );
3129  if (Pstream::master() && decompose)
3130  {
3131  runTime.caseName() = proc0CaseName;
3132  }
3133 
3134 
3135  // Load fields, do all distribution (mesh and fields - but not
3136  // lagrangian fields; these are done later)
3137  autoPtr<mapDistributePolyMesh> distMap = redistributeAndWrite
3138  (
3139  baseRunTime,
3140  tolDim,
3141  haveMesh,
3142  meshSubDir,
3143  true, // read fields
3144  decompose, // decompose, i.e. read from undecomposed case
3145  false, // no reconstruction
3146  overwrite,
3147  proc0CaseName,
3148  nDestProcs,
3149  finalDecomp,
3150  masterInstDir,
3151  mesh
3152  );
3153 
3154 
3155  // Redistribute any clouds
3156  redistributeLagrangian
3157  (
3158  lagrangianReconstructorPtr,
3159  mesh,
3160  nOldCells,
3161  distMap(),
3162  clouds
3163  );
3164 
3165 
3166  // Copy any uniform data
3167  const fileName uniformDir("uniform");
3168  if (isDir(baseRunTime.timePath()/uniformDir))
3169  {
3170  Info<< "Detected additional non-decomposed files in "
3171  << baseRunTime.timePath()/uniformDir << endl;
3172  cp
3173  (
3174  baseRunTime.timePath()/uniformDir,
3175  runTime.timePath()/uniformDir
3176  );
3177  }
3178  }
3179  }
3180 
3181 
3182  Info<< "End\n" << endl;
3183 
3184  return 0;
3185 }
3186 
3187 
3188 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::TimePaths::findTimes
static instantList findTimes(const fileName &directory, const word &constantName="constant")
Search a given directory for valid time directories.
Definition: TimePaths.C:140
Foam::IOobject::NO_WRITE
Definition: IOobject.H:130
Foam::autoPtr::New
static autoPtr< T > New(Args &&... args)
Construct autoPtr of T with forwarding arguments.
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:71
volFields.H
runTime
engineTime & runTime
Definition: createEngineTime.H:13
Foam::mapDistributeBase::subMap
const labelListList & subMap() const
From subsetted data back to original data.
Definition: mapDistributeBase.H:289
Foam::autoPtr::reset
void reset(T *p=nullptr) noexcept
Delete managed object and set to new given pointer.
Definition: autoPtrI.H:109
Foam::fvc::reconstruct
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> reconstruct(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcReconstruct.C:56
Foam::cloud::prefix
static const word prefix
The prefix to local: lagrangian.
Definition: cloud.H:87
Foam::UPstream::commsTypes::blocking
Foam::parFvFieldReconstructor::reconstructFvVolumeInternalFields
label reconstructFvVolumeInternalFields(const IOobjectList &objects, const wordRes &selectedFields=wordRes()) const
Read, reconstruct and write all/selected volume internal fields.
Foam::TimePaths::globalCaseName
const fileName & globalCaseName() const
Return global case name.
Definition: TimePathsI.H:48
Foam::Tensor< scalar >
Foam::mapDistributeBase::constructHasFlip
bool constructHasFlip() const
Does constructMap include a sign.
Definition: mapDistributeBase.H:325
Foam::volTensorField
GeometricField< tensor, fvPatchField, volMesh > volTensorField
Definition: volFieldsFwd.H:66
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:104
faceProcAddressing
PtrList< labelIOList > & faceProcAddressing
Definition: checkFaceAddressingComp.H:9
Foam::SymmTensor< scalar >
meshPtr
Foam::autoPtr< Foam::fvMesh > meshPtr(nullptr)
Foam::UPstream::commsTypes::nonBlocking
Foam::IOobject::AUTO_WRITE
Definition: IOobject.H:129
Foam::boundBox::mag
scalar mag() const
The magnitude of the bounding box span.
Definition: boundBoxI.H:133
Foam::dimless
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
Definition: dimensionSets.H:50
Foam::UPstream::masterNo
static constexpr int masterNo() noexcept
Process index of the master (always 0)
Definition: UPstream.H:452
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::faceMap
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
Definition: blockMeshMergeTopological.C:94
cloudName
const word cloudName(propsDict.get< word >("cloud"))
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:62
regionProperties.H
Foam::fileOperation
An encapsulation of filesystem-related operations.
Definition: fileOperation.H:77
Foam::parLagrangianRedistributor::redistributeStoredFields
label redistributeStoredFields(const mapDistributeBase &map, passivePositionParticleCloud &cloud) const
Redistribute and write stored lagrangian fields.
Foam::fileName
A class for handling file names.
Definition: fileName.H:69
Foam::IOField
A primitive field of type <T> with automated input and output.
Definition: foamVtkLagrangianWriter.H:61
Foam::fvMesh::write
virtual bool write(const bool valid=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:1027
Foam::returnReduce
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Definition: PstreamReduceOps.H:94
Foam::bitSet
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Definition: bitSet.H:63
processorFvPatchField.H
Foam::polyBoundaryMesh
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Definition: polyBoundaryMesh.H:62
Foam::polyMesh::POINTS_MOVED
Definition: polyMesh.H:93
Foam::polyMesh::defaultRegion
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:318
Foam::argList::getOrDefault
T getOrDefault(const word &optName, const T &deflt) const
Get a value from the named option if present, or return default.
Definition: argListI.H:280
Foam::processorFvPatchField
This boundary condition enables processor communication across patches.
Definition: processorFvPatchField.H:66
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::fvMeshSubset
Given the original mesh and the list of selected cells, it creates the mesh consisting only of the de...
Definition: fvMeshSubset.H:73
topoSet.H
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::DynamicList< word >
Foam::globalMeshData::processorPatches
const labelList & processorPatches() const
Return list of processor patch labels.
Definition: globalMeshData.H:381
Foam::dictionary::found
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Definition: dictionary.C:364
unmappedPassivePositionParticleCloud.H
Foam::OPstream
Output inter-processor communications stream.
Definition: OPstream.H:52
Foam::surfaceInterpolation::geometry
virtual const fvGeometryScheme & geometry() const
Return reference to geometry calculation scheme.
Definition: surfaceInterpolation.C:80
Foam::Pstream::scatterList
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
Definition: gatherScatterList.C:215
Foam::Time::writeFormat
IOstream::streamFormat writeFormat() const
The write stream format.
Definition: Time.H:387
Foam::mapDistributePolyMesh::distributeCellData
void distributeCellData(List< T > &lst) const
Distribute list of cell data.
Definition: mapDistributePolyMesh.H:256
addOverwriteOption.H
Foam::primitiveMesh::nFaces
label nFaces() const
Number of mesh faces.
Definition: primitiveMeshI.H:90
globalIndex.H
IOmapDistributePolyMesh.H
Foam::fileName::name
static std::string name(const std::string &str)
Return basename (part beyond last /), including its extension.
Definition: fileNameI.H:209
Foam::polyMesh::meshSubDir
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:321
Foam::UPstream::parRun
static bool & parRun()
Test if this a parallel run, or allow modify access.
Definition: UPstream.H:434
Foam::regionProperties::names
wordList names() const
The region names. Sorted by region type.
Definition: regionProperties.C:89
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::argList::addNote
static void addNote(const string &note)
Add extra notes for the usage information.
Definition: argList.C:426
Foam::Time::functionObjects
const functionObjectList & functionObjects() const
Return the list of function objects.
Definition: Time.H:499
Foam::timeSelector::select
instantList select(const instantList &times) const
Select a list of Time values that are within the ranges.
Definition: timeSelector.C:95
Foam::polyMesh::facesInstance
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:852
Foam::mapDistributePolyMesh::transfer
void transfer(mapDistributePolyMesh &map)
Transfer the contents of the argument and annul the argument.
Definition: mapDistributePolyMesh.C:193
Foam::UPstream::waitRequests
static void waitRequests(const label start=0)
Wait until all requests (from start onwards) have finished.
Definition: UPstream.C:234
Foam::UPstream::master
static bool master(const label communicator=worldComm)
Am I the master process.
Definition: UPstream.H:458
Foam::combineReduce
void combineReduce(const List< UPstream::commsStruct > &comms, T &Value, const CombineOp &cop, const int tag, const label comm)
Definition: PstreamCombineReduceOps.H:54
Foam::processorPolyPatch::neighbProcNo
int neighbProcNo() const
Return neighbour processor number.
Definition: processorPolyPatch.H:274
Foam::Pstream::scatter
static void scatter(const List< commsStruct > &comms, T &Value, const int tag, const label comm)
Scatter data. Distribute without modification. Reverse of gather.
Definition: gatherScatter.C:150
Foam::argList
Extract command arguments and options from the supplied argc and argv parameters.
Definition: argList.H:123
Foam::unmappedPassivePositionParticleCloud
passivePositionParticleCloud but with autoMap and writing disabled. Only used for its objectRegistry ...
Definition: unmappedPassivePositionParticleCloud.H:52
Foam::isFile
bool isFile(const fileName &name, const bool checkGzip=true, const bool followLink=true)
Does the name exist as a FILE in the file system?
Definition: MSwindows.C:658
Foam::sigFpe::unset
static void unset(bool verbose=false)
Deactivate SIGFPE signal handler and NaN memory initialisation.
Definition: sigFpe.C:204
Foam::fileHandler
const fileOperation & fileHandler()
Get current file handler.
Definition: fileOperation.C:1354
Foam::rm
bool rm(const fileName &file)
Remove a file (or its gz equivalent), returning true if successful.
Definition: MSwindows.C:994
IOobjectList.H
Foam::polyMesh::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:444
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:350
surfaceFields.H
Foam::surfaceFields.
Foam::PtrList::set
const T * set(const label i) const
Return const pointer to element (can be nullptr),.
Definition: PtrList.H:138
Foam::UPstream::defaultCommsType
static commsTypes defaultCommsType
Default commsType.
Definition: UPstream.H:283
Foam::CompactIOField
A Field of objects of type <T> with automated input and output using a compact storage....
Definition: CompactIOField.H:53
decompositionMethod.H
Foam::dictionary::get
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:81
Foam::HashSet< word >
parLagrangianRedistributor.H
Foam::meshRefinement::removeFiles
static void removeFiles(const polyMesh &)
Helper: remove all relevant files from mesh instance.
Definition: meshRefinement.C:3416
Foam::decompositionModel
MeshObject wrapper of decompositionMethod.
Definition: decompositionModel.H:57
Foam::polyBoundaryMesh::start
label start() const
The start label of the boundary faces in the polyMesh face list.
Definition: polyBoundaryMesh.C:638
Foam::argList::rootPath
const fileName & rootPath() const
Return root path.
Definition: argListI.H:63
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::polyBoundaryMesh::names
wordList names() const
Return a list of patch names.
Definition: polyBoundaryMesh.C:593
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
foamDlOpenLibs.H
Foam::mapDistributePolyMesh::cellMap
const mapDistribute & cellMap() const
Cell distribute map.
Definition: mapDistributePolyMesh.H:223
Foam::fvMesh::readUpdate
virtual readUpdateState readUpdate()
Update the mesh based on the mesh files saved in time.
Definition: fvMesh.C:643
Foam::parLagrangianRedistributor::findClouds
static void findClouds(const fvMesh &, wordList &cloudNames, List< wordList > &objectNames)
Find all clouds (on all processors) and for each cloud all.
Foam::polyMesh::TOPO_PATCH_CHANGE
Definition: polyMesh.H:95
Foam::parFvFieldReconstructor::reconstructFvSurfaceFields
label reconstructFvSurfaceFields(const IOobjectList &objects, const wordRes &selectedFields=wordRes()) const
Read, reconstruct and write all/selected surface fields.
Foam::polyMesh::pointsInstance
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:846
regionName
Foam::word regionName
Definition: createNamedDynamicFvMesh.H:1
Foam::decompositionMethod::parallelAware
virtual bool parallelAware() const =0
Is method parallel aware?
Foam::argList::noFunctionObjects
static void noFunctionObjects(bool addWithOption=false)
Remove '-noFunctionObjects' option and ignore any occurrences.
Definition: argList.C:467
Foam::UPstream::subProcs
static rangeType subProcs(const label communicator=worldComm)
Range of process indices for sub-processes.
Definition: UPstream.H:516
Foam::primitiveMesh::nCells
label nCells() const
Number of mesh cells.
Definition: primitiveMeshI.H:96
Foam::Field< scalar >
Foam::IOobject::writeOpt
writeOption writeOpt() const
The write option.
Definition: IOobjectI.H:177
Foam::volSymmTensorField
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:65
Foam::functionObjectList::off
void off()
Switch the function objects off.
Definition: functionObjectList.C:652
Foam::Info
messageStream Info
Information stream (uses stdout - output is on the master only)
Foam::mapDistributeBase::distribute
static void distribute(const Pstream::commsTypes commsType, const List< labelPair > &schedule, const label constructSize, const labelListList &subMap, const bool subHasFlip, const labelListList &constructMap, const bool constructHasFlip, List< T > &, const negateOp &negOp, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Distribute data. Note:schedule only used for.
Definition: mapDistributeBaseTemplates.C:122
Foam::name
word name(const complex &c)
Return string representation of complex.
Definition: complex.C:76
Foam::parLagrangianRedistributor
Lagrangian field redistributor.
Definition: parLagrangianRedistributor.H:61
Foam::Time::timePath
fileName timePath() const
Return current time path.
Definition: Time.H:375
Foam::fvMeshTools::newMesh
static autoPtr< fvMesh > newMesh(const IOobject &io, const bool masterOnlyReading)
Read mesh or create dummy mesh (0 cells, >0 patches). Works in two.
Definition: fvMeshTools.C:424
Foam::polyMesh::removeFiles
void removeFiles(const fileName &instanceDir) const
Remove all files from mesh instance.
Definition: polyMesh.C:1325
Foam::DynamicList::append
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
Definition: DynamicListI.H:474
Foam::mapDistributePolyMesh::faceMap
const mapDistribute & faceMap() const
Face distribute map.
Definition: mapDistributePolyMesh.H:217
argList.H
Foam::mapDistribute
Class containing processor-to-processor mapping information.
Definition: mapDistribute.H:163
Foam::decompositionModel::decomposer
decompositionMethod & decomposer() const
Definition: decompositionModel.H:122
init
mesh init(true)
Foam::IOobject::READ_IF_PRESENT
Definition: IOobject.H:122
Foam::mapDistributePolyMesh::patchMap
const mapDistribute & patchMap() const
Patch distribute map.
Definition: mapDistributePolyMesh.H:229
Foam::decompositionMethod::nDomains
static label nDomains(const dictionary &decompDict)
Return the numberOfSubdomains entry from the dictionary.
Definition: decompositionMethod.C:62
Foam::TimePaths::times
instantList times() const
Search the case for valid time directories.
Definition: TimePaths.C:149
addRegionOption.H
Foam::parLagrangianRedistributor::redistributeLagrangianPositions
autoPtr< mapDistributeBase > redistributeLagrangianPositions(passivePositionParticleCloud &cloud) const
Redistribute and write lagrangian positions.
Foam::CompactIOList< face, label >
Foam::andOp
Definition: ops.H:233
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:43
Foam::processorPolyPatch
Neighbour processor patch.
Definition: processorPolyPatch.H:58
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:62
Foam::polyMesh::TOPO_CHANGE
Definition: polyMesh.H:94
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Foam::regionProperties
Simple class to hold region information for coupled region simulations.
Definition: regionProperties.H:60
Foam::pow
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Definition: dimensionedScalar.C:75
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
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::hexRef8Data
Various for reading/decomposing/reconstructing/distributing refinement data.
Definition: hexRef8Data.H:60
Foam::UPstream::commsTypeNames
static const Enum< commsTypes > commsTypeNames
Names of the communication types.
Definition: UPstream.H:77
Foam::timeSelector::selectIfPresent
static instantList selectIfPresent(Time &runTime, const argList &args)
Definition: timeSelector.C:272
Foam::mapDistributePolyMesh::oldPatchSizes
const labelList & oldPatchSizes() const
List of the old patch sizes.
Definition: mapDistributePolyMesh.H:193
Foam::UPstream::commsTypes::scheduled
forAllIters
#define forAllIters(container, iter)
Iterate across all elements in the container object.
Definition: stdFoam.H:223
timeName
word timeName
Definition: getTimeIndex.H:3
Foam::argList::path
fileName path() const
Return the full path to the (processor local) case.
Definition: argListI.H:81
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:121
Foam::primitiveMesh::nBoundaryFaces
label nBoundaryFaces() const
Number of boundary faces (== nFaces - nInternalFaces)
Definition: primitiveMeshI.H:84
Foam::Time::controlDict
const dictionary & controlDict() const
Return read access to the controlDict dictionary.
Definition: Time.H:364
Foam::eqOp
Definition: ops.H:71
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
meshRefinement.H
Foam::decompositionMethod
Abstract base class for domain decomposition.
Definition: decompositionMethod.H:51
Foam::mapDistributeBase::subHasFlip
bool subHasFlip() const
Does subMap include a sign.
Definition: mapDistributeBase.H:313
hexRef8Data.H
Foam::GeoMesh
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:48
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:83
fvMesh.H
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::flatOutput
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
Definition: FlatOutput.H:217
cloudNames
const wordList cloudNames(cloudFields.sortedToc())
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::globalIndex
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:68
Foam::ListOps::uniqueEqOp
List helper to append y unique elements onto the end of x.
Definition: ListOps.H:577
PstreamReduceOps.H
Inter-processor communication reduction functions.
Foam::SphericalTensor< scalar >
Foam::IOobjectList
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:55
Foam::Time::findInstance
word findInstance(const fileName &dir, const word &name=word::null, const IOobject::readOption rOpt=IOobject::MUST_READ, const word &stopInstance=word::null) const
Definition: Time.C:797
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::readFields
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
Definition: ReadFieldsTemplates.C:312
Foam::polyMesh::readUpdateState
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:90
Foam::HashTable
A HashTable similar to std::unordered_map.
Definition: HashTable.H:105
Foam::polyBoundaryMesh::nNonProcessor
label nNonProcessor() const
The number of patches before the first processor patch.
Definition: polyBoundaryMesh.C:573
Foam::volVectorField
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:62
Foam::polyMesh::bounds
const boundBox & bounds() const
Return mesh bounding box.
Definition: polyMesh.H:450
Foam::fvMesh::boundary
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:679
Foam::mapDistributeBase::constructSize
label constructSize() const
Constructed data size.
Definition: mapDistributeBase.H:277
Foam::decompositionModel::New
static const decompositionModel & New(const polyMesh &mesh, const fileName &decompDictFile="")
Read (optionally from absolute path) and register on mesh.
Definition: decompositionModel.C:116
Foam::argList::addBoolOption
static void addBoolOption(const word &optName, const string &usage="", bool advanced=false)
Add a bool option to validOptions with usage information.
Definition: argList.C:338
Time.H
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: HashPtrTable.H:53
Foam::fvMeshSubset::interpolate
static tmp< GeometricField< Type, fvPatchField, volMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &, const fvMesh &sMesh, const labelUList &patchMap, const labelUList &cellMap, const labelUList &faceMap, const bool allowUnmapped=false)
Map volume field. Optionally allow unmapped faces not to produce.
Foam::TimePaths::system
const word & system() const
Return system name.
Definition: TimePathsI.H:94
Foam::IOstreamOption::ASCII
"ascii" (normal default)
Definition: IOstreamOption.H:72
Foam::UPstream::msgType
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:541
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:381
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.
Definition: gatherScatterList.C:52
Foam::UPstream::allProcs
static rangeType allProcs(const label communicator=worldComm)
Range of process indices for all processes.
Definition: UPstream.H:509
Foam::Time::caseName
const fileName & caseName() const
Return case name.
Definition: TimePathsI.H:54
Foam::parFvFieldReconstructor::reconstructFvVolumeFields
label reconstructFvVolumeFields(const IOobjectList &objects, const wordRes &selectedFields=wordRes()) const
Read, reconstruct and write all/selected volume fields.
Foam::UPstream::nRequests
static label nRequests()
Get number of outstanding requests.
Definition: UPstream.C:224
fvMeshDistribute.H
Foam::cp
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy the source to the destination (recursively if necessary).
Definition: MSwindows.C:802
Foam::UPstream::myProcNo
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:464
Foam::nl
constexpr char nl
Definition: Ostream.H:385
Foam::parLagrangianRedistributor::readFields
static label readFields(const passivePositionParticleCloud &cloud, const IOobjectList &objects, const wordRes &selectedFields=wordRes())
Read and store all fields of a cloud.
Foam::Time::path
fileName path() const
Return path.
Definition: Time.H:358
Foam::IOstream::defaultPrecision
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:333
fvMeshTools.H
Foam::rmDir
bool rmDir(const fileName &directory, const bool silent=false)
Remove a directory and its contents (optionally silencing warnings)
Definition: MSwindows.C:1018
Foam::lagrangianReconstructor
Reconstructor for lagrangian positions and fields.
Definition: lagrangianReconstructor.H:56
Foam::timeSelector::addOptions
static void addOptions(const bool constant=true, const bool withZero=false)
Add timeSelector options to argList::validOptions.
Definition: timeSelector.C:108
Foam::Vector< scalar >
Foam::Time::rootPath
const fileName & rootPath() const
Return root path.
Definition: TimePathsI.H:42
Foam::List< labelList >
Foam::IOmapDistributePolyMesh
IOmapDistributePolyMesh is derived from mapDistributePolyMesh and IOobject to give the mapDistributeP...
Definition: IOmapDistributePolyMesh.H:53
Foam::decompositionMethod::decompose
virtual labelList decompose(const pointField &points, const scalarField &pointWeights) const
Return for every coordinate the wanted processor number.
Definition: decompositionMethod.H:244
Foam::Time::setTime
virtual void setTime(const Time &t)
Reset the time and time-index to those of the given time.
Definition: Time.C:1003
Foam::UList< label >
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
Foam::primitiveMesh::nPoints
label nPoints() const
Number of mesh points.
Definition: primitiveMeshI.H:37
Foam::identity
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
Definition: labelList.C:38
Foam::IOList< label >
Foam::wordRes
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:51
Foam::argList::globalPath
fileName globalPath() const
Return the full path to the global case.
Definition: argListI.H:87
timeSelector.H
createTime.H
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:65
Foam::autoPtr::clear
void clear() noexcept
Same as reset(nullptr)
Definition: autoPtr.H:181
Foam::boundBox
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:63
Foam::mapDistributeBase
Class containing processor-to-processor mapping information.
Definition: mapDistributeBase.H:103
Foam::mapDistribute::reverseDistribute
void reverseDistribute(const label constructSize, List< T > &, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Reverse distribute data using default commsType.
Definition: mapDistributeTemplates.C:182
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:275
Foam::IPstream
Input inter-processor communications stream.
Definition: IPstream.H:52
Foam::mapDistributePolyMesh::pointMap
const mapDistribute & pointMap() const
Point distribute map.
Definition: mapDistributePolyMesh.H:211
decompositionModel.H
cyclicACMIFvPatch.H
Foam::mapDistributePolyMesh::nOldFaces
label nOldFaces() const
Number of faces in mesh before distribution.
Definition: mapDistributePolyMesh.H:181
Foam::instant
An instant of time. Contains the time value and name.
Definition: instant.H:52
Foam::volSphericalTensorField
GeometricField< sphericalTensor, fvPatchField, volMesh > volSphericalTensorField
Definition: volFieldsFwd.H:64
Foam::mapDistributePolyMesh
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Definition: mapDistributePolyMesh.H:66
Foam::polyMesh::globalData
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1295
Foam::mapDistributePolyMesh::nOldPoints
label nOldPoints() const
Number of points in mesh before distribution.
Definition: mapDistributePolyMesh.H:175
Foam::fvGeometryScheme::New
static tmp< fvGeometryScheme > New(const fvMesh &mesh, const dictionary &dict, const word &defaultScheme)
Return new tmp interpolation scheme.
Definition: fvGeometryScheme.C:44
Foam::IOobject::objectPath
fileName objectPath() const
The complete path + object name.
Definition: IOobjectI.H:209
Foam::mapDistributePolyMesh::nOldCells
label nOldCells() const
Number of cells in mesh before distribution.
Definition: mapDistributePolyMesh.H:187
Foam::List::setSize
void setSize(const label newSize)
Alias for resize(const label)
Definition: ListI.H:146
Foam::TimePaths::constant
const word & constant() const
Return constant name.
Definition: TimePathsI.H:88
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::mapDistributeBase::constructMap
const labelListList & constructMap() const
From subsetted data to new reconstructed data.
Definition: mapDistributeBase.H:301
Foam::IOobjectList::lookupClass
IOobjectList lookupClass(const char *clsName) const
The list of IOobjects with the given headerClassName.
Definition: IOobjectList.C:290
Foam::argList::addOption
static void addOption(const word &optName, const string &param="", const string &usage="", bool advanced=false)
Add an option to validOptions with usage information.
Definition: argList.C:349
Foam::IOobject::NO_READ
Definition: IOobject.H:123
args
Foam::argList args(argc, argv)
Foam::mkDir
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
Definition: MSwindows.C:507
parFvFieldReconstructor.H
Foam::mapDistributePolyMesh::distributePointData
void distributePointData(List< T > &lst) const
Distribute list of point data.
Definition: mapDistributePolyMesh.H:242
basicFvGeometryScheme.H
Foam::topoSet::removeFiles
static void removeFiles(const polyMesh &)
Helper: remove all sets files from mesh instance.
Definition: topoSet.C:641
DebugVar
#define DebugVar(var)
Report a variable name and value.
Definition: messageStream.H:381
zeroGradientFvPatchFields.H
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:303
Foam::fvMeshDistribute
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
Definition: fvMeshDistribute.H:73
pointFields.H
Foam::polyMesh::setInstance
void setInstance(const fileName &instance, const IOobject::writeOption wOpt=IOobject::AUTO_WRITE)
Set the instance for mesh files.
Definition: polyMeshIO.C:36
fields
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:97
Foam::UPstream::nProcs
static label nProcs(const label communicator=worldComm)
Number of processes in parallel run, and 1 for serial run.
Definition: UPstream.H:446
Foam::parFvFieldReconstructor
Finite volume reconstructor for volume and surface fields.
Definition: parFvFieldReconstructor.H:61
Foam::isDir
bool isDir(const fileName &name, const bool followLink=true)
Does the name exist as a DIRECTORY in the file system?
Definition: MSwindows.C:643
Foam::loadOrCreateMesh
autoPtr< fvMesh > loadOrCreateMesh(const IOobject &io)
Load (if it exists) or create zero cell mesh given an IOobject:
Foam::argList::found
bool found(const word &optName) const
Return true if the named option is found.
Definition: argListI.H:151
loadOrCreateMesh.H
Load or create (0 size) a mesh. Used in distributing meshes to a larger number of processors.
Foam::flipOp
Functor to negate primitives. Dummy for most other types.
Definition: flipOp.H:53
Foam::IOobject::MUST_READ
Definition: IOobject.H:120
Foam::flipLabelOp
Negate integer values.
Definition: flipOp.H:85
Foam::polyMesh::tetBasePtIs
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:892