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-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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
27Application
28 redistributePar
29
30Group
31 grpParallelUtilities
32
33Description
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
43Usage
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 "fieldsDistributor.H"
82#include "decompositionMethod.H"
83#include "decompositionModel.H"
84#include "timeSelector.H"
85#include "PstreamReduceOps.H"
86#include "volFields.H"
87#include "surfaceFields.H"
89#include "IOobjectList.H"
90#include "globalIndex.H"
91#include "loadOrCreateMesh.H"
94#include "topoSet.H"
95#include "regionProperties.H"
96
99#include "hexRef8Data.H"
100#include "meshRefinement.H"
101#include "pointFields.H"
102
103#include "faMeshSubset.H"
104#include "faMeshTools.H"
105#include "faMeshDistributor.H"
107
109
110#include "cyclicACMIFvPatch.H"
114
115using namespace Foam;
116
117// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
118
119const int debug(::Foam::debug::debugSwitch("redistributePar", 0));
120
121void createTimeDirs(const fileName& path)
122{
123 // Get current set of local processor's time directories. Uses
124 // fileHandler
125 const instantList localTimeDirs(Time::findTimes(path, "constant"));
126
127 instantList masterTimeDirs;
128 if (Pstream::master())
129 {
130 //const bool oldParRun = Pstream::parRun(false);
131 //timeDirs = Time::findTimes(path, "constant");
132 //Pstream::parRun(oldParRun); // Restore parallel state
133 masterTimeDirs = localTimeDirs;
134 }
135 Pstream::scatter(masterTimeDirs);
136 //DebugVar(masterTimeDirs);
137 //DebugVar(localTimeDirs);
138
139 // Sync any cached times (e.g. masterUncollatedFileOperation::times_)
140 // since only master would have done the findTimes
141 for (const instant& t : masterTimeDirs)
142 {
143 if (!localTimeDirs.found(t))
144 {
145 const fileName timePath(path/t.name());
146
147 //Pout<< "Time:" << t << nl
148 // << " raw :" << timePath << nl
149 // << endl;
150 Foam::mkDir(timePath);
151 }
152 }
153
154 // Just to make sure remove all state and re-scan
155 fileHandler().flush();
156 (void)Time::findTimes(path, "constant");
157}
158
159
160void copyUniform
161(
162 const fileOperation& fh,
163 const bool decompose,
164 const bool reconstruct,
165 const word& readTimeName,
166 const objectRegistry& readDb,
167 const objectRegistry& writeDb
168)
169{
170 // Detect uniform/ at original database + time
171 const IOobject readIO("uniform", readTimeName, readDb);
172 const fileName readPath
173 (
174 fh.dirPath
175 (
176 false, // local directory
177 readIO,
178 false // do not search in time
179 )
180 );
181 //if (Pstream::master() && !readPath.empty())
182 if (!readPath.empty())
183 {
184 Info<< "Detected additional non-decomposed files in "
185 << readPath << endl;
186
187 // readPath: searching is the same for all file handlers. Typical:
188 // <case>/0.1/uniform (parent dir, decompose mode)
189 // <case>/processor1/0.1/uniform (redistribute/reconstruct mode)
190 // <case>/processors2/0.1/uniform ,,
191 // writePath:
192 // uncollated : <case>/0.1/uniform (reconstruct mode). Should only
193 // be done by master
194 // uncollated : <case>/processorXXX/0.1/uniform. Should be done by all.
195 // collated : <case>/processors2/0.1/uniform. Should be done by
196 // local master only.
197
198 // See what local directory
199 const IOobject writeIO("uniform", writeDb.time().timeName(), writeDb);
200 const fileName writePath
201 (
202 fh.objectPath
203 (
204 writeIO,
206 )
207 );
208 // Do we already have this directory?
209 const fileName currentPath(fh.dirPath(false, writeIO, false));
210
211 if (::debug)
212 {
213 Pout<< " readPath :" << readPath << endl;
214 Pout<< " writePath :" << writePath << endl;
215 Pout<< " currentPath:" << currentPath << endl;
216 }
217
218 if (readPath == writePath)
219 {
220 return;
221 }
222
223 if (currentPath.empty())
224 {
225 if (decompose)
226 {
227 // All processors copy to destination
228 fh.cp(readPath, writePath);
229 }
230 else if (reconstruct)
231 {
232 // Only master
233 if (Pstream::master())
234 {
235 const bool oldParRun = Pstream::parRun(false);
236 fh.cp(readPath, writePath);
237 Pstream::parRun(oldParRun);
238 }
239 }
240 else
241 {
242 // Redistribute. If same destination path do only on master,
243 // if different path do on all processors. For now check
244 // if collated file handler only. tbd.
245 if (isA<fileOperations::collatedFileOperation>(fh))
246 {
247 // Collated
248 if (Pstream::master())
249 {
250 const bool oldParRun = Pstream::parRun(false);
251 fh.cp(readPath, writePath);
252 Pstream::parRun(oldParRun);
253 }
254 }
255 else
256 {
257 // Assume uncollated
258 fh.cp(readPath, writePath);
259 }
260 }
261 }
262 }
263}
264
265
266void printMeshData(const polyMesh& mesh)
267{
268 // Collect all data on master
269
270 labelListList patchNeiProcNo(Pstream::nProcs());
271 labelListList patchSize(Pstream::nProcs());
272 const labelList& pPatches = mesh.globalData().processorPatches();
273 patchNeiProcNo[Pstream::myProcNo()].setSize(pPatches.size());
274 patchSize[Pstream::myProcNo()].setSize(pPatches.size());
275 forAll(pPatches, i)
276 {
277 const processorPolyPatch& ppp = refCast<const processorPolyPatch>
278 (
279 mesh.boundaryMesh()[pPatches[i]]
280 );
281 patchNeiProcNo[Pstream::myProcNo()][i] = ppp.neighbProcNo();
282 patchSize[Pstream::myProcNo()][i] = ppp.size();
283 }
284 Pstream::gatherList(patchNeiProcNo);
285 Pstream::gatherList(patchSize);
286
287
288 // Print stats
289
290 const globalIndex globalCells(mesh.nCells());
291 const globalIndex globalBoundaryFaces(mesh.nBoundaryFaces());
292
293 label maxProcCells = 0;
294 label maxProcFaces = 0;
295 label totProcFaces = 0;
296 label maxProcPatches = 0;
297 label totProcPatches = 0;
298
299 for (const int proci : Pstream::allProcs())
300 {
301 const label nLocalCells = globalCells.localSize(proci);
302 const label nBndFaces = globalBoundaryFaces.localSize(proci);
303
304 Info<< nl
305 << "Processor " << proci;
306
307 if (!nLocalCells)
308 {
309 Info<< " (empty)" << endl;
310 continue;
311 }
312 else
313 {
314 Info<< nl
315 << " Number of cells = " << nLocalCells << endl;
316 }
317
318 label nProcFaces = 0;
319 const labelList& nei = patchNeiProcNo[proci];
320
321 forAll(patchNeiProcNo[proci], i)
322 {
323 Info<< " Number of faces shared with processor "
324 << patchNeiProcNo[proci][i] << " = "
325 << patchSize[proci][i] << nl;
326
327 nProcFaces += patchSize[proci][i];
328 }
329
330 {
331 Info<< " Number of processor patches = " << nei.size() << nl
332 << " Number of processor faces = " << nProcFaces << nl
333 << " Number of boundary faces = "
334 << nBndFaces-nProcFaces << endl;
335 }
336
337 maxProcCells = max(maxProcCells, nLocalCells);
338 totProcFaces += nProcFaces;
339 totProcPatches += nei.size();
340 maxProcFaces = max(maxProcFaces, nProcFaces);
341 maxProcPatches = max(maxProcPatches, nei.size());
342 }
343
344 // Summary stats
345
346 Info<< nl
347 << "Number of processor faces = " << (totProcFaces/2) << nl
348 << "Max number of cells = " << maxProcCells;
349
350 if (maxProcCells != globalCells.totalSize())
351 {
352 scalar avgValue = scalar(globalCells.totalSize())/Pstream::nProcs();
353
354 Info<< " (" << 100.0*(maxProcCells-avgValue)/avgValue
355 << "% above average " << avgValue << ')';
356 }
357 Info<< nl;
358
359 Info<< "Max number of processor patches = " << maxProcPatches;
360 if (totProcPatches)
361 {
362 scalar avgValue = scalar(totProcPatches)/Pstream::nProcs();
363
364 Info<< " (" << 100.0*(maxProcPatches-avgValue)/avgValue
365 << "% above average " << avgValue << ')';
366 }
367 Info<< nl;
368
369 Info<< "Max number of faces between processors = " << maxProcFaces;
370 if (totProcFaces)
371 {
372 scalar avgValue = scalar(totProcFaces)/Pstream::nProcs();
373
374 Info<< " (" << 100.0*(maxProcFaces-avgValue)/avgValue
375 << "% above average " << avgValue << ')';
376 }
377 Info<< nl << endl;
378}
379
380
381// Debugging: write volScalarField with decomposition for post processing.
382void writeDecomposition
383(
384 const word& name,
385 const fvMesh& mesh,
386 const labelUList& decomp
387)
388{
389 // Write the decomposition as labelList for use with 'manual'
390 // decomposition method.
391 labelIOList cellDecomposition
392 (
394 (
395 "cellDecomposition",
396 mesh.facesInstance(), // mesh read from facesInstance
397 mesh,
398 IOobject::NO_READ,
399 IOobject::NO_WRITE,
400 false
401 ),
402 decomp
403 );
404 cellDecomposition.write();
405
406 Info<< "Writing wanted cell distribution to volScalarField " << name
407 << " for postprocessing purposes." << nl << endl;
408
409 volScalarField procCells
410 (
412 (
413 name,
414 mesh.time().timeName(),
415 mesh,
416 IOobject::NO_READ,
417 IOobject::NO_WRITE,
418 false // do not register
419 ),
420 mesh,
421 dimensionedScalar(name, dimless, -1),
423 );
424
425 forAll(procCells, celli)
426 {
427 procCells[celli] = decomp[celli];
428 }
429
430 procCells.correctBoundaryConditions();
431 procCells.write();
432}
433
434
435void determineDecomposition
436(
437 const Time& baseRunTime,
438 const fileName& decompDictFile, // optional location for decomposeParDict
439 const bool decompose, // decompose, i.e. read from undecomposed case
440 const fileName& proc0CaseName,
441 const fvMesh& mesh,
442 const bool writeCellDist,
443
444 label& nDestProcs,
445 labelList& decomp
446)
447{
448 // Read decomposeParDict (on all processors)
449 const decompositionModel& method = decompositionModel::New
450 (
451 mesh,
452 decompDictFile
453 );
454
455 decompositionMethod& decomposer = method.decomposer();
456
457 if (!decomposer.parallelAware())
458 {
460 << "You have selected decomposition method \""
461 << decomposer.type() << "\n"
462 << " which does not synchronise decomposition across"
463 " processor patches.\n"
464 " You might want to select a decomposition method"
465 " that is aware of this. Continuing...." << endl;
466 }
467
468 Time& tm = const_cast<Time&>(mesh.time());
469
470 const bool oldProcCase = tm.processorCase();
471 if (Pstream::master() && decompose)
472 {
473 Info<< "Setting caseName to " << baseRunTime.caseName()
474 << " to read decomposeParDict" << endl;
475 tm.caseName() = baseRunTime.caseName();
476 tm.processorCase(false);
477 }
478
479 scalarField cellWeights;
480 if (method.found("weightField"))
481 {
482 word weightName = method.get<word>("weightField");
483
484 volScalarField weights
485 (
487 (
488 weightName,
489 tm.timeName(),
490 mesh,
491 IOobject::MUST_READ,
492 IOobject::NO_WRITE
493 ),
494 mesh
495 );
496 cellWeights = weights.internalField();
497 }
498
499 nDestProcs = decomposer.nDomains();
500 decomp = decomposer.decompose(mesh, cellWeights);
501
502 if (Pstream::master() && decompose)
503 {
504 Info<< "Restoring caseName" << endl;
505 tm.caseName() = proc0CaseName;
506 tm.processorCase(oldProcCase);
507 }
508
509 // Dump decomposition to volScalarField
510 if (writeCellDist)
511 {
512 // Note: on master make sure to write to processor0
513 if (decompose)
514 {
515 if (Pstream::master())
516 {
517 const bool oldParRun = Pstream::parRun(false);
518
519 Info<< "Setting caseName to " << baseRunTime.caseName()
520 << " to write undecomposed cellDist" << endl;
521
522 tm.caseName() = baseRunTime.caseName();
523 tm.processorCase(false);
524 writeDecomposition("cellDist", mesh, decomp);
525 Info<< "Restoring caseName" << endl;
526 tm.caseName() = proc0CaseName;
527 tm.processorCase(oldProcCase);
528
529 Pstream::parRun(oldParRun);
530 }
531 }
532 else
533 {
534 writeDecomposition("cellDist", mesh, decomp);
535 }
536 }
537}
538
539
540// Variant of GeometricField::correctBoundaryConditions that only
541// evaluates selected patch fields
542template<class CoupledPatchType, class GeoField>
543void correctCoupledBoundaryConditions(fvMesh& mesh)
544{
545 for (GeoField& fld : mesh.sorted<GeoField>())
546 {
547 fld.boundaryFieldRef().template evaluateCoupled<CoupledPatchType>();
548 }
549}
550
551
552// Inplace redistribute mesh and any fields
553autoPtr<mapDistributePolyMesh> redistributeAndWrite
554(
555 autoPtr<fileOperation>&& writeHandler,
556 const Time& baseRunTime,
557 const fileName& proc0CaseName,
558
559 // Controls
560 const bool doReadFields,
561 const bool decompose, // decompose, i.e. read from undecomposed case
562 const bool reconstruct,
563 const bool overwrite,
564
565 // Decomposition information
566 const label nDestProcs,
567 const labelList& decomp,
568
569 // Mesh information
570 const boolList& volMeshOnProc,
571 const fileName& volMeshInstance,
572 fvMesh& mesh
573)
574{
575 Time& runTime = const_cast<Time&>(mesh.time());
576 const bool oldProcCase = runTime.processorCase();
577
579 //Info<< "Before distribution:" << endl;
580 //printMeshData(mesh);
581
582 // Storage of fields
583
584 PtrList<volScalarField> volScalarFields;
585 PtrList<volVectorField> volVectorFields;
586 PtrList<volSphericalTensorField> volSphereTensorFields;
587 PtrList<volSymmTensorField> volSymmTensorFields;
588 PtrList<volTensorField> volTensorFields;
589
590 PtrList<surfaceScalarField> surfScalarFields;
591 PtrList<surfaceVectorField> surfVectorFields;
592 PtrList<surfaceSphericalTensorField> surfSphereTensorFields;
593 PtrList<surfaceSymmTensorField> surfSymmTensorFields;
594 PtrList<surfaceTensorField> surfTensorFields;
595
601
602 PtrList<pointScalarField> pointScalarFields;
603 PtrList<pointVectorField> pointVectorFields;
604 PtrList<pointTensorField> pointTensorFields;
605 PtrList<pointSphericalTensorField> pointSphTensorFields;
606 PtrList<pointSymmTensorField> pointSymmTensorFields;
607
608 // Self-contained pointMesh for reading pointFields
609 const pointMesh oldPointMesh(mesh);
610
611 // Track how many (if any) pointFields are read/mapped
612 label nPointFields = 0;
613
614 parPointFieldDistributor pointDistributor
615 (
616 oldPointMesh, // source mesh
617 false, // savePoints=false (ie, delay until later)
618 false // Do not write
619 );
620
621
622 if (doReadFields)
623 {
624 // Create 0 sized mesh to do all the generation of zero sized
625 // fields on processors that have zero sized meshes. Note that this is
626 // only necessary on master but since polyMesh construction with
627 // Pstream::parRun does parallel comms we have to do it on all
628 // processors
629 autoPtr<fvMeshSubset> subsetterPtr;
630
631 // Missing a volume mesh somewhere?
632 if (volMeshOnProc.found(false))
633 {
634 // A zero-sized mesh with boundaries.
635 // This is used to create zero-sized fields.
636 subsetterPtr.reset(new fvMeshSubset(mesh, zero{}));
637 }
638
639
640 // Get original objects (before incrementing time!)
641 if (Pstream::master() && decompose)
642 {
643 runTime.caseName() = baseRunTime.caseName();
644 runTime.processorCase(false);
645 }
646 IOobjectList objects(mesh, runTime.timeName());
647 if (Pstream::master() && decompose)
648 {
649 runTime.caseName() = proc0CaseName;
650 runTime.processorCase(oldProcCase);
651 }
652
653 Info<< "From time " << runTime.timeName()
654 << " mesh:" << mesh.objectRegistry::objectRelPath()
655 << " have objects:" << objects.names() << endl;
656
657 // We don't want to map the decomposition (mapping already tested when
658 // mapping the cell centre field)
659 auto iter = objects.find("cellDist");
660 if (iter.found())
661 {
662 objects.erase(iter);
663 }
664
665
666 if (Pstream::master() && decompose)
667 {
668 runTime.caseName() = baseRunTime.caseName();
669 runTime.processorCase(false);
670 }
671
672 // Field reading
673
674 #undef doFieldReading
675 #define doFieldReading(Storage) \
676 { \
677 fieldsDistributor::readFields \
678 ( \
679 volMeshOnProc, mesh, subsetterPtr, objects, Storage \
680 ); \
681 }
682
683 // volField
684 doFieldReading(volScalarFields);
685 doFieldReading(volVectorFields);
686 doFieldReading(volSphereTensorFields);
687 doFieldReading(volSymmTensorFields);
688 doFieldReading(volTensorFields);
689
690 // surfaceField
691 doFieldReading(surfScalarFields);
692 doFieldReading(surfVectorFields);
693 doFieldReading(surfSphereTensorFields);
694 doFieldReading(surfSymmTensorFields);
695 doFieldReading(surfTensorFields);
696
697 // Dimensioned internal fields
698 doFieldReading(dimScalarFields);
699 doFieldReading(dimVectorFields);
700 doFieldReading(dimSphereTensorFields);
701 doFieldReading(dimSymmTensorFields);
702 doFieldReading(dimTensorFields);
703
704 // pointFields
705 nPointFields = 0;
706
707 #undef doFieldReading
708 #define doFieldReading(Storage) \
709 { \
710 fieldsDistributor::readFields \
711 ( \
712 volMeshOnProc, oldPointMesh, subsetterPtr, objects, Storage, \
713 true /* (deregister field) */ \
714 ); \
715 nPointFields += Storage.size(); \
716 }
717
718 doFieldReading(pointScalarFields);
719 doFieldReading(pointVectorFields);
720 doFieldReading(pointSphTensorFields);
721 doFieldReading(pointSymmTensorFields);
722 doFieldReading(pointTensorFields);
723 #undef doFieldReading
724
725
726 // Done reading
727
728 if (Pstream::master() && decompose)
729 {
730 runTime.caseName() = proc0CaseName;
731 runTime.processorCase(oldProcCase);
732 }
733 }
734
735 // Save pointMesh information before any topology changes occur!
736 if (nPointFields)
737 {
738 pointDistributor.saveMeshPoints();
739 }
740
741
742 // Mesh distribution engine
743 fvMeshDistribute distributor(mesh);
744
745 // Do all the distribution of mesh and fields
746 autoPtr<mapDistributePolyMesh> distMap = distributor.distribute(decomp);
747
748 // Print some statistics
749 Info<< "After distribution:" << endl;
750 printMeshData(mesh);
751
752 // Get other side of processor boundaries
753 do
754 {
755 #undef doCorrectCoupled
756 #define doCorrectCoupled(FieldType) \
757 correctCoupledBoundaryConditions<processorFvPatch, FieldType>(mesh);
758
759 doCorrectCoupled(volScalarField);
760 doCorrectCoupled(volVectorField);
761 doCorrectCoupled(volSphericalTensorField);
762 doCorrectCoupled(volSymmTensorField);
763 doCorrectCoupled(volTensorField);
764 #undef doCorrectCoupled
765 }
766 while (false);
767
768 // No update surface fields
769
770
771 // Map pointFields
772 if (nPointFields)
773 {
774 // Construct new pointMesh from distributed mesh
775 const pointMesh& newPointMesh = pointMesh::New(mesh);
776
777 pointDistributor.resetTarget(newPointMesh, distMap());
778
779 pointDistributor.distributeAndStore(pointScalarFields);
780 pointDistributor.distributeAndStore(pointVectorFields);
781 pointDistributor.distributeAndStore(pointSphTensorFields);
782 pointDistributor.distributeAndStore(pointSymmTensorFields);
783 pointDistributor.distributeAndStore(pointTensorFields);
784 }
785
786
787 // Set the minimum write precision
788 IOstream::defaultPrecision(max(10u, IOstream::defaultPrecision()));
789
790
791 if (!overwrite)
792 {
793 ++runTime;
794 mesh.setInstance(runTime.timeName());
795 }
796 else
797 {
798 mesh.setInstance(volMeshInstance);
799 }
800
801
802 // Register mapDistributePolyMesh for automatic writing...
804 (
806 (
807 "procAddressing",
808 mesh.facesInstance(),
809 polyMesh::meshSubDir,
810 mesh.thisDb(),
811 IOobject::NO_READ,
812 IOobject::AUTO_WRITE
813 ),
814 distMap()
815 );
816
817
818 if (reconstruct)
819 {
820 if (Pstream::master())
821 {
822 Info<< "Setting caseName to " << baseRunTime.caseName()
823 << " to write reconstructed mesh (and fields)." << endl;
824 runTime.caseName() = baseRunTime.caseName();
825 const bool oldProcCase(runTime.processorCase(false));
826
827 mesh.write();
828 topoSet::removeFiles(mesh);
829
830 // Now we've written all. Reset caseName on master
831 Info<< "Restoring caseName" << endl;
832 runTime.caseName() = proc0CaseName;
833 runTime.processorCase(oldProcCase);
834 }
835 }
836 else
837 {
838 autoPtr<fileOperation> defaultHandler;
839 if (writeHandler)
840 {
841 defaultHandler = fileHandler(std::move(writeHandler));
842 }
843
844 mesh.write();
845
846 if (defaultHandler)
847 {
848 writeHandler = fileHandler(std::move(defaultHandler));
849 }
850 topoSet::removeFiles(mesh);
851 }
852 Info<< "Written redistributed mesh to "
853 << mesh.facesInstance() << nl << endl;
854
855
856 if (decompose || reconstruct)
857 {
858 // Decompose (1 -> N) or reconstruct (N -> 1)
859 // so {boundary,cell,face,point}ProcAddressing have meaning
860 fvMeshTools::writeProcAddressing
861 (
862 mesh,
863 distMap(),
864 decompose,
865 std::move(writeHandler)
866 );
867 }
868 else
869 {
870 // Redistribute (N -> M)
871 // {boundary,cell,face,point}ProcAddressing would be incorrect
872 // - can either remove or redistribute previous
874 }
875
876
877 // Refinement data
878 {
879
880 // Read refinement data
881 if (Pstream::master() && decompose)
882 {
883 runTime.caseName() = baseRunTime.caseName();
884 runTime.processorCase(false);
885 }
887 (
888 "dummy",
889 mesh.facesInstance(),
890 polyMesh::meshSubDir,
891 mesh,
892 IOobject::READ_IF_PRESENT,
893 IOobject::NO_WRITE,
894 false
895 );
896
897 hexRef8Data refData(io);
898 if (Pstream::master() && decompose)
899 {
900 runTime.caseName() = proc0CaseName;
901 runTime.processorCase(oldProcCase);
902 }
903
904 // Make sure all processors have valid data (since only some will
905 // read)
906 refData.sync(io);
907
908 // Distribute
909 refData.distribute(distMap());
910
911
912 // Now we've read refinement data we can remove it
913 meshRefinement::removeFiles(mesh);
914
915 if (reconstruct)
916 {
917 if (Pstream::master())
918 {
919 const bool oldParRun = Pstream::parRun(false);
920
921 Info<< "Setting caseName to " << baseRunTime.caseName()
922 << " to write reconstructed refinement data." << endl;
923 runTime.caseName() = baseRunTime.caseName();
924 const bool oldProcCase(runTime.processorCase(false));
925
926 refData.write();
927
928 // Now we've written all. Reset caseName on master
929 Info<< "Restoring caseName" << endl;
930 runTime.caseName() = proc0CaseName;
931 runTime.processorCase(oldProcCase);
932
933 Pstream::parRun(oldParRun);
934 }
935 }
936 else
937 {
938 refData.write();
939 }
940 }
941
943 //{
944 // // Read sets
945 // if (Pstream::master() && decompose)
946 // {
947 // runTime.caseName() = baseRunTime.caseName();
948 // runTime.processorCase(false);
949 // }
950 // IOobjectList objects(mesh, mesh.facesInstance(), "polyMesh/sets");
951 //
952 // PtrList<cellSet> cellSets;
953 // ReadFields(objects, cellSets);
954 //
955 // if (Pstream::master() && decompose)
956 // {
957 // runTime.caseName() = proc0CaseName;
958 // runTime.processorCase(oldProcCase);
959 // }
960 //
961 // forAll(cellSets, i)
962 // {
963 // cellSets[i].distribute(distMap());
964 // }
965 //
966 // if (reconstruct)
967 // {
968 // if (Pstream::master())
969 // {
970 // Info<< "Setting caseName to " << baseRunTime.caseName()
971 // << " to write reconstructed refinement data." << endl;
972 // runTime.caseName() = baseRunTime.caseName();
973 // const bool oldProcCase(runTime.processorCase(false));
974 //
975 // forAll(cellSets, i)
976 // {
977 // cellSets[i].distribute(distMap());
978 // }
979 //
980 // // Now we've written all. Reset caseName on master
981 // Info<< "Restoring caseName" << endl;
982 // runTime.caseName() = proc0CaseName;
983 // runTime.processorCase(oldProcCase);
984 // }
985 // }
986 // else
987 // {
988 // forAll(cellSets, i)
989 // {
990 // cellSets[i].distribute(distMap());
991 // }
992 // }
993 //}
994
995
996 return distMap;
997}
998
999
1000/*---------------------------------------------------------------------------*\
1001 main
1002\*---------------------------------------------------------------------------*/
1003
1004int main(int argc, char *argv[])
1005{
1006 argList::addNote
1007 (
1008 "Redistribute decomposed mesh and fields according"
1009 " to the decomposeParDict settings.\n"
1010 "Optionally run in decompose/reconstruct mode"
1011 );
1012
1013 argList::noFunctionObjects(); // Never use function objects
1014
1015 // enable -constant ... if someone really wants it
1016 // enable -zeroTime to prevent accidentally trashing the initial fields
1017 timeSelector::addOptions(true, true);
1018
1019 #include "addAllRegionOptions.H"
1020
1021 #include "addOverwriteOption.H"
1022 argList::addBoolOption("decompose", "Decompose case");
1023 argList::addBoolOption("reconstruct", "Reconstruct case");
1024 argList::addVerboseOption("Additional verbosity");
1025 argList::addDryRunOption
1026 (
1027 "Test without writing the decomposition. "
1028 "Changes -cellDist to only write volScalarField."
1029 );
1030 argList::addVerboseOption("Additional verbosity");
1031 argList::addBoolOption
1032 (
1033 "cellDist",
1034 "Write cell distribution as a labelList - for use with 'manual' "
1035 "decomposition method or as a volScalarField for post-processing."
1036 );
1037 argList::addBoolOption
1038 (
1039 "newTimes",
1040 "Only reconstruct new times (i.e. that do not exist already)"
1041 );
1042 argList::addVerboseOption
1043 (
1044 "Additional verbosity. (Can be used multiple times)"
1045 );
1046 argList::addBoolOption
1047 (
1048 "no-finite-area",
1049 "Suppress finiteArea mesh/field handling",
1050 true // Advanced option
1051 );
1052
1053 // Handle arguments
1054 // ~~~~~~~~~~~~~~~~
1055 // (replacement for setRootCase that does not abort)
1056
1057 argList args(argc, argv);
1058
1059
1060 // As much as possible avoid synchronised operation. To be looked at more
1061 // closely for the three scenarios:
1062 // - decompose - reads on master (and from parent directory) and sends
1063 // dictionary to slaves
1064 // - distribute - reads on potentially a different number of processors
1065 // than it writes to
1066 // - reconstruct - reads parallel, write on master only and to parent
1067 // directory
1068 autoPtr<fileOperation> writeHandler;
1069 if
1070 (
1071 fileHandler().type()
1072 != fileOperations::uncollatedFileOperation::typeName
1073 )
1074 {
1075 // Install 'uncollated' as fileHandler. Save old one in writeHandler.
1076 writeHandler = fileHandler(fileOperation::NewUncollated());
1077 }
1078
1079 // Switch off parallel synchronisation of cached time directories
1080 fileHandler().distributed(true);
1081
1082 // File handler to be used for writing
1083 const fileOperation& fh
1084 (
1085 writeHandler
1086 ? writeHandler()
1087 : fileHandler()
1088 );
1089
1090 // Make sure to call findTimes on all processors to force caching of
1091 // time directories
1092 (void)fileHandler().findTimes(args.path(), "constant");
1093
1094 // Need this line since we don't include "setRootCase.H"
1095 #include "foamDlOpenLibs.H"
1096
1097 const bool reconstruct = args.found("reconstruct");
1098 const bool writeCellDist = args.found("cellDist");
1099 const bool dryrun = args.dryRun();
1100 const bool newTimes = args.found("newTimes");
1101 const int optVerbose = args.verbose();
1102
1103 const bool doFiniteArea = !args.found("no-finite-area");
1104 bool decompose = args.found("decompose");
1105 bool overwrite = args.found("overwrite");
1106
1107 if (optVerbose)
1108 {
1109 // Report on output
1110 faMeshDistributor::verbose_ = 1;
1111 parPointFieldDistributor::verbose_ = 1;
1112 }
1113
1114 // Disable NaN setting and floating point error trapping. This is to avoid
1115 // any issues inside the field redistribution inside fvMeshDistribute
1116 // which temporarily moves processor faces into existing patches. These
1117 // will now not have correct values. After all bits have been assembled
1118 // the processor fields will be restored but until then there might
1119 // be uninitialised values which might upset some patch field constructors.
1120 // Workaround by disabling floating point error trapping. TBD: have
1121 // actual field redistribution instead of subsetting inside
1122 // fvMeshDistribute.
1123 Foam::sigFpe::unset(true);
1124
1125 const wordRes selectedFields;
1126 const wordRes selectedLagrangianFields;
1127
1128
1129 if (decompose)
1130 {
1131 Info<< "Decomposing case (like decomposePar)" << nl << endl;
1132 if (reconstruct)
1133 {
1135 << "Cannot specify both -decompose and -reconstruct"
1136 << exit(FatalError);
1137 }
1138 }
1139 else if (reconstruct)
1140 {
1141 Info<< "Reconstructing case (like reconstructParMesh)" << nl << endl;
1142 }
1143
1144 if ((decompose || reconstruct) && !overwrite)
1145 {
1146 overwrite = true;
1148 << "Working in -decompose or -reconstruct mode:"
1149 " automatically implies -overwrite" << nl << endl;
1150 }
1151
1152 if (!Pstream::parRun())
1153 {
1155 << ": This utility can only be run parallel"
1156 << exit(FatalError);
1157 }
1158
1159
1160 if (!Foam::isDir(args.rootPath()))
1161 {
1163 << ": cannot open root directory " << args.rootPath()
1164 << exit(FatalError);
1165 }
1166
1167 // Detect if running data-distributed (multiple roots)
1168 bool nfs = true;
1169 {
1170 List<fileName> roots(1, args.rootPath());
1171 Pstream::combineAllGather(roots, ListOps::uniqueEqOp<fileName>());
1172 nfs = (roots.size() == 1);
1173 }
1174
1175 if (!nfs)
1176 {
1177 Info<< "Detected multiple roots i.e. non-nfs running"
1178 << nl << endl;
1179 }
1180
1181 // Check if we have processor directories. Ideally would like to
1182 // use fileHandler().dirPath here but we don't have runTime yet and
1183 // want to delay constructing runTime until we've synced all time
1184 // directories...
1185 const fileName procDir(fileHandler().filePath(args.path()));
1186 if (Foam::isDir(procDir))
1187 {
1188 if (decompose)
1189 {
1190 Info<< "Removing existing processor directory" << procDir << endl;
1191 fileHandler().rmDir(procDir);
1192 }
1193 }
1194 else
1195 {
1196 // Directory does not exist. If this happens on master -> decompose mode
1197 if (Pstream::master())
1198 {
1199 decompose = true;
1200 Info<< "No processor directories; switching on decompose mode"
1201 << nl << endl;
1202 }
1203 }
1204 // If master changed to decompose mode make sure all nodes know about it
1205 Pstream::scatter(decompose);
1206
1207
1208 // If running distributed we have problem of new processors not finding
1209 // a system/controlDict. However if we switch on the master-only reading
1210 // the problem becomes that the time directories are differing sizes and
1211 // e.g. latestTime will pick up a different time (which causes createTime.H
1212 // to abort). So for now make sure to have master times on all
1213 // processors
1214 if (!reconstruct)
1215 {
1216 Info<< "Creating time directories on all processors" << nl << endl;
1217 createTimeDirs(args.path());
1218 }
1219
1220 // Construct time
1221 // ~~~~~~~~~~~~~~
1222
1223 #include "createTime.H"
1224 runTime.functionObjects().off(); // Extra safety?
1225
1226
1227 // Save local processor0 casename
1228 const fileName proc0CaseName = runTime.caseName();
1229 const bool oldProcCase = runTime.processorCase();
1230
1231
1232 // Construct undecomposed Time
1233 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
1234 // This will read the same controlDict but might have a different
1235 // set of times so enforce same times
1236
1237 if (!nfs)
1238 {
1239 Info<< "Creating time directories for undecomposed Time"
1240 << " on all processors" << nl << endl;
1241 createTimeDirs(args.globalPath());
1242 }
1243
1244
1245 Info<< "Create undecomposed database" << nl << endl;
1246 Time baseRunTime
1247 (
1248 runTime.controlDict(),
1249 runTime.rootPath(),
1250 runTime.globalCaseName(),
1251 runTime.system(),
1252 runTime.constant(),
1253 false // enableFunctionObjects
1254 );
1255
1256
1257 wordHashSet masterTimeDirSet;
1258 if (newTimes)
1259 {
1260 instantList baseTimeDirs(baseRunTime.times());
1261 for (const instant& t : baseTimeDirs)
1262 {
1263 masterTimeDirSet.insert(t.name());
1264 }
1265 }
1266
1267
1268 // Allow override of decomposeParDict location
1269 const fileName decompDictFile =
1270 args.getOrDefault<fileName>("decomposeParDict", "");
1271
1272 // Get region names
1273 #include "getAllRegionOptions.H"
1274
1275 if (regionNames.size() == 1 && regionNames[0] != polyMesh::defaultRegion)
1276 {
1277 Info<< "Using region: " << regionNames[0] << nl << endl;
1278 }
1279
1280
1281 // Demand driven lagrangian mapper
1282 autoPtr<parLagrangianDistributor> lagrangianDistributorPtr;
1283
1284
1285 if (reconstruct)
1286 {
1287 // use the times list from the master processor
1288 // and select a subset based on the command-line options
1289 instantList timeDirs = timeSelector::select(runTime.times(), args);
1290 Pstream::scatter(timeDirs);
1291
1292 if (timeDirs.empty())
1293 {
1295 << "No times selected"
1296 << exit(FatalError);
1297 }
1298
1299
1300 // Pass1 : reconstruct mesh and addressing
1301 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1302
1303
1304 Info<< nl
1305 << "Reconstructing mesh and addressing" << nl << endl;
1306
1307 forAll(regionNames, regioni)
1308 {
1309 const word& regionName = regionNames[regioni];
1310 const word& regionDir = polyMesh::regionName(regionName);
1311
1312 const fileName volMeshSubDir(regionDir/polyMesh::meshSubDir);
1313 const fileName areaMeshSubDir(regionDir/faMesh::meshSubDir);
1314
1315 Info<< nl
1316 << "Reconstructing mesh " << regionDir << nl << endl;
1317
1318 bool areaMeshDetected = false;
1319
1320 // Loop over all times
1321 forAll(timeDirs, timeI)
1322 {
1323 // Set time for global database
1324 runTime.setTime(timeDirs[timeI], timeI);
1325 baseRunTime.setTime(timeDirs[timeI], timeI);
1326
1327 Info<< "Time = " << runTime.timeName() << endl << endl;
1328
1329 // Where meshes are
1330 fileName volMeshInstance;
1331 fileName areaMeshInstance;
1332
1333 volMeshInstance = runTime.findInstance
1334 (
1335 volMeshSubDir,
1336 "faces",
1337 IOobject::READ_IF_PRESENT
1338 );
1339
1340 if (doFiniteArea)
1341 {
1342 areaMeshInstance = runTime.findInstance
1343 (
1344 areaMeshSubDir,
1345 "faceLabels",
1346 IOobject::READ_IF_PRESENT
1347 );
1348 }
1349
1350 Pstream::broadcasts
1351 (
1352 UPstream::worldComm,
1353 volMeshInstance,
1354 areaMeshInstance
1355 );
1356
1357
1358 // Check processors have meshes
1359 // - check for 'faces' file (polyMesh)
1360 // - check for 'faceLabels' file (faMesh)
1361 boolList volMeshOnProc;
1362 boolList areaMeshOnProc;
1363
1364 volMeshOnProc = haveMeshFile
1365 (
1366 runTime,
1367 volMeshInstance/volMeshSubDir,
1368 "faces"
1369 );
1370
1371 if (doFiniteArea)
1372 {
1373 areaMeshOnProc = haveMeshFile
1374 (
1375 runTime,
1376 areaMeshInstance/areaMeshSubDir,
1377 "faceLabels"
1378 );
1379
1380 areaMeshDetected = areaMeshOnProc.found(true);
1381 }
1382
1383
1384 // Addressing back to reconstructed mesh as xxxProcAddressing.
1385 // - all processors have consistent faceProcAddressing
1386 // - processors without a mesh don't need faceProcAddressing
1387
1388
1389 // Note: filePath searches up on processors that don't have
1390 // processor if instance = constant so explicitly check
1391 // found filename.
1392 bool haveVolAddressing = false;
1393 if (volMeshOnProc[Pstream::myProcNo()])
1394 {
1395 // Read faces (just to know their size)
1396 faceCompactIOList faces
1397 (
1398 IOobject
1399 (
1400 "faces",
1401 volMeshInstance,
1402 volMeshSubDir,
1403 runTime,
1404 IOobject::MUST_READ
1405 )
1406 );
1407
1408 // Check faceProcAddressing
1409 labelIOList faceProcAddressing
1410 (
1411 IOobject
1412 (
1413 "faceProcAddressing",
1414 volMeshInstance,
1415 volMeshSubDir,
1416 runTime,
1417 IOobject::READ_IF_PRESENT
1418 )
1419 );
1420
1421 haveVolAddressing =
1422 (
1423 faceProcAddressing.headerOk()
1424 && faceProcAddressing.size() == faces.size()
1425 );
1426 }
1427 else
1428 {
1429 // Have no mesh. Don't need addressing
1430 haveVolAddressing = true;
1431 }
1432
1433 bool haveAreaAddressing = false;
1434 if (areaMeshOnProc[Pstream::myProcNo()])
1435 {
1436 // Read faces (just to know their size)
1437 labelIOList faceLabels
1438 (
1439 IOobject
1440 (
1441 "faceLabels",
1442 areaMeshInstance,
1443 areaMeshSubDir,
1444 runTime,
1445 IOobject::MUST_READ
1446 )
1447 );
1448
1449 // Check faceProcAddressing
1450 labelIOList faceProcAddressing
1451 (
1452 IOobject
1453 (
1454 "faceProcAddressing",
1455 areaMeshInstance,
1456 areaMeshSubDir,
1457 runTime,
1458 IOobject::READ_IF_PRESENT
1459 )
1460 );
1461
1462 haveAreaAddressing =
1463 (
1464 faceProcAddressing.headerOk()
1465 && faceProcAddressing.size() == faceLabels.size()
1466 );
1467 }
1468 else if (areaMeshDetected)
1469 {
1470 // Have no mesh. Don't need addressing
1471 haveAreaAddressing = true;
1472 }
1473
1474
1475 // Additionally check for master faces being readable. Could
1476 // do even more checks, e.g. global number of cells same
1477 // as cellProcAddressing
1478
1479 bool volMeshHaveUndecomposed = false;
1480 bool areaMeshHaveUndecomposed = false;
1481
1482 if (Pstream::master())
1483 {
1484 Info<< "Checking " << baseRunTime.caseName()
1485 << " for undecomposed volume and area meshes..."
1486 << endl;
1487
1488 const bool oldParRun = Pstream::parRun(false);
1489
1490 // Volume
1491 {
1492 faceCompactIOList facesIO
1493 (
1494 IOobject
1495 (
1496 "faces",
1497 volMeshInstance,
1498 volMeshSubDir,
1499 baseRunTime,
1500 IOobject::NO_READ
1501 ),
1502 label(0)
1503 );
1504 volMeshHaveUndecomposed = facesIO.headerOk();
1505 }
1506
1507 // Area
1508 if (doFiniteArea)
1509 {
1510 labelIOList labelsIO
1511 (
1512 IOobject
1513 (
1514 "faceLabels",
1515 areaMeshInstance,
1516 areaMeshSubDir,
1517 baseRunTime,
1518 IOobject::NO_READ
1519 )
1520 );
1521 areaMeshHaveUndecomposed = labelsIO.headerOk();
1522 }
1523
1524 Pstream::parRun(oldParRun); // Restore parallel state
1525 }
1526
1527 Pstream::broadcasts
1528 (
1529 UPstream::worldComm,
1530 volMeshHaveUndecomposed,
1531 areaMeshHaveUndecomposed
1532 );
1533
1534 // Report
1535 {
1536 Info<< " volume mesh ["
1537 << volMeshHaveUndecomposed << "] : "
1538 << volMeshInstance << nl
1539 << " area mesh ["
1540 << areaMeshHaveUndecomposed << "] : "
1541 << areaMeshInstance << nl
1542 << endl;
1543 }
1544
1545
1546 if
1547 (
1548 !volMeshHaveUndecomposed
1549 || !returnReduce(haveVolAddressing, andOp<bool>())
1550 )
1551 {
1552 Info<< "No undecomposed mesh. Creating from: "
1553 << volMeshInstance << endl;
1554
1555 if (areaMeshHaveUndecomposed)
1556 {
1557 areaMeshHaveUndecomposed = false;
1558 Info<< "Also ignore any undecomposed area mesh"
1559 << endl;
1560 }
1561
1562 autoPtr<fvMesh> volMeshPtr = fvMeshTools::loadOrCreateMesh
1563 (
1564 IOobject
1565 (
1566 regionName,
1567 volMeshInstance,
1568 runTime,
1570 ),
1571 decompose
1572 );
1573 fvMeshTools::setBasicGeometry(volMeshPtr());
1574 fvMesh& mesh = volMeshPtr();
1575
1576
1577 Info<< nl << "Reconstructing mesh" << nl << endl;
1578
1579 // Reconstruct (1 processor)
1580 const label nDestProcs(1);
1581 const labelList finalDecomp(mesh.nCells(), Zero);
1582
1583 redistributeAndWrite
1584 (
1585 std::move(writeHandler),
1586 baseRunTime,
1587 proc0CaseName,
1588
1589 // Controls
1590 false, // do not read fields
1591 false, // do not read undecomposed case on proc0
1592 true, // write redistributed files to proc0
1593 overwrite,
1594
1595 // Decomposition information
1596 nDestProcs,
1597 finalDecomp,
1598
1599 // For finite-volume
1600 volMeshOnProc,
1601 volMeshInstance,
1602 mesh
1603 );
1604 }
1605
1606
1607 // Similarly for finiteArea
1608 // - may or may not have undecomposed mesh
1609 // - may or may not have decomposed meshes
1610
1611 if
1612 (
1613 areaMeshOnProc.found(true) // ie, areaMeshDetected
1614 &&
1615 (
1616 !areaMeshHaveUndecomposed
1617 || !returnReduce(haveAreaAddressing, andOp<bool>())
1618 )
1619 )
1620 {
1621 Info<< "Loading area mesh from "
1622 << areaMeshInstance << endl;
1623
1624 Info<< " getting volume mesh support" << endl;
1625
1626 autoPtr<fvMesh> baseMeshPtr = fvMeshTools::newMesh
1627 (
1628 IOobject
1629 (
1630 regionName,
1631 baseRunTime.timeName(),
1632 baseRunTime,
1633 IOobject::MUST_READ
1634 ),
1635 true // read on master only
1636 );
1637 fvMeshTools::setBasicGeometry(baseMeshPtr());
1638
1639 autoPtr<fvMesh> volMeshPtr = fvMeshTools::loadOrCreateMesh
1640 (
1641 IOobject
1642 (
1643 regionName,
1644 baseMeshPtr().facesInstance(),
1645 runTime,
1647 ),
1648 decompose
1649 );
1650 fvMeshTools::setBasicGeometry(volMeshPtr());
1651 fvMesh& mesh = volMeshPtr();
1652
1653 // Read volume proc addressing back to base mesh
1655 (
1656 fvMeshTools::readProcAddressing(mesh, baseMeshPtr)
1657 );
1658
1659
1660 autoPtr<faMesh> areaMeshPtr = faMeshTools::loadOrCreateMesh
1661 (
1662 IOobject
1663 (
1664 regionName,
1665 areaMeshInstance,
1666 runTime,
1668 ),
1669 mesh, // <- The referenced polyMesh (from above)
1670 decompose
1671 );
1672 faMesh& areaMesh = areaMeshPtr();
1673
1674 faMeshTools::forceDemandDriven(areaMesh);
1675 faMeshTools::unregisterMesh(areaMesh);
1676
1677 autoPtr<faMesh> areaBaseMeshPtr;
1678
1679 // Reconstruct using polyMesh distribute map
1680 mapDistributePolyMesh faDistMap
1681 (
1682 faMeshDistributor::distribute
1683 (
1684 areaMesh,
1685 distMap(), // The polyMesh distMap
1686 baseMeshPtr(), // Target polyMesh
1687 areaBaseMeshPtr
1688 )
1689 );
1690
1691 faMeshTools::forceDemandDriven(areaBaseMeshPtr());
1692 faMeshTools::unregisterMesh(areaBaseMeshPtr());
1693
1694
1695 if (Pstream::master())
1696 {
1697 Info<< "Setting caseName to " << baseRunTime.caseName()
1698 << " to write reconstructed area mesh." << endl;
1699 runTime.caseName() = baseRunTime.caseName();
1700 const bool oldProcCase(runTime.processorCase(false));
1701
1702 areaBaseMeshPtr().write();
1703
1704 // Now we've written all. Reset caseName on master
1705 Info<< "Restoring caseName" << endl;
1706 runTime.caseName() = proc0CaseName;
1707 runTime.processorCase(oldProcCase);
1708 }
1709
1710 // Update for the reconstructed procAddressing
1711 faMeshTools::writeProcAddressing
1712 (
1713 areaBaseMeshPtr(), // Reconstruct location
1714 faDistMap,
1715 false, // decompose=false
1716 std::move(writeHandler),
1717 areaMeshPtr.get() // procMesh
1718 );
1719 }
1720 }
1721
1722 // Make sure all is finished writing until re-reading in pass2
1723 // below
1724 fileHandler().flush();
1725
1726
1727 // Pass2 : read mesh and addressing and reconstruct fields
1728 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1729
1730 Info<< nl
1731 << "Reconstructing fields" << nl << endl;
1732
1733 runTime.setTime(timeDirs[0], 0);
1734 baseRunTime.setTime(timeDirs[0], 0);
1735 Info<< "Time = " << runTime.timeName() << endl << endl;
1736
1737
1738 // Read undecomposed mesh on master and 'empty' mesh
1739 // (zero faces, point, cells but valid patches and zones) on slaves.
1740 // This is a bit of tricky code and hidden inside fvMeshTools for
1741 // now.
1742 Info<< "Reading undecomposed mesh (on master)" << endl;
1743
1744 autoPtr<fvMesh> baseMeshPtr = fvMeshTools::newMesh
1745 (
1746 IOobject
1747 (
1748 regionName,
1749 baseRunTime.timeName(),
1750 baseRunTime,
1751 IOobject::MUST_READ
1752 ),
1753 true // read on master only
1754 );
1755
1756 fvMeshTools::setBasicGeometry(baseMeshPtr());
1757
1758 Info<< "Reading local, decomposed mesh" << endl;
1759 autoPtr<fvMesh> volMeshPtr = fvMeshTools::loadOrCreateMesh
1760 (
1761 IOobject
1762 (
1763 regionName,
1764 baseMeshPtr().facesInstance(),
1765 runTime,
1767 ),
1768 decompose
1769 );
1770 fvMesh& mesh = volMeshPtr();
1771
1772
1773 // Similarly for finiteArea
1774 autoPtr<faMesh> areaBaseMeshPtr;
1775 autoPtr<faMesh> areaMeshPtr;
1776 autoPtr<faMeshDistributor> faDistributor;
1777 mapDistributePolyMesh areaDistMap;
1778
1779 if (areaMeshDetected)
1780 {
1781 areaBaseMeshPtr = faMeshTools::newMesh
1782 (
1783 IOobject
1784 (
1785 regionName,
1786 baseRunTime.timeName(),
1787 baseRunTime,
1788 IOobject::MUST_READ
1789 ),
1790 baseMeshPtr(),
1791 true // read on master only
1792 );
1793
1794 areaMeshPtr = faMeshTools::loadOrCreateMesh
1795 (
1796 IOobject
1797 (
1798 regionName,
1799 areaBaseMeshPtr().facesInstance(),
1800 runTime,
1801 IOobject::MUST_READ
1802 ),
1803 mesh,
1804 decompose
1805 );
1806
1807 areaDistMap =
1808 faMeshTools::readProcAddressing
1809 (
1810 areaMeshPtr(),
1811 areaBaseMeshPtr
1812 );
1813
1814 faMeshTools::forceDemandDriven(areaMeshPtr());
1815
1816 // Create an appropriate field distributor
1817 faDistributor.reset
1818 (
1820 (
1821 areaMeshPtr(), // source
1822 areaBaseMeshPtr(), // target
1823 areaDistMap,
1824 Pstream::master() // only write on master
1825 )
1826 );
1827 // Report some messages. Tbd.
1828 faMeshDistributor::verbose_ = 1;
1829 }
1830
1831
1832 if (writeHandler && Pstream::master())
1833 {
1834 // Remove any left-over empty processor directories created
1835 // by loadOrCreateMesh to get around e.g. collated start-up
1836 // problems. Should not happen in reconstruct mode ...
1837 const bool oldParRun = Pstream::parRun(false);
1838 removeEmptyDirs(mesh.time().path());
1839 Pstream::parRun(oldParRun);
1840 }
1841
1842
1843 // Read addressing back to base mesh
1845 distMap = fvMeshTools::readProcAddressing(mesh, baseMeshPtr);
1846
1847 // Construct field mapper
1848 auto fvDistributorPtr =
1850 (
1851 mesh, // source
1852 baseMeshPtr(), // target
1853 distMap(),
1854 UPstream::master() // Write reconstructed on master
1855 );
1856
1857 // Construct point field mapper
1858 const auto& basePointMesh = pointMesh::New(baseMeshPtr());
1859 const auto& procPointMesh = pointMesh::New(mesh);
1860
1861 auto pointFieldDistributorPtr =
1863 (
1864 procPointMesh, // source
1865 basePointMesh, // target
1866 distMap(),
1867 false, // delay
1868 UPstream::master() // Write reconstructed on master
1869 );
1870
1871
1872 // Since we start from Times[0] and not runTime.timeName() we
1873 // might overlook point motion in the first timestep
1874 // (since mesh.readUpdate() below will not be triggered). Instead
1875 // detect points by hand
1876 if (mesh.pointsInstance() != mesh.facesInstance())
1877 {
1878 Info<< " Detected initial mesh motion;"
1879 << " reconstructing points" << nl
1880 << endl;
1881 fvDistributorPtr().reconstructPoints();
1882 }
1883
1884
1885 // Loop over all times
1886 forAll(timeDirs, timeI)
1887 {
1888 if (newTimes && masterTimeDirSet.found(timeDirs[timeI].name()))
1889 {
1890 Info<< "Skipping time " << timeDirs[timeI].name()
1891 << nl << endl;
1892 continue;
1893 }
1894
1895 // Set time for global database
1896 runTime.setTime(timeDirs[timeI], timeI);
1897 baseRunTime.setTime(timeDirs[timeI], timeI);
1898
1899 Info<< "Time = " << runTime.timeName() << endl << endl;
1900
1901
1902 // Check if any new meshes need to be read.
1903 polyMesh::readUpdateState procStat = mesh.readUpdate();
1904
1905 if (procStat == polyMesh::POINTS_MOVED)
1906 {
1907 Info<< " Detected mesh motion; reconstructing points"
1908 << nl << endl;
1909 fvDistributorPtr().reconstructPoints();
1910 }
1911 else if
1912 (
1913 procStat == polyMesh::TOPO_CHANGE
1914 || procStat == polyMesh::TOPO_PATCH_CHANGE
1915 )
1916 {
1917 Info<< " Detected topology change;"
1918 << " reconstructing addressing" << nl << endl;
1919
1920 if (baseMeshPtr)
1921 {
1922 // Cannot do a baseMesh::readUpdate() since not all
1923 // processors will have mesh files. So instead just
1924 // recreate baseMesh
1925 baseMeshPtr.clear();
1926 baseMeshPtr = fvMeshTools::newMesh
1927 (
1928 IOobject
1929 (
1930 regionName,
1931 baseRunTime.timeName(),
1932 baseRunTime,
1933 IOobject::MUST_READ
1934 ),
1935 true // read on master only
1936 );
1937 }
1938
1939 // Re-read procAddressing
1940 distMap =
1941 fvMeshTools::readProcAddressing(mesh, baseMeshPtr);
1942
1943 // Reset field mappers
1944
1945 fvDistributorPtr.reset
1946 (
1948 (
1949 mesh, // source
1950 baseMeshPtr(), // target
1951 distMap(),
1952 UPstream::master() // Write reconstruct on master
1953 )
1954 );
1955
1956 // Construct point field mapper
1957 const auto& basePointMesh = pointMesh::New(baseMeshPtr());
1958 const auto& procPointMesh = pointMesh::New(mesh);
1959
1960 pointFieldDistributorPtr.reset
1961 (
1963 (
1964 procPointMesh, // source
1965 basePointMesh, // target
1966 distMap(),
1967 false, // delay until later
1968 UPstream::master() // Write reconstruct on master
1969 )
1970 );
1971
1972 lagrangianDistributorPtr.reset();
1973
1974 if (areaMeshPtr)
1975 {
1976 Info<< " Discarding finite-area addressing"
1977 << " (TODO)" << nl << endl;
1978
1979 areaBaseMeshPtr.reset();
1980 areaMeshPtr.reset();
1981 faDistributor.reset();
1982 areaDistMap.clear();
1983 }
1984 }
1985
1986
1987 // Get list of objects
1988 IOobjectList objects(mesh, runTime.timeName());
1989
1990
1991 // Mesh fields (vol, surface, volInternal)
1992 fvDistributorPtr()
1993 .distributeAllFields(objects, selectedFields);
1994
1995 // pointfields
1996 // - distribute and write (verbose)
1997 pointFieldDistributorPtr()
1998 .distributeAllFields(objects, selectedFields);
1999
2000
2001 // Clouds (note: might not be present on all processors)
2003 (
2004 lagrangianDistributorPtr,
2005 baseMeshPtr(),
2006 mesh,
2007 distMap(),
2008 selectedLagrangianFields
2009 );
2010
2011 if (faDistributor)
2012 {
2013 faDistributor()
2014 .distributeAllFields(objects, selectedFields);
2015 }
2016
2017 // If there are any "uniform" directories copy them from
2018 // the master processor
2019 copyUniform
2020 (
2021 fh,
2022 decompose,
2023 reconstruct,
2024 mesh.time().timeName(),
2025 mesh,
2026 baseMeshPtr()
2027 );
2028 // Non-region specific. Note: should do outside region loop
2029 // but would then have to replicate the whole time loop ...
2030 copyUniform
2031 (
2032 fh,
2033 decompose,
2034 reconstruct,
2035 mesh.time().timeName(),
2036 mesh.time(), // runTime
2037 baseMeshPtr().time() // baseRunTime
2038 );
2039 }
2040 }
2041 }
2042 else
2043 {
2044 // decompose or redistribution mode.
2045 // decompose : master : read from parent dir
2046 // slave : dummy mesh
2047 // redistribute : all read mesh or dummy mesh
2048
2049 // Time coming from processor0 (or undecomposed if no processor0)
2050 scalar masterTime;
2051 if (decompose)
2052 {
2053 // Use base time. This is to handle e.g. startTime = latestTime
2054 // which will not do anything if there are no processor directories
2055 masterTime = timeSelector::selectIfPresent
2056 (
2057 baseRunTime,
2058 args
2059 )[0].value();
2060 }
2061 else
2062 {
2063 masterTime = timeSelector::selectIfPresent
2064 (
2065 runTime,
2066 args
2067 )[0].value();
2068 }
2069 Pstream::scatter(masterTime);
2070 Info<< "Setting time to that of master or undecomposed case : "
2071 << masterTime << endl;
2072 runTime.setTime(masterTime, 0);
2073 baseRunTime.setTime(masterTime, 0);
2074
2075 // Save old time name (since might be incremented)
2076 const word oldTimeName(runTime.timeName());
2077
2078 forAll(regionNames, regioni)
2079 {
2080 const word& regionName = regionNames[regioni];
2081 const word& regionDir = polyMesh::regionName(regionName);
2082
2083 const fileName volMeshSubDir(regionDir/polyMesh::meshSubDir);
2084 const fileName areaMeshSubDir(regionDir/faMesh::meshSubDir);
2085
2086 Info<< nl << nl
2087 << (decompose ? "Decomposing" : "Redistributing")
2088 << " mesh " << regionDir << nl << endl;
2089
2090
2091 // Get time instance directory
2092 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
2093 // At this point we should be able to read at least a mesh on
2094 // processor0. Note the changing of the processor0 casename to
2095 // enforce it to read/write from the undecomposed case
2096
2097 fileName volMeshMasterInstance;
2098 fileName areaMeshMasterInstance;
2099
2100 // Assume to be true
2101 bool volMeshHaveUndecomposed = true;
2102 bool areaMeshHaveUndecomposed = doFiniteArea;
2103
2104 if (Pstream::master())
2105 {
2106 if (decompose)
2107 {
2108 Info<< "Checking undecomposed mesh in case: "
2109 << baseRunTime.caseName() << endl;
2110 runTime.caseName() = baseRunTime.caseName();
2111 runTime.processorCase(false);
2112 }
2113
2114 const bool oldParRun = Pstream::parRun(false);
2115 volMeshMasterInstance = runTime.findInstance
2116 (
2117 volMeshSubDir,
2118 "faces",
2119 IOobject::READ_IF_PRESENT
2120 );
2121
2122 if (doFiniteArea)
2123 {
2124 areaMeshMasterInstance = runTime.findInstance
2125 (
2126 areaMeshSubDir,
2127 "faceLabels",
2128 IOobject::READ_IF_PRESENT
2129 );
2130
2131 // Note: findInstance returns "constant" even if not found,
2132 // so recheck now for a false positive.
2133
2134 if ("constant" == areaMeshMasterInstance)
2135 {
2136 const boolList areaMeshOnProc
2137 (
2138 haveMeshFile
2139 (
2140 runTime,
2141 areaMeshMasterInstance/areaMeshSubDir,
2142 "faceLabels",
2143 false // verbose=false
2144 )
2145 );
2146
2147 if (areaMeshOnProc.empty() || !areaMeshOnProc[0])
2148 {
2149 areaMeshHaveUndecomposed = false;
2150 }
2151 }
2152 }
2153
2154 Pstream::parRun(oldParRun); // Restore parallel state
2155
2156 if (decompose)
2157 {
2158 Info<< " volume mesh ["
2159 << volMeshHaveUndecomposed << "] : "
2160 << volMeshMasterInstance << nl
2161 << " area mesh ["
2162 << areaMeshHaveUndecomposed << "] : "
2163 << areaMeshMasterInstance << nl
2164 << nl << nl;
2165
2166 // Restoring caseName
2167 runTime.caseName() = proc0CaseName;
2168 runTime.processorCase(oldProcCase);
2169 }
2170 }
2171
2172 Pstream::broadcasts
2173 (
2174 UPstream::worldComm,
2175 volMeshHaveUndecomposed,
2176 areaMeshHaveUndecomposed,
2177 volMeshMasterInstance,
2178 areaMeshMasterInstance
2179 );
2180
2181 // Check processors have meshes
2182 // - check for 'faces' file (polyMesh)
2183 // - check for 'faceLabels' file (faMesh)
2184 boolList volMeshOnProc;
2185 boolList areaMeshOnProc;
2186
2187 volMeshOnProc = haveMeshFile
2188 (
2189 runTime,
2190 volMeshMasterInstance/volMeshSubDir,
2191 "faces"
2192 );
2193
2194 if (doFiniteArea)
2195 {
2196 areaMeshOnProc = haveMeshFile
2197 (
2198 runTime,
2199 areaMeshMasterInstance/areaMeshSubDir,
2200 "faceLabels"
2201 );
2202 }
2203
2204 // Prior to loadOrCreateMesh, note which meshes already exist
2205 // for the current file handler.
2206 // - where mesh would be written if it didn't exist already.
2207 fileNameList volMeshDir(Pstream::nProcs());
2208 {
2209 volMeshDir[Pstream::myProcNo()] =
2210 (
2212 (
2213 IOobject
2214 (
2215 "faces",
2216 volMeshMasterInstance/volMeshSubDir,
2217 runTime
2218 ),
2220 ).path()
2221 );
2222
2223 Pstream::allGatherList(volMeshDir);
2224
2225 if (optVerbose && Pstream::master())
2226 {
2227 Info<< "Per processor faces dirs:" << nl
2228 << '(' << nl;
2229
2230 for (const int proci : Pstream::allProcs())
2231 {
2232 Info<< " "
2233 << runTime.relativePath(volMeshDir[proci]);
2234
2235 if (!volMeshOnProc[proci])
2236 {
2237 Info<< " [missing]";
2238 }
2239 Info<< nl;
2240 }
2241 Info<< ')' << nl << endl;
2242 }
2243 }
2244
2245 fileNameList areaMeshDir(Pstream::nProcs());
2246 if (doFiniteArea)
2247 {
2248 areaMeshDir[Pstream::myProcNo()] =
2249 (
2251 (
2252 IOobject
2253 (
2254 "faceLabels",
2255 areaMeshMasterInstance/areaMeshSubDir,
2256 runTime
2257 ),
2259 ).path()
2260 );
2261
2262 Pstream::allGatherList(areaMeshDir);
2263
2264 if (optVerbose && Pstream::master())
2265 {
2266 Info<< "Per processor faceLabels dirs:" << nl
2267 << '(' << nl;
2268
2269 for (const int proci : Pstream::allProcs())
2270 {
2271 Info<< " "
2272 << runTime.relativePath(areaMeshDir[proci]);
2273
2274 if (!areaMeshOnProc[proci])
2275 {
2276 Info<< " [missing]";
2277 }
2278 Info<< nl;
2279 }
2280 Info<< ')' << nl << endl;
2281 }
2282 }
2283
2284
2285 // Load mesh (or create dummy one)
2286 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2287
2288 if (Pstream::master() && decompose)
2289 {
2290 Info<< "Setting caseName to " << baseRunTime.caseName()
2291 << " to read undecomposed mesh" << endl;
2292 runTime.caseName() = baseRunTime.caseName();
2293 runTime.processorCase(false);
2294 }
2295
2296 // Volume mesh
2297 autoPtr<fvMesh> volMeshPtr = fvMeshTools::loadOrCreateMesh
2298 (
2299 IOobject
2300 (
2301 regionName,
2302 volMeshMasterInstance,
2303 runTime,
2305 ),
2306 decompose
2307 );
2308 fvMesh& mesh = volMeshPtr();
2309
2310
2311 // Area mesh
2312
2313 autoPtr<faMesh> areaMeshPtr;
2314
2315 // Decomposing: must have an undecomposed mesh
2316 // Redistributing: have any proc mesh
2317 if
2318 (
2319 doFiniteArea
2320 &&
2321 (
2322 decompose
2323 ? areaMeshHaveUndecomposed
2324 : areaMeshOnProc.found(true)
2325 )
2326 )
2327 {
2328 areaMeshPtr = faMeshTools::loadOrCreateMesh
2329 (
2330 IOobject
2331 (
2332 regionName,
2333 areaMeshMasterInstance,
2334 runTime,
2336 ),
2337 mesh, // <- The referenced polyMesh (from above)
2338 decompose
2339 );
2340
2341 faMeshTools::forceDemandDriven(*areaMeshPtr);
2342 faMeshTools::unregisterMesh(*areaMeshPtr);
2343 }
2344
2345
2346 if (writeHandler)
2347 {
2348 // Remove any left-over empty processor directories created
2349 // by loadOrCreateMesh to get around the collated start-up
2350 // problems
2351 if (Pstream::master()) //fileHandler().comm()))
2352 {
2353 const auto myProci = UPstream::myProcNo(); //comm()
2354 const auto& procs = UPstream::procID
2355 (
2356 UPstream::worldComm
2357 );
2358 const bool oldParRun = Pstream::parRun(false);
2359 for (const auto proci : procs)
2360 {
2361 if
2362 (
2363 !volMeshOnProc[proci]
2364 && volMeshDir[proci] != volMeshDir[myProci]
2365 )
2366 {
2367 Info<< "Deleting mesh dir:"
2368 << volMeshDir[proci] << endl;
2369 Foam::rmDir(volMeshDir[proci]);
2370 }
2371
2372 if
2373 (
2374 !areaMeshOnProc[proci]
2375 && areaMeshDir[proci] != areaMeshDir[myProci]
2376 )
2377 {
2378 Info<< "Deleting mesh dir:"
2379 << areaMeshDir[proci] << endl;
2380 Foam::rmDir(areaMeshDir[proci]);
2381 }
2382 }
2383
2384 // Remove empty directory
2385 removeEmptyDirs(mesh.time().path());
2386
2387 Pstream::parRun(oldParRun);
2388 }
2389 }
2390
2391
2392 if (Pstream::master() && decompose)
2393 {
2394 Info<< "Restoring caseName" << endl;
2395 runTime.caseName() = proc0CaseName;
2396 runTime.processorCase(oldProcCase);
2397 }
2398
2399 const label nOldCells = mesh.nCells();
2400
2401 // const label nOldAreaFaces =
2402 // (areaMeshPtr ? areaMeshPtr().nFaces() : 0);
2403 //
2404 //Pout<< "Loaded mesh : nCells:" << nOldCells
2405 // << " nPatches:" << mesh.boundaryMesh().size() << endl;
2406 //Pout<< "Loaded area mesh : nFaces:" << nOldAreaFaces << endl;
2407
2408 // Determine decomposition
2409 // ~~~~~~~~~~~~~~~~~~~~~~~
2410
2411 label nDestProcs;
2412 labelList finalDecomp;
2413 determineDecomposition
2414 (
2415 baseRunTime,
2416 decompDictFile,
2417 decompose,
2418 proc0CaseName,
2419 mesh,
2420 writeCellDist,
2421
2422 nDestProcs,
2423 finalDecomp
2424 );
2425
2426
2427 if (dryrun)
2428 {
2429 if (!Pstream::master())
2430 {
2431 if (areaMeshPtr && !areaMeshOnProc[Pstream::myProcNo()])
2432 {
2433 // Remove dummy mesh created by loadOrCreateMesh
2434 const bool oldParRun = Pstream::parRun(false);
2435 areaMeshPtr->removeFiles();
2436 Pstream::parRun(oldParRun); // Restore parallel state
2437 }
2438
2439 if (!volMeshOnProc[Pstream::myProcNo()])
2440 {
2441 // Remove dummy mesh created by loadOrCreateMesh
2442 const bool oldParRun = Pstream::parRun(false);
2443 mesh.removeFiles();
2444 Foam::rmDir(mesh.objectRegistry::objectPath());
2445 Pstream::parRun(oldParRun); // Restore parallel state
2446 }
2447 }
2448 continue;
2449 }
2450
2451 // Area fields first. Read and deregister
2453 if (areaMeshPtr)
2454 {
2455 areaFields.read
2456 (
2457 baseRunTime,
2458 proc0CaseName,
2459 decompose,
2460
2461 areaMeshOnProc,
2462 areaMeshMasterInstance,
2463 (*areaMeshPtr)
2464 );
2465 }
2466
2467
2468 // Detect lagrangian fields
2469 if (Pstream::master() && decompose)
2470 {
2471 runTime.caseName() = baseRunTime.caseName();
2472 runTime.processorCase(false);
2473 }
2474
2475 // Read lagrangian fields and store on cloud (objectRegistry)
2477 (
2478 readLagrangian
2479 (
2480 mesh,
2481 selectedLagrangianFields
2482 )
2483 );
2484
2485 if (Pstream::master() && decompose)
2486 {
2487 runTime.caseName() = proc0CaseName;
2488 runTime.processorCase(oldProcCase);
2489 }
2490
2491 // Load fields, do all distribution (mesh and fields)
2492 // - but not lagrangian fields; these are done later
2493 autoPtr<mapDistributePolyMesh> distMap = redistributeAndWrite
2494 (
2495 std::move(writeHandler),
2496 baseRunTime,
2497 proc0CaseName,
2498
2499 // Controls
2500 true, // read fields
2501 decompose, // decompose, i.e. read from undecomposed case
2502 false, // no reconstruction
2503 overwrite,
2504
2505 // Decomposition information
2506 nDestProcs,
2507 finalDecomp,
2508
2509 // For finite volume
2510 volMeshOnProc,
2511 volMeshMasterInstance,
2512 mesh
2513 );
2514
2515
2516 // Redistribute any clouds
2518 (
2519 lagrangianDistributorPtr,
2520 mesh,
2521 nOldCells,
2522 distMap(),
2523 clouds
2524 );
2525
2526
2527 // Redistribute area fields
2528
2529 mapDistributePolyMesh faDistMap;
2530 autoPtr<faMesh> areaProcMeshPtr;
2531
2532 if (areaMeshPtr)
2533 {
2534 faDistMap = faMeshDistributor::distribute
2535 (
2536 areaMeshPtr(),
2537 distMap(),
2538 areaProcMeshPtr
2539 );
2540
2541 // Force recreation of everything that might vaguely
2542 // be used by patches:
2543
2544 faMeshTools::forceDemandDriven(areaProcMeshPtr());
2545
2546
2547 if (reconstruct)
2548 {
2549 if (Pstream::master())
2550 {
2551 Info<< "Setting caseName to " << baseRunTime.caseName()
2552 << " to write reconstructed mesh (and fields)."
2553 << endl;
2554 runTime.caseName() = baseRunTime.caseName();
2555 const bool oldProcCase(runTime.processorCase(false));
2556 //const bool oldParRun = Pstream::parRun(false);
2557
2558 areaProcMeshPtr->write();
2559
2560 // Now we've written all. Reset caseName on master
2561 Info<< "Restoring caseName" << endl;
2562 runTime.caseName() = proc0CaseName;
2563 runTime.processorCase(oldProcCase);
2564 }
2565 }
2566 else
2567 {
2568 autoPtr<fileOperation> defaultHandler;
2569 if (writeHandler)
2570 {
2571 defaultHandler = fileHandler(std::move(writeHandler));
2572 }
2573
2575 (
2576 IOobject
2577 (
2578 "procAddressing",
2579 areaProcMeshPtr->facesInstance(),
2580 faMesh::meshSubDir,
2581 areaProcMeshPtr->thisDb(),
2582 IOobject::NO_READ,
2583 IOobject::NO_WRITE,
2584 false
2585 ),
2586 faDistMap
2587 ).write();
2588
2589 areaProcMeshPtr->write();
2590
2591 if (defaultHandler)
2592 {
2593 writeHandler = fileHandler(std::move(defaultHandler));
2594 }
2595
2596 if (decompose)
2597 {
2598 faMeshTools::writeProcAddressing
2599 (
2600 areaProcMeshPtr(),
2601 faDistMap,
2602 decompose,
2603 std::move(writeHandler)
2604 );
2605 }
2606 }
2607
2608 Info<< "Written redistributed mesh to "
2609 << areaProcMeshPtr->facesInstance() << nl << endl;
2610
2611 faMeshDistributor distributor
2612 (
2613 areaMeshPtr(), // source
2614 areaProcMeshPtr(), // target
2615 faDistMap
2616 );
2617
2618 areaFields.redistributeAndWrite(distributor, true);
2619 }
2620
2621 // Copy region-specific uniform
2622 // (e.g. solid/uniform/cumulativeContErr)
2623 copyUniform
2624 (
2625 fh,
2626 decompose,
2627 reconstruct,
2628 oldTimeName, // provided read time
2629 mesh, // read location is mesh (but oldTimeName)
2630 mesh // write location is mesh
2631 );
2632 }
2633
2634 // Copy non-region specific uniform (e.g. uniform/time)
2635 copyUniform
2636 (
2637 fh,
2638 decompose,
2639 reconstruct,
2640 oldTimeName, // provided read time
2641 ( // read location
2642 decompose
2643 ? baseRunTime
2644 : runTime
2645 ),
2646 runTime // writing location
2647 );
2648 }
2649
2650
2651 Info<< "End\n" << endl;
2652
2653 return 0;
2654}
2655
2656
2657// ************************************************************************* //
Inter-processor communication reduction functions.
Y[inertIndex] max(0.0)
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition: HashSet.H:191
bool found(const Key &key) const
Return true if hashed entry is found in table.
Definition: HashTableI.H:100
A IOmapDistributePolyMesh wrapper for using referenced external data.
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:59
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
instantList times() const
Search the case for valid time directories.
Definition: TimePaths.C:149
bool processorCase() const noexcept
Return true if this is a processor case.
Definition: TimePathsI.H:36
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
const fileName & caseName() const
Return case name.
Definition: TimePathsI.H:62
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
virtual void setTime(const Time &t)
Reset the time and time-index to those of the given time.
Definition: Time.C:997
bool found(const T &val, label pos=0) const
True if the value if found in the list.
Definition: UListI.H:265
bool empty() const noexcept
True if the UList is empty (ie, size() is zero)
Definition: UListI.H:427
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Finite area area (element) fields.
Mesh data needed to do the Finite Area discretisation.
Definition: areaFaMesh.H:56
Extract command arguments and options from the supplied argc and argv parameters.
Definition: argList.H:124
int dryRun() const noexcept
Return the dry-run flag.
Definition: argListI.H:116
int verbose() const noexcept
Return the verbose flag.
Definition: argListI.H:128
const fileName & rootPath() const noexcept
Return root path.
Definition: argListI.H:63
bool found(const word &optName) const
Return true if the named option is found.
Definition: argListI.H:178
fileName path() const
Return the full path to the (processor local) case.
Definition: argListI.H:81
T getOrDefault(const word &optName, const T &deflt) const
Get a value from the named option if present, or return default.
Definition: argListI.H:307
fileName globalPath() const
Return the full path to the global case.
Definition: argListI.H:87
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
void clear() noexcept
Same as reset(nullptr)
Definition: autoPtr.H:176
T * get() noexcept
Return pointer to managed object without nullptr checking.
Definition: autoPtr.H:152
void reset(autoPtr< T > &&other) noexcept
Delete managed object and set to new given pointer.
Definition: autoPtrI.H:117
Abstract base class for domain decomposition.
static label nDomains(const dictionary &decompDict, const word &regionName="")
Return region-specific or top-level numberOfSubdomains entry.
virtual labelList decompose(const pointField &points, const scalarField &pointWeights) const
Return the wanted processor number for every coordinate.
virtual bool parallelAware() const =0
Is method parallel aware?
MeshObject wrapper of decompositionMethod.
decompositionMethod & decomposer() const
Return demand-driven decomposition method.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Definition: dictionaryI.H:87
Holds a reference to the original mesh (the baseMesh) and optionally to a subset of that mesh (the su...
Finite area mesh (used for 2-D non-Euclidian finite area method) defined using a patch of faces on a ...
Definition: faMesh.H:100
A class for handling file names.
Definition: fileName.H:76
static std::string path(const std::string &str)
Return directory path name (part before last /)
Definition: fileNameI.H:176
An encapsulation of filesystem-related operations.
Definition: fileOperation.H:69
bool distributed() const noexcept
Distributed roots (parallel run)
virtual void flush() const
Forcibly wait until all output done. Flush any cached data.
virtual bool cp(const fileName &src, const fileName &dst, const bool followLink=true) const =0
Copy, recursively if necessary, the source to the destination.
virtual bool rmDir(const fileName &dir, const bool silent=false) const =0
Remove a directory and its contents.
virtual fileName dirPath(const bool checkGlobal, const IOobject &io, const bool search=true) const =0
Search for a directory. checkGlobal : also check undecomposed.
virtual instantList findTimes(const fileName &, const word &) const
Get sorted list of times.
virtual fileName objectPath(const IOobject &io, const word &typeName) const
Generate disk file name for object. Opposite of filePath.
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
Holds a reference to the original mesh (the baseMesh) and optionally to a subset of that mesh (the su...
Definition: fvMeshSubset.H:80
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:68
Various for reading/decomposing/reconstructing/distributing refinement data.
Definition: hexRef8Data.H:61
An instant of time. Contains the time value and name. Uses Foam::Time when formatting the name.
Definition: instant.H:56
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
void clear()
Reset to zero size, only retaining communicator(s)
Registry of regIOobjects.
const Time & time() const noexcept
Return time registry.
Simple container to manage read/write, redistribute finiteArea fields.
Finite volume reconstructor for volume and surface fields.
Distributor/redistributor for point fields, uses a two (or three) stage construction.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:55
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:91
Neighbour processor patch.
int neighbProcNo() const
Return neighbour processor number.
virtual bool write(const bool valid=true) const
Write using setting from DB.
static const triad unset
Definition: triad.H:97
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:54
A class for handling words, derived from Foam::string.
Definition: word.H:68
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Definition: zero.H:63
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
dynamicFvMesh & mesh
engineTime & runTime
Foam::word regionName(Foam::polyMesh::defaultRegion)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
const word & regionDir
label nPointFields
wordList regionNames
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
Miscellaneous file handling for meshes.
#define WarningInFunction
Report a warning using Foam::Warning.
int debugSwitch(const char *name, const int deflt=0)
Lookup debug switch or add default value.
Definition: debug.C:225
int debug
Static debugging option.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > reconstruct(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
type
Types of root.
Definition: Roots.H:55
Namespace for OpenFOAM.
const fileOperation & fileHandler()
Get current file handler.
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:515
void reconstructLagrangian(autoPtr< parLagrangianDistributor > &distributorPtr, const fvMesh &baseMesh, const fvMesh &mesh, const mapDistributePolyMesh &distMap, const wordRes &selectedFields)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
bool rmDir(const fileName &directory, const bool silent=false)
Remove a directory and its contents (optionally silencing warnings)
Definition: MSwindows.C:1036
void removeEmptyDirs(const fileName &path)
Remove empty directories from bottom up.
void removeProcAddressing(const faMesh &mesh)
Remove procAddressing.
boolList haveMeshFile(const Time &runTime, const fileName &meshPath, const word &meshFile="faces", const bool verbose=true)
Check for availability of specified mesh file (default: "faces")
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
bool isDir(const fileName &name, const bool followLink=true)
Does the name exist as a DIRECTORY in the file system?
Definition: MSwindows.C:651
void redistributeLagrangian(autoPtr< parLagrangianDistributor > &distributorPtr, const fvMesh &mesh, const label nOldCells, const mapDistributePolyMesh &distMap, PtrList< unmappedPassivePositionParticleCloud > &clouds)
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Reading, reconstruct, redistribution of lagrangian fields.
Foam::argList args(argc, argv)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
List helper to append y unique elements onto the end of x.
Definition: ListOps.H:570
static const char *const typeName
The type name used in ensight case files.
Foam::surfaceFields.