domainDecomposition.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2024 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "domainDecomposition.H"
27 #include "decompositionMethod.H"
28 #include "IOobjectList.H"
29 #include "cellSet.H"
30 #include "faceSet.H"
31 #include "fvMeshStitcher.H"
32 #include "pointSet.H"
33 #include "hexRef8Data.H"
34 #include "cyclicFvPatch.H"
35 #include "processorCyclicFvPatch.H"
36 #include "nonConformalFvPatch.H"
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 namespace Foam
41 {
43 }
44 
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47 
48 Foam::label Foam::domainDecomposition::compareInstances
49 (
50  const fileName& a,
51  const fileName& b
52 ) const
53 {
54  const word& constant = runTimes_.completeTime().constant();
55 
56  if (a == constant && b == constant) return 0;
57 
58  if (a == constant) return +1;
59 
60  if (b == constant) return -1;
61 
62  const scalar aValue = instant(a).value();
63  const scalar bValue = instant(b).value();
64 
65  if (aValue < bValue) return +1;
66 
67  if (aValue > bValue) return -1;
68 
69  return 0;
70 }
71 
72 
73 void Foam::domainDecomposition::validateComplete() const
74 {
75  if (!completeMesh_.valid())
76  {
78  << "Complete data requested but complete mesh has not been "
79  << "generated or read" << exit(FatalError);
80  }
81 }
82 
83 
84 void Foam::domainDecomposition::validateProcs() const
85 {
86  if (!procMeshes_.set(0))
87  {
89  << "Decomposed data requested but decomposed mesh has not been "
90  << "generated or read" << exit(FatalError);
91  }
92 }
93 
94 
95 void Foam::domainDecomposition::readComplete()
96 {
97  completeMesh_.reset
98  (
99  new fvMesh
100  (
101  IOobject
102  (
103  regionName_,
104  runTimes_.completeTime().name(),
105  runTimes_.completeTime(),
108  ),
109  false
110  )
111  );
112 
113  completeMesh_->postConstruct(false, fvMesh::stitchType::none);
114 }
115 
116 
117 void Foam::domainDecomposition::readProcs()
118 {
119  for (label proci = 0; proci < nProcs(); proci++)
120  {
121  procMeshes_.set
122  (
123  proci,
124  new fvMesh
125  (
126  IOobject
127  (
128  regionName_,
129  runTimes_.procTimes()[proci].name(),
130  runTimes_.procTimes()[proci],
133  ),
134  false
135  )
136  );
137 
138  procMeshes_[proci].postConstruct(false, fvMesh::stitchType::none);
139  }
140 }
141 
142 
143 void Foam::domainDecomposition::readCompleteAddressing()
144 {
145  cellProc_ =
147  (
148  IOobject
149  (
150  "cellProc",
151  completeMesh().facesInstance(),
152  completeMesh().meshSubDir,
153  completeMesh(),
156  false
157  )
158  );
159 }
160 
161 
162 void Foam::domainDecomposition::readProcsAddressing()
163 {
164  for (label proci = 0; proci < nProcs(); proci++)
165  {
166  const fvMesh& procMesh = procMeshes_[proci];
167 
168  procPointAddressing_[proci] =
170  (
171  IOobject
172  (
173  "pointProcAddressing",
174  procMesh.facesInstance(),
175  procMesh.meshSubDir,
176  procMesh,
179  false
180  )
181  );
182 
183  procFaceAddressing_[proci] =
185  (
186  IOobject
187  (
188  "faceProcAddressing",
189  procMesh.facesInstance(),
190  procMesh.meshSubDir,
191  procMesh,
194  false
195  )
196  );
197 
198  procCellAddressing_[proci] =
200  (
201  IOobject
202  (
203  "cellProcAddressing",
204  procMesh.facesInstance(),
205  procMesh.meshSubDir,
206  procMesh,
209  false
210  )
211  );
212  }
213 }
214 
215 
216 void Foam::domainDecomposition::readAddressing()
217 {
218  readCompleteAddressing();
219  readProcsAddressing();
220  procFaceAddressingBf_.clear();
221 }
222 
223 
224 Foam::fvMesh::readUpdateState Foam::domainDecomposition::readUpdate()
225 {
226  validateComplete();
227  validateProcs();
228 
229  // Do read-update on all meshes
231  completeMesh_->readUpdate(fvMesh::stitchType::none);
232 
233  forAll(runTimes_.procTimes(), proci)
234  {
235  fvMesh::readUpdateState procStat =
236  procMeshes_[proci].readUpdate(fvMesh::stitchType::none);
237 
238  stat = procStat > stat ? procStat : stat;
239  }
240 
241  return stat;
242 }
243 
244 
245 void Foam::domainDecomposition::writeCompleteAddressing() const
246 {
247  labelIOList cellProc
248  (
249  IOobject
250  (
251  "cellProc",
252  completeMesh().facesInstance(),
253  completeMesh().meshSubDir,
254  completeMesh(),
257  ),
258  cellProc_
259  );
260 
261  cellProc.write();
262 }
263 
264 
265 void Foam::domainDecomposition::writeProcsAddressing() const
266 {
267  for (label proci = 0; proci < nProcs(); proci++)
268  {
269  const fvMesh& procMesh = procMeshes_[proci];
270 
271  labelIOList pointProcAddressing
272  (
273  IOobject
274  (
275  "pointProcAddressing",
276  procMesh.facesInstance(),
277  procMesh.meshSubDir,
278  procMesh,
281  ),
282  procPointAddressing_[proci]
283  );
284  pointProcAddressing.write();
285 
286  labelIOList faceProcAddressing
287  (
288  IOobject
289  (
290  "faceProcAddressing",
291  procMesh.facesInstance(),
292  procMesh.meshSubDir,
293  procMesh,
296  ),
297  procFaceAddressing_[proci]
298  );
299  faceProcAddressing.write();
300 
301  labelIOList cellProcAddressing
302  (
303  IOobject
304  (
305  "cellProcAddressing",
306  procMesh.facesInstance(),
307  procMesh.meshSubDir,
308  procMesh,
311  ),
312  procCellAddressing_[proci]
313  );
314  cellProcAddressing.write();
315  }
316 }
317 
318 
319 void Foam::domainDecomposition::writeAddressing() const
320 {
321  writeCompleteAddressing();
322  writeProcsAddressing();
323 }
324 
325 
326 void Foam::domainDecomposition::writeProcPoints(const fileName& inst)
327 {
328  IOobject completePointsIo
329  (
330  "points",
331  inst,
333  completeMesh(),
336  false
337  );
338 
339  if (!completePointsIo.headerOk()) return;
340 
341  const pointIOField completePoints(completePointsIo);
342 
343  for (label proci = 0; proci < nProcs(); proci++)
344  {
345  pointIOField procPoints
346  (
347  IOobject
348  (
349  "points",
350  inst,
352  procMeshes()[proci],
355  false
356  ),
357  pointField
358  (
359  completePoints,
360  procPointAddressing_[proci]
361  )
362  );
363 
364  procPoints.write();
365  }
366 }
367 
368 
369 void Foam::domainDecomposition::writeCompletePoints(const fileName& inst)
370 {
371  pointIOField completePoints
372  (
373  IOobject
374  (
375  "points",
376  inst,
378  completeMesh(),
381  false
382  ),
383  pointField(completeMesh().nPoints())
384  );
385 
386  for (label proci = 0; proci < nProcs(); proci++)
387  {
388  IOobject procPointsIo
389  (
390  "points",
391  inst,
393  procMeshes()[proci],
396  false
397  );
398 
399  if (!procPointsIo.headerOk()) return;
400 
401  completePoints.rmap
402  (
403  pointIOField(procPointsIo),
404  procPointAddressing_[proci]
405  );
406  }
407 
408  completePoints.write();
409 }
410 
411 
412 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
413 
415 (
416  const processorRunTimes& runTimes,
417  const word& regionName,
418  const multiDomainDecomposition& regionMeshes
419 )
420 :
421  runTimes_(runTimes),
422  regionName_(regionName),
423  completeMesh_(nullptr),
424  procMeshes_(nProcs()),
425  regionMeshes_(regionMeshes),
426  cellProc_(),
427  procPointAddressing_(nProcs()),
428  procFaceAddressing_(nProcs()),
429  procCellAddressing_(nProcs()),
430  procFaceAddressingBf_()
431 {}
432 
433 
434 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
435 
437 {}
438 
439 
440 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
441 
443 {
444  readComplete();
445 
447  (
448  "cellProc",
449  completeMesh().facesInstance(),
451  completeMesh()
452  );
453  IOobject procFaceIo
454  (
455  "faces",
456  completeMesh().facesInstance(),
457  completeMesh().meshDir(),
458  runTimes_.procTimes()[0]
459  );
460  typeIOobject<labelIOList> procAddrIo
461  (
462  "cellProcAddressing",
463  completeMesh().facesInstance(),
464  completeMesh().meshDir(),
465  runTimes_.procTimes()[0]
466  );
467 
468  const bool load = addrIo.headerOk() && procFaceIo.headerOk();
469 
470  if (load)
471  {
472  readProcs();
473 
474  if (procAddrIo.headerOk())
475  {
476  readAddressing();
477  }
478  else
479  {
480  readCompleteAddressing();
481 
483  << nl << " Processor meshes exist but have no addressing."
484  << nl << nl << " This could be because the processor meshes "
485  << "have changed. Decomposing the" << nl << " mesh would "
486  << "overwrite that change. If you are sure that this is "
487  << "appropriate," << nl << " then delete the "
488  << fileName("processor*")/procFaceIo.relativePath().c_str()
489  << " directories and re-run this" << nl << " command."
490  << exit(FatalError);
491  }
492 
493  decomposePoints();
494  }
495  else
496  {
497  if
498  (
499  completeMesh().facesInstance()
500  != runTimes_.completeTime().name()
501  && completeMesh().facesInstance()
502  != runTimes_.completeTime().constant()
503  )
504  {
506  << "Cannot begin mesh decomposition at time "
507  << fileName(runTimes_.completeTime().name()) << nl
508  << "The mesh at this instant is that of an earlier"
509  << " time " << completeMesh().facesInstance() << nl
510  << "Decomposition must start from this earlier time"
511  << exit(FatalError);
512  }
513 
514  decompose();
515  }
516 
517  return !load;
518 }
519 
520 
522 {
523  completeMesh_->stitcher().connect(false, false, true);
524 }
525 
526 
528 {
529  if (!completeConformal())
530  {
531  procFaceAddressingBf_.clear();
532 
533  forAll(procMeshes_, proci)
534  {
535  procMeshes_[proci].conform();
536  }
537 
538  unconform();
539  }
540 }
541 
542 
544 (
545  const bool decomposed,
546  const bool doSets
547 )
548 {
549  writeProcs(doSets);
550 
551  if (decomposed)
552  {
553  writeProcPoints(completeMesh().facesInstance());
554  }
555 }
556 
557 
559 {
560  readProcs();
561 
562  IOobject faceIo
563  (
564  "faces",
565  procMeshes()[0].facesInstance(),
566  procMeshes()[0].meshDir(),
567  runTimes_.completeTime()
568  );
570  (
571  "cellProc",
572  procMeshes()[0].facesInstance(),
573  procMeshes()[0].meshDir(),
574  runTimes_.completeTime()
575  );
576  typeIOobject<labelIOList> procAddrIo
577  (
578  "cellProcAddressing",
579  procMeshes()[0].facesInstance(),
581  procMeshes()[0]
582  );
583 
584  const bool load = faceIo.headerOk() && procAddrIo.headerOk();
585 
586  if (load)
587  {
588  typeIOobject<pointIOField> completePointsIo
589  (
590  "points",
591  procMeshes()[0].pointsInstance(),
592  procMeshes()[0].meshDir(),
593  runTimes_.completeTime()
594  );
595 
596  readComplete();
597 
598  if (addrIo.headerOk())
599  {
600  readAddressing();
601  }
602  else
603  {
604  readProcsAddressing();
605 
607  << nl << " A complete mesh exists but has no "
608  << addrIo.name() << " addressing." << nl << nl << " This "
609  << "could be because the complete mesh has changed. "
610  << "Reconstructing the" << nl << " mesh would overwrite "
611  << "that change. If you are sure that this is appropriate,"
612  << nl << " then delete the " << faceIo.relativePath()
613  << " directory and re-run this command." << nl << nl
614  << " Or, it could be because the complete and processor "
615  << "meshes were decomposed" << nl << " by a version of "
616  << "OpenFOAM that pre-dates the automatic generation of "
617  << nl << " " << addrIo.name() << " addressing. This will be "
618  << "assumed and the " << addrIo.name() << " addressing will"
619  << nl << " be re-built" << nl << endl;
620 
621  cellProc_ = labelList(completeMesh().nCells(), -1);
622 
623  for (label proci = 0; proci < nProcs(); proci++)
624  {
626  (
627  cellProc_,
628  procCellAddressing_[proci]
629  ) = proci;
630  }
631 
632  writeCompleteAddressing();
633  }
634 
635  reconstructPoints();
636  }
637  else
638  {
639  if
640  (
641  procMeshes()[0].facesInstance()
642  != runTimes_.procTimes()[0].name()
643  && procMeshes()[0].facesInstance()
644  != runTimes_.procTimes()[0].constant()
645  )
646  {
648  << "Cannot begin mesh reconstruction at time "
649  << fileName(runTimes_.procTimes()[0].name()) << nl
650  << "The mesh at this instant is that of an earlier"
651  << " time " << procMeshes()[0].facesInstance() << nl
652  << "Reconstruction must start from this earlier time"
653  << exit(FatalError);
654  }
655 
656  reconstruct();
657  }
658 
659  return !load;
660 }
661 
662 
664 {
665  forAll(procMeshes_, proci)
666  {
667  procMeshes_[proci].stitcher().connect(false, false, true);
668  }
669 }
670 
671 
673 {
674  if (!procsConformal())
675  {
676  procFaceAddressingBf_.clear();
677 
678  completeMesh_->conform();
679 
680  unconform();
681  }
682 }
683 
684 
686 (
687  const bool reconstructed,
688  const bool doSets
689 )
690 {
691  writeComplete(doSets);
692 
693  if (reconstructed)
694  {
695  writeCompletePoints(procMeshes()[0].facesInstance());
696  }
697 }
698 
699 
701 {
702  const fvMesh::readUpdateState stat = readUpdate();
703 
704  // Topology changes
705  {
706  const label facesCompare =
707  compareInstances
708  (
709  completeMesh().facesInstance(),
710  procMeshes_[0].facesInstance()
711  );
712 
713  // If the complete mesh has newer topology then we need to decompose
714  if (facesCompare == -1)
715  {
716  decompose();
717  }
718 
719  // If there has been matching topology change then reload the addressing
720  if (facesCompare == 0 && stat >= fvMesh::TOPO_CHANGE)
721  {
722  readAddressing();
723  }
724 
725  // The processor meshes should not have newer topology when decomposing
726  if (facesCompare == +1)
727  {
729  << "Cannot decompose at time "
730  << procMeshes_[0].facesInstance()
731  << " because the processor mesh topology has evolved further"
732  << " than the complete mesh topology." << exit(FatalError);
733  }
734  }
735 
736  // Geometry changes
737  {
738  const label pointsCompare =
739  compareInstances
740  (
741  completeMesh().pointsInstance(),
742  procMeshes_[0].pointsInstance()
743  );
744 
745  // If the complete mesh has newer geometry then we need to decompose
746  // the points
747  if (pointsCompare == -1)
748  {
749  decomposePoints();
750  }
751 
752  // The processor meshes should not have newer geometry when decomposing
753  if (pointsCompare == +1)
754  {
756  << "Cannot decompose at time "
757  << procMeshes_[0].pointsInstance()
758  << " because the processor mesh geometry has evolved further"
759  << " than the complete mesh geometry." << exit(FatalError);
760  }
761  }
762 
763  return stat;
764 }
765 
766 
768 (
769  const fvMesh::readUpdateState stat
770 )
771 {
772  if (completeMesh_->stitcher().stitches() && stat != fvMesh::UNCHANGED)
773  {
774  procFaceAddressingBf_.clear();
775 
776  completeMesh_->stitcher().connect(false, false, true);
777  }
778 }
779 
780 
782 (
783  const fvMesh::readUpdateState stat
784 )
785 {
786  if (completeMesh_->stitcher().stitches() && stat != fvMesh::UNCHANGED)
787  {
788  forAll(procMeshes_, proci)
789  {
790  procMeshes_[proci].conform();
791  }
792 
793  unconform();
794  }
795 }
796 
797 
799 {
800  const fvMesh::readUpdateState stat = readUpdate();
801 
802  // Topology changes
803  {
804  const label facesCompare =
805  compareInstances
806  (
807  completeMesh().facesInstance(),
808  procMeshes_[0].facesInstance()
809  );
810 
811  // The complete mesh should not have newer topology when reconstructing
812  if (facesCompare == -1)
813  {
815  << "Cannot reconstruct at time "
816  << completeMesh().facesInstance()
817  << " because the complete mesh topology has evolved further"
818  << " than the processor mesh topology." << exit(FatalError);
819  }
820 
821  // If there has been matching topology change then reload the addressing
822  if (facesCompare == 0 && stat >= fvMesh::TOPO_CHANGE)
823  {
824  readAddressing();
825  }
826 
827  // If the processor meshes have newer topology then we need to
828  // reconstruct
829  if (facesCompare == +1)
830  {
831  reconstruct();
832  }
833  }
834 
835  // Geometry changes
836  {
837  const label pointsCompare =
838  compareInstances
839  (
840  completeMesh().pointsInstance(),
841  procMeshes_[0].pointsInstance()
842  );
843 
844  // The complete mesh should not have newer geometry when reconstructing
845  if (pointsCompare == -1)
846  {
848  << "Cannot reconstruct at time "
849  << completeMesh().pointsInstance()
850  << " because the complete mesh geometry has evolved further"
851  << " than the processor mesh geometry." << exit(FatalError);
852  }
853 
854  // If the processor meshes have newer geometry then we need to
855  // reconstruct the points
856  if (pointsCompare == +1)
857  {
858  reconstructPoints();
859  }
860  }
861 
862  return stat;
863 }
864 
865 
867 (
868  const fvMesh::readUpdateState stat
869 )
870 {
871  if (completeMesh_->stitcher().stitches() && stat != fvMesh::UNCHANGED)
872  {
873  procFaceAddressingBf_.clear();
874 
875  forAll(procMeshes_, proci)
876  {
877  procMeshes_[proci].stitcher().connect(false, false, true);
878  }
879  }
880 }
881 
882 
884 (
885  const fvMesh::readUpdateState stat
886 )
887 {
888  if (completeMesh_->stitcher().stitches() && stat != fvMesh::UNCHANGED)
889  {
890  completeMesh_->conform();
891 
892  unconform();
893  }
894 }
895 
896 
899 {
900  validateComplete();
901  validateProcs();
902 
903  if (procFaceAddressingBf_.empty())
904  {
905  // Get any non-conformal proc-face addressing
906  List<List<DynamicList<label>>> nonConformalProcFaceAddressingBf =
907  this->nonConformalProcFaceAddressingBf();
908 
909  // Build finite volume face addressing boundary fields
910  procFaceAddressingBf_.resize(nProcs());
911  forAll(procMeshes_, proci)
912  {
913  const fvMesh& procMesh = procMeshes_[proci];
914 
915  procFaceAddressingBf_.set
916  (
917  proci,
919  (
920  procMesh.boundary(),
922  calculatedFvsPatchLabelField::typeName
923  )
924  );
925 
926  forAll(procMesh.boundary(), procPatchi)
927  {
928  const fvPatch& fvp = procMesh.boundary()[procPatchi];
929 
930  if (isA<nonConformalFvPatch>(fvp))
931  {
932  procFaceAddressingBf_[proci][procPatchi] =
933  nonConformalProcFaceAddressingBf[proci][procPatchi];
934  }
935  else if (isA<processorCyclicFvPatch>(fvp))
936  {
937  const label completePatchi =
938  refCast<const processorCyclicFvPatch>(fvp)
939  .referPatchIndex();
940 
941  procFaceAddressingBf_[proci][procPatchi] =
942  mag(fvp.patchSlice(procFaceAddressing_[proci]))
943  - completeMesh().boundaryMesh()[completePatchi].start();
944  }
945  else if (isA<processorFvPatch>(fvp))
946  {
947  procFaceAddressingBf_[proci][procPatchi] =
948  fvp.patchSlice(procFaceAddressing_[proci]);
949  }
950  else
951  {
952  procFaceAddressingBf_[proci][procPatchi] =
953  mag(fvp.patchSlice(procFaceAddressing_[proci]))
954  - completeMesh().boundaryMesh()[procPatchi].start();
955  }
956  }
957  }
958  }
959 
960  return procFaceAddressingBf_;
961 }
962 
963 
964 void Foam::domainDecomposition::writeComplete(const bool doSets) const
965 {
966  const bool topologyWrite =
967  static_cast<const faceCompactIOList&>(completeMesh().faces())
968  .writeOpt() == IOobject::AUTO_WRITE;
969 
970  // Set the precision of the points data to be min 10
972 
973  // Write the complete mesh
974  completeMesh().write();
975 
976  // Everything written below is topological data, so quit here if not
977  // writing topology
978  if (!topologyWrite) return;
979 
980  // Read, reconstruct and write sets
981  if (doSets)
982  {
985  HashPtrTable<pointSet> pointSets;
986 
987  for (label proci = 0; proci < nProcs(); proci++)
988  {
989  const fvMesh& procMesh = procMeshes_[proci];
990 
991  const labelList& cellMap = procCellAddressing()[proci];
992  const labelList& faceMap = procFaceAddressing()[proci];
993  const labelList& pointMap = procPointAddressing()[proci];
994 
995  // Scan the contents of the sets directory
996  IOobjectList setObjects
997  (
998  procMesh,
999  runTimes_.procTimes()[proci].name(),
1000  polyMesh::meshSubDir/"sets"
1001  );
1002  IOobjectList cellSetObjects
1003  (
1004  setObjects.lookupClass(cellSet::typeName)
1005  );
1006  IOobjectList faceSetObjects
1007  (
1008  setObjects.lookupClass(faceSet::typeName)
1009  );
1010  IOobjectList pointSetObjects
1011  (
1012  setObjects.lookupClass(pointSet::typeName)
1013  );
1014 
1015  // Read and reconstruct the sets
1016  forAllConstIter(IOobjectList, cellSetObjects, iter)
1017  {
1018  const cellSet procSet(*iter());
1019 
1020  if (!cellSets.found(iter.key()))
1021  {
1022  cellSets.insert
1023  (
1024  iter.key(),
1025  new cellSet
1026  (
1027  completeMesh(),
1028  iter.key(),
1029  procSet.size()
1030  )
1031  );
1032  }
1033 
1034  cellSet& cSet = *cellSets[iter.key()];
1035 
1036  cSet.instance() = runTimes_.completeTime().name();
1037 
1038  forAllConstIter(cellSet, procSet, iter)
1039  {
1040  cSet.insert(cellMap[iter.key()]);
1041  }
1042  }
1043  forAllConstIter(IOobjectList, faceSetObjects, iter)
1044  {
1045  const faceSet procSet(*iter());
1046 
1047  if (!faceSets.found(iter.key()))
1048  {
1049  faceSets.insert
1050  (
1051  iter.key(),
1052  new faceSet
1053  (
1054  completeMesh(),
1055  iter.key(),
1056  procSet.size()
1057  )
1058  );
1059  }
1060 
1061  faceSet& cSet = *faceSets[iter.key()];
1062 
1063  cSet.instance() = runTimes_.completeTime().name();
1064 
1065  forAllConstIter(faceSet, procSet, iter)
1066  {
1067  cSet.insert(mag(faceMap[iter.key()]) - 1);
1068  }
1069  }
1070  forAllConstIter(IOobjectList, pointSetObjects, iter)
1071  {
1072  const pointSet procSet(*iter());
1073 
1074  if (!pointSets.found(iter.key()))
1075  {
1076  pointSets.insert
1077  (
1078  iter.key(),
1079  new pointSet
1080  (
1081  completeMesh(),
1082  iter.key(),
1083  procSet.size()
1084  )
1085  );
1086  }
1087 
1088  pointSet& cSet = *pointSets[iter.key()];
1089 
1090  cSet.instance() = runTimes_.completeTime().name();
1091 
1092  forAllConstIter(pointSet, procSet, iter)
1093  {
1094  cSet.insert(pointMap[iter.key()]);
1095  }
1096  }
1097  }
1098 
1099  // Write the sets
1101  {
1102  iter()->write();
1103  }
1105  {
1106  iter()->write();
1107  }
1108  forAllConstIter(HashPtrTable<pointSet>, pointSets, iter)
1109  {
1110  iter()->write();
1111  }
1112  }
1113 
1114  // Read, decompose, and write refinement data (if any)
1115  UPtrList<const labelList> cellMaps(nProcs());
1116  UPtrList<const labelList> pointMaps(nProcs());
1117  PtrList<const hexRef8Data> refinementDatas(nProcs());
1118  for (label proci = 0; proci < nProcs(); proci++)
1119  {
1120  const fvMesh& procMesh = procMeshes_[proci];
1121 
1122  cellMaps.set(proci, &procCellAddressing_[proci]);
1123  pointMaps.set(proci, &procPointAddressing_[proci]);
1124  refinementDatas.set
1125  (
1126  proci,
1127  new hexRef8Data
1128  (
1129  IOobject
1130  (
1131  "dummy",
1132  completeMesh().facesInstance(),
1134  procMesh,
1137  false
1138  )
1139  )
1140  );
1141  }
1142  hexRef8Data
1143  (
1144  IOobject
1145  (
1146  "dummy",
1147  completeMesh().facesInstance(),
1149  completeMesh(),
1152  false
1153  ),
1154  cellMaps,
1155  pointMaps,
1156  refinementDatas
1157  ).write();
1158 
1159  // Write decomposition addressing
1160  writeAddressing();
1161 }
1162 
1163 
1164 void Foam::domainDecomposition::writeProcs(const bool doSets) const
1165 {
1166  const bool topologyWrite =
1167  static_cast<const faceCompactIOList&>(procMeshes()[0].faces())
1168  .writeOpt() == IOobject::AUTO_WRITE;
1169 
1170  // Write out the meshes
1171  for (label proci = 0; proci < nProcs(); proci++)
1172  {
1173  const fvMesh& procMesh = procMeshes_[proci];
1174 
1175  // Set the precision of the points data to be min 10
1177 
1178  // Write the processor mesh
1179  procMesh.write();
1180  }
1181 
1182  // Everything written below is topological data, so quit here if not
1183  // writing topology
1184  if (!topologyWrite) return;
1185 
1186  // Read, decompose, and write any sets
1187  if (doSets)
1188  {
1189  // Scan the contents of the sets directory
1190  IOobjectList setObjects
1191  (
1192  completeMesh(),
1193  completeMesh().facesInstance(),
1194  polyMesh::meshSubDir/"sets"
1195  );
1196  IOobjectList cellSetObjects
1197  (
1198  setObjects.lookupClass(cellSet::typeName)
1199  );
1200  IOobjectList faceSetObjects
1201  (
1202  setObjects.lookupClass(faceSet::typeName)
1203  );
1204  IOobjectList pointSetObjects
1205  (
1206  setObjects.lookupClass(pointSet::typeName)
1207  );
1208 
1209  // Read the sets
1211  forAllConstIter(IOobjectList, cellSetObjects, iter)
1212  {
1213  cellSets.append(new cellSet(*iter()));
1214  }
1216  forAllConstIter(IOobjectList, faceSetObjects, iter)
1217  {
1218  faceSets.append(new faceSet(*iter()));
1219  }
1220  PtrList<const pointSet> pointSets;
1221  forAllConstIter(IOobjectList, pointSetObjects, iter)
1222  {
1223  pointSets.append(new pointSet(*iter()));
1224  }
1225 
1226  // Decompose and write sets into the processor mesh directories
1227  for (label proci = 0; proci < nProcs(); proci++)
1228  {
1229  const fvMesh& procMesh = procMeshes_[proci];
1230 
1231  forAll(cellSets, i)
1232  {
1233  const cellSet& cs = cellSets[i];
1234  cellSet set(procMesh, cs.name(), cs.size()/nProcs());
1235  forAll(procCellAddressing_[proci], i)
1236  {
1237  if (cs.found(procCellAddressing_[proci][i]))
1238  {
1239  set.insert(i);
1240  }
1241  }
1242  set.write();
1243  }
1244  forAll(faceSets, i)
1245  {
1246  const faceSet& cs = faceSets[i];
1247  faceSet set(procMesh, cs.name(), cs.size()/nProcs());
1248  forAll(procFaceAddressing_[proci], i)
1249  {
1250  if (cs.found(mag(procFaceAddressing_[proci][i]) - 1))
1251  {
1252  set.insert(i);
1253  }
1254  }
1255  set.write();
1256  }
1257  forAll(pointSets, i)
1258  {
1259  const pointSet& cs = pointSets[i];
1260  pointSet set(procMesh, cs.name(), cs.size()/nProcs());
1261  forAll(procPointAddressing_[proci], i)
1262  {
1263  if (cs.found(procPointAddressing_[proci][i]))
1264  {
1265  set.insert(i);
1266  }
1267  }
1268  set.write();
1269  }
1270  }
1271  }
1272 
1273  // Read, decompose, and write refinement data (if any)
1275  (
1276  IOobject
1277  (
1278  "dummy",
1279  completeMesh().facesInstance(),
1281  completeMesh(),
1284  false
1285  )
1286  );
1287  for (label proci = 0; proci < nProcs(); proci++)
1288  {
1289  const fvMesh& procMesh = procMeshes_[proci];
1290 
1291  hexRef8Data
1292  (
1293  IOobject
1294  (
1295  "dummy",
1296  completeMesh().facesInstance(),
1298  procMesh,
1301  false
1302  ),
1304  procCellAddressing_[proci],
1305  procPointAddressing_[proci]
1306  ).write();
1307  }
1308 
1309  // Write decomposition addressing
1310  writeAddressing();
1311 }
1312 
1313 
1314 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Generic GeometricBoundaryField class.
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
A HashTable specialisation for hashing pointers.
Definition: HashPtrTable.H:67
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:109
bool set(const Key &, const T &newElmt)
Set a new hashedEntry, overwriting existing entries.
Definition: HashTableI.H:91
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:138
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:53
IOobjectList lookupClass(const word &className) const
Return the list for all IOobjects of a given class.
Definition: IOobjectList.C:190
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
fileName relativePath() const
Return the path relative to the case directory.
Definition: IOobject.C:393
fileName & instance() const
Return the instance directory, constant, system, <time> etc.
Definition: IOobject.C:355
bool headerOk()
Read header of local object without type-checking.
const word & name() const
Return name.
Definition: IOobject.H:310
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:458
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
bool set(const label) const
Is element set.
Definition: PtrListI.H:62
void append(T *)
Append an element at the end of the list.
Definition: PtrListI.H:39
static const word & constant()
Return constant name.
Definition: TimePaths.H:122
A List with indirect addressing.
Definition: UIndirectList.H:60
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: UPtrList.H:66
bool set(const label) const
Is element set.
Definition: UPtrListI.H:87
A collection of cell labels.
Definition: cellSet.H:51
Automatic domain decomposition class for finite-volume meshes.
bool readReconstruct()
Read in the processor meshes. Read the complete mesh if it is.
void postReadReconstruct()
Post-read-construction steps for the meshes after read-reconstruct.
void postReadDecompose()
Post-read-construction steps for the meshes after read-decompose.
bool readDecompose()
Read in the complete mesh. Read the processor meshes if they are.
void unconformReadReconstruct()
Synchronise non-conformal structure after read-reconstruct.
void unconformReadUpdateReconstruct(const fvMesh::readUpdateState stat)
Synchronise non-conformal structure after read-update-reconstruct.
domainDecomposition(const processorRunTimes &procRunTimes, const word &regionName, const multiDomainDecomposition &regionMeshes=NullObjectRef< multiDomainDecomposition >())
Construct from processor run times and region name.
void postReadUpdateReconstruct(const fvMesh::readUpdateState stat)
Complete read-update-reconstruct.
void writeProcs(const bool doSets) const
Write the decomposed meshes and associated data.
void writeReadReconstruct(const bool reconstructed, const bool doSets)
Write following read-reconstruct.
void writeReadDecompose(const bool decomposed, const bool doSets)
Write following read-decompose.
void writeComplete(const bool doSets) const
Write the decomposed meshes and associated data.
void postReadUpdateDecompose(const fvMesh::readUpdateState stat)
Complete read-update-decompose.
const PtrList< surfaceLabelField::Boundary > & procFaceAddressingBf() const
Access the labels of faces for each processor (see notes above)
void unconformReadDecompose()
Synchronise non-conformal structure after read-decompose.
fvMesh::readUpdateState readUpdateReconstruct()
Read-update for reconstruction.
void unconformReadUpdateDecompose(const fvMesh::readUpdateState stat)
Synchronise non-conformal structure after read-update-decompose.
fvMesh::readUpdateState readUpdateDecompose()
Read-update for decomposition.
virtual ~domainDecomposition()
Destructor.
A list of face labels.
Definition: faceSet.H:51
A class for handling file names.
Definition: fileName.H:82
word name() const
Return file name (part beyond last /)
Definition: fileName.C:195
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:857
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: fvMesh.H:108
@ TOPO_CHANGE
Definition: fvMesh.H:112
virtual bool write(const bool write=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:1737
const word & name() const
Return reference to name.
Definition: fvMesh.H:432
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
const List< T >::subList patchSlice(const List< T > &l) const
Slice list to patch.
Definition: fvPatch.H:172
Various for reading/decomposing/reconstructing/distributing refinement data.
Definition: hexRef8Data.H:58
bool write() const
Write.
Definition: hexRef8Data.C:372
A set of point labels.
Definition: pointSet.H:51
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:270
const Time & completeTime() const
Access the complete run time.
Transfers refinement levels such that slow transition between levels is maintained....
virtual bool write(const bool write=true) const
Write using setting from DB.
Templated form of IOobject providing type information for file reading and header type checking.
Definition: IOobject.H:531
bool headerOk()
Read header (uses typeGlobalFile to find file) and check.
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label nPoints
volScalarField & b
Definition: createFields.H:25
#define WarningInFunction
Report a warning using Foam::Warning.
tmp< VolField< typename outerProduct< vector, Type >::type > > reconstruct(const SurfaceField< Type > &ssf)
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
List< label > labelList
A List of labels.
Definition: labelList.H:56
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
const word & regionName(const solver &region)
Definition: solver.H:209
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
dimensioned< scalar > mag(const dimensioned< Type > &)
vectorIOField pointIOField
pointIOField is a vectorIOField.
Definition: pointIOField.H:42
defineTypeNameAndDebug(combustionModel, 0)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
IOList< label > labelIOList
Label container classes.
Definition: labelIOList.H:42
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
error FatalError
static const char nl
Definition: Ostream.H:266