decomposePar.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) 2016-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 Application
28  decomposePar
29 
30 Group
31  grpParallelUtilities
32 
33 Description
34  Automatically decomposes a mesh and fields of a case for parallel
35  execution of OpenFOAM.
36 
37 Usage
38  \b decomposePar [OPTIONS]
39 
40  Options:
41  - \par -allRegions
42  Decompose all regions in regionProperties. Does not check for
43  existence of processor*.
44 
45  - \par -case <dir>
46  Specify case directory to use (instead of the cwd).
47 
48  - \par -cellDist
49  Write the cell distribution as a labelList, for use with 'manual'
50  decomposition method and as a volScalarField for visualization.
51 
52  - \par -constant
53  Include the 'constant/' dir in the times list.
54 
55  - \par -copyUniform
56  Copy any \a uniform directories too.
57 
58  - \par -copyZero
59  Copy \a 0 directory to processor* rather than decompose the fields.
60 
61  - \par -debug-switch <name=val>
62  Specify the value of a registered debug switch. Default is 1
63  if the value is omitted. (Can be used multiple times)
64 
65  - \par -decomposeParDict <file>
66  Use specified file for decomposePar dictionary.
67 
68  - \par -dry-run
69  Test without writing the decomposition. Changes -cellDist to
70  only write volScalarField.
71 
72  - \par -fields
73  Use existing geometry decomposition and convert fields only.
74 
75  - \par fileHandler <handler>
76  Override the file handler type.
77 
78  - \par -force
79  Remove any existing \a processor subdirectories before decomposing the
80  geometry.
81 
82  - \par -ifRequired
83  Only decompose the geometry if the number of domains has changed from a
84  previous decomposition. No \a processor subdirectories will be removed
85  unless the \a -force option is also specified. This option can be used
86  to avoid redundant geometry decomposition (eg, in scripts), but should
87  be used with caution when the underlying (serial) geometry or the
88  decomposition method etc. have been changed between decompositions.
89 
90  - \par -info-switch <name=val>
91  Specify the value of a registered info switch. Default is 1
92  if the value is omitted. (Can be used multiple times)
93 
94  - \par -latestTime
95  Select the latest time.
96 
97  - \par -lib <name>
98  Additional library or library list to load (can be used multiple times).
99 
100  - \par -noFunctionObjects
101  Do not execute function objects.
102 
103  - \par -noSets
104  Skip decomposing cellSets, faceSets, pointSets.
105 
106  - \par -noZero
107  Exclude the \a 0 dir from the times list.
108 
109  - \par -opt-switch <name=val>
110  Specify the value of a registered optimisation switch (int/bool).
111  Default is 1 if the value is omitted. (Can be used multiple times)
112 
113  - \par -region <regionName>
114  Decompose named region. Does not check for existence of processor*.
115 
116  - \par -time <ranges>
117  Override controlDict settings and decompose selected times. Does not
118  re-decompose the mesh i.e. does not handle moving mesh or changing
119  mesh cases. Eg, ':10,20 40:70 1000:', 'none', etc.
120 
121  - \par -verbose
122  Additional verbosity.
123 
124  - \par -doc
125  Display documentation in browser.
126 
127  - \par -doc-source
128  Display source code in browser.
129 
130  - \par -help
131  Display short help and exit.
132 
133  - \par -help-man
134  Display full help (manpage format) and exit.
135 
136  - \par -help-notes
137  Display help notes (description) and exit.
138 
139  - \par -help-full
140  Display full help and exit.
141 
142 \*---------------------------------------------------------------------------*/
143 
144 #include "OSspecific.H"
145 #include "fvCFD.H"
146 #include "IOobjectList.H"
147 #include "domainDecomposition.H"
149 #include "labelIOField.H"
150 #include "labelFieldIOField.H"
151 #include "scalarIOField.H"
152 #include "scalarFieldIOField.H"
153 #include "vectorIOField.H"
154 #include "vectorFieldIOField.H"
155 #include "sphericalTensorIOField.H"
157 #include "symmTensorIOField.H"
158 #include "symmTensorFieldIOField.H"
159 #include "tensorIOField.H"
160 #include "tensorFieldIOField.H"
161 #include "pointFields.H"
162 #include "regionProperties.H"
163 
164 #include "readFields.H"
165 #include "dimFieldDecomposer.H"
166 #include "fvFieldDecomposer.H"
167 #include "pointFieldDecomposer.H"
169 #include "decompositionModel.H"
170 
171 #include "faCFD.H"
172 #include "emptyFaPatch.H"
173 #include "faMeshDecomposition.H"
174 #include "faFieldDecomposer.H"
175 
176 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
177 
178 namespace Foam
179 {
180 
181 const labelIOList& procAddressing
182 (
183  const PtrList<fvMesh>& procMeshList,
184  const label proci,
185  const word& name,
186  PtrList<labelIOList>& procAddressingList
187 )
188 {
189  const fvMesh& procMesh = procMeshList[proci];
190 
191  if (!procAddressingList.set(proci))
192  {
193  procAddressingList.set
194  (
195  proci,
196  new labelIOList
197  (
198  IOobject
199  (
200  name,
201  procMesh.facesInstance(),
202  procMesh.meshSubDir,
203  procMesh,
206  false
207  )
208  )
209  );
210  }
211  return procAddressingList[proci];
212 }
213 
214 
215 void decomposeUniform
216 (
217  const bool copyUniform,
218  const domainDecomposition& mesh,
219  const Time& processorDb,
220  const word& regionDir = word::null
221 )
222 {
223  const Time& runTime = mesh.time();
224 
225  // Any uniform data to copy/link?
226  const fileName uniformDir(regionDir/"uniform");
227 
228  if (fileHandler().isDir(runTime.timePath()/uniformDir))
229  {
230  Info<< "Detected additional non-decomposed files in "
231  << runTime.timePath()/uniformDir
232  << endl;
233 
234  const fileName timePath =
235  fileHandler().filePath(processorDb.timePath());
236 
237  // If no fields have been decomposed the destination
238  // directory will not have been created so make sure.
239  mkDir(timePath);
240 
241  if (copyUniform || mesh.distributed())
242  {
243  if (!fileHandler().exists(timePath/uniformDir))
244  {
245  fileHandler().cp
246  (
247  runTime.timePath()/uniformDir,
248  timePath/uniformDir
249  );
250  }
251  }
252  else
253  {
254  // Link with relative paths
255  string parentPath = string("..")/"..";
256 
257  if (regionDir != word::null)
258  {
259  parentPath = parentPath/"..";
260  }
261 
262  fileName currentDir(cwd());
263  chDir(timePath);
264 
265  if (!fileHandler().exists(uniformDir))
266  {
267  fileHandler().ln
268  (
269  parentPath/runTime.timeName()/uniformDir,
270  uniformDir
271  );
272  }
273  chDir(currentDir);
274  }
275  }
276 }
277 
278 }
279 
280 
281 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
282 
283 int main(int argc, char *argv[])
284 {
285  argList::addNote
286  (
287  "Decompose a mesh and fields of a case for parallel execution"
288  );
289 
290  argList::noParallel();
291  argList::addOption
292  (
293  "decomposeParDict",
294  "file",
295  "Use specified file for decomposePar dictionary"
296  );
297  #include "addRegionOption.H"
298  argList::addBoolOption
299  (
300  "allRegions",
301  "Operate on all regions in regionProperties"
302  );
303  argList::addBoolOption
304  (
305  "dry-run",
306  "Test without writing the decomposition. "
307  "Changes -cellDist to only write volScalarField."
308  );
309  argList::addBoolOption
310  (
311  "verbose",
312  "Additional verbosity"
313  );
314  argList::addBoolOption
315  (
316  "cellDist",
317  "Write cell distribution as a labelList - for use with 'manual' "
318  "decomposition method and as a volScalarField for visualization."
319  );
320  argList::addBoolOption
321  (
322  "copyZero",
323  "Copy 0/ directory to processor*/ rather than decompose the fields"
324  );
325  argList::addBoolOption
326  (
327  "copyUniform",
328  "Copy any uniform/ directories too"
329  );
330  argList::addBoolOption
331  (
332  "fields",
333  "Use existing geometry decomposition and convert fields only"
334  );
335  argList::addBoolOption
336  (
337  "noSets",
338  "Skip decomposing cellSets, faceSets, pointSets"
339  );
340  argList::addBoolOption
341  (
342  "force",
343  "Remove existing processor*/ subdirs before decomposing the geometry"
344  );
345  argList::addBoolOption
346  (
347  "ifRequired",
348  "Only decompose geometry if the number of domains has changed"
349  );
350 
351  // Allow explicit -constant, have zero from time range
352  timeSelector::addOptions(true, false); // constant(true), zero(false)
353 
354  #include "setRootCase.H"
355 
356  const bool dryrun = args.found("dry-run");
357  const bool optRegion = args.found("region");
358  const bool allRegions = args.found("allRegions");
359  const bool writeCellDist = args.found("cellDist");
360  const bool verbose = args.found("verbose");
361 
362  // Most of these are ignored for dry-run (not triggered anywhere)
363  const bool copyZero = args.found("copyZero");
364  const bool copyUniform = args.found("copyUniform");
365  const bool decomposeSets = !args.found("noSets");
366  const bool decomposeIfRequired = args.found("ifRequired");
367 
368  bool decomposeFieldsOnly = args.found("fields");
369  bool forceOverwrite = args.found("force");
370 
371 
372  // Set time from database
373  #include "createTime.H"
374 
375  // Allow override of time (unless dry-run)
376  instantList times;
377  if (dryrun)
378  {
379  Info<< "\ndry-run: ignoring -copy*, -fields, -force, time selection"
380  << nl;
381  }
382  else
383  {
384  times = timeSelector::selectIfPresent(runTime, args);
385  }
386 
387  // Allow override of decomposeParDict location
388  fileName decompDictFile(args.get<fileName>("decomposeParDict", ""));
389  if (!decompDictFile.empty() && !decompDictFile.isAbsolute())
390  {
391  decompDictFile = runTime.globalPath()/decompDictFile;
392  }
393 
394  // Get all region names
395  wordList regionNames;
396  if (allRegions)
397  {
398  regionNames = regionProperties(runTime).names();
399 
400  Info<< "Decomposing all regions in regionProperties" << nl
401  << " " << flatOutput(regionNames) << nl << endl;
402  }
403  else
404  {
405  regionNames.resize(1);
406  regionNames.first() =
407  args.getOrDefault<word>("region", fvMesh::defaultRegion);
408  }
409 
410  forAll(regionNames, regioni)
411  {
412  const word& regionName = regionNames[regioni];
413  const word& regionDir =
414  (regionName == fvMesh::defaultRegion ? word::null : regionName);
415 
416  if (dryrun)
417  {
418  Info<< "dry-run: decomposing mesh " << regionName << nl << nl
419  << "Create mesh..." << flush;
420 
421  domainDecompositionDryRun decompTest
422  (
423  IOobject
424  (
425  regionName,
426  runTime.timeName(),
427  runTime,
428  IOobject::MUST_READ,
429  IOobject::NO_WRITE,
430  false
431  ),
432  decompDictFile
433  );
434 
435  decompTest.execute(writeCellDist, verbose);
436  continue;
437  }
438 
439  Info<< "\n\nDecomposing mesh " << regionName << nl << endl;
440 
441  // Determine the existing processor count directly
442  label nProcs = fileHandler().nProcs(runTime.path(), regionDir);
443 
444  // Get requested numberOfSubdomains directly from the dictionary.
445  // Note: have no mesh yet so cannot use decompositionModel::New
446  const label nDomains = decompositionMethod::nDomains
447  (
448  IOdictionary
449  (
450  IOobject::selectIO
451  (
452  IOobject
453  (
454  decompositionModel::canonicalName,
455  runTime.time().system(),
456  regionDir, // region (if non-default)
457  runTime,
458  IOobject::MUST_READ,
459  IOobject::NO_WRITE,
460  false
461  ),
462  decompDictFile
463  )
464  )
465  );
466 
467  // Give file handler a chance to determine the output directory
468  const_cast<fileOperation&>(fileHandler()).setNProcs(nDomains);
469 
470  if (decomposeFieldsOnly)
471  {
472  // Sanity check on previously decomposed case
473  if (nProcs != nDomains)
474  {
476  << "Specified -fields, but the case was decomposed with "
477  << nProcs << " domains"
478  << nl
479  << "instead of " << nDomains
480  << " domains as specified in decomposeParDict" << nl
481  << exit(FatalError);
482  }
483  }
484  else if (nProcs)
485  {
486  bool procDirsProblem = true;
487 
488  if (decomposeIfRequired && nProcs == nDomains)
489  {
490  // We can reuse the decomposition
491  decomposeFieldsOnly = true;
492  procDirsProblem = false;
493  forceOverwrite = false;
494 
495  Info<< "Using existing processor directories" << nl;
496  }
497 
498  if (allRegions || optRegion)
499  {
500  procDirsProblem = false;
501  forceOverwrite = false;
502  }
503 
504  if (forceOverwrite)
505  {
506  Info<< "Removing " << nProcs
507  << " existing processor directories" << endl;
508 
509  // Remove existing processors directory
510  fileNameList dirs
511  (
513  (
514  runTime.path(),
515  fileName::Type::DIRECTORY
516  )
517  );
518  forAllReverse(dirs, diri)
519  {
520  const fileName& d = dirs[diri];
521 
522  label proci = -1;
523 
524  if
525  (
526  d.starts_with("processor")
527  &&
528  (
529  // Collated is "processors"
530  d[9] == 's'
531 
532  // Uncollated has integer(s) after 'processor'
533  || Foam::read(d.substr(9), proci)
534  )
535  )
536  {
537  if (fileHandler().exists(d))
538  {
539  fileHandler().rmDir(d);
540  }
541  }
542  }
543 
544  procDirsProblem = false;
545  }
546 
547  if (procDirsProblem)
548  {
550  << "Case is already decomposed with " << nProcs
551  << " domains, use the -force option or manually" << nl
552  << "remove processor directories before decomposing. e.g.,"
553  << nl
554  << " rm -rf " << runTime.path().c_str() << "/processor*"
555  << nl
556  << exit(FatalError);
557  }
558  }
559 
560  Info<< "Create mesh" << endl;
561  domainDecomposition mesh
562  (
563  IOobject
564  (
565  regionName,
566  runTime.timeName(),
567  runTime,
568  IOobject::NO_READ,
569  IOobject::NO_WRITE,
570  false
571  ),
572  decompDictFile
573  );
574 
575  // Decompose the mesh
576  if (!decomposeFieldsOnly)
577  {
578  mesh.decomposeMesh();
579 
580  mesh.writeDecomposition(decomposeSets);
581 
582  if (writeCellDist)
583  {
584  const labelList& procIds = mesh.cellToProc();
585 
586  // Write decomposition as volScalarField for visualization
587  volScalarField cellDist
588  (
589  IOobject
590  (
591  "cellDist",
592  runTime.timeName(),
593  mesh,
594  IOobject::NO_READ,
595  IOobject::AUTO_WRITE
596  ),
597  mesh,
598  dimensionedScalar("cellDist", dimless, -1),
599  zeroGradientFvPatchScalarField::typeName
600  );
601 
602  forAll(procIds, celli)
603  {
604  cellDist[celli] = procIds[celli];
605  }
606 
607  cellDist.correctBoundaryConditions();
608  cellDist.write();
609 
610  Info<< nl << "Wrote decomposition as volScalarField to "
611  << cellDist.name() << " for visualization."
612  << endl;
613 
614  // Write decomposition as labelList for use with 'manual'
615  // decomposition method.
616  labelIOList cellDecomposition
617  (
618  IOobject
619  (
620  "cellDecomposition",
621  mesh.facesInstance(),
622  mesh,
623  IOobject::NO_READ,
624  IOobject::NO_WRITE,
625  false
626  ),
627  procIds
628  );
629  cellDecomposition.write();
630 
631  Info<< nl << "Wrote decomposition to "
632  << cellDecomposition.objectPath()
633  << " for use in manual decomposition." << endl;
634  }
635 
636  fileHandler().flush();
637  }
638 
639 
640  if (copyZero)
641  {
642  // Copy the 0 directory into each of the processor directories
643  fileName prevTimePath;
644  for (label proci = 0; proci < mesh.nProcs(); ++proci)
645  {
646  Time processorDb
647  (
648  Time::controlDictName,
649  args.rootPath(),
650  args.caseName()/("processor" + Foam::name(proci))
651  );
652  processorDb.setTime(runTime);
653 
654  if (fileHandler().isDir(runTime.timePath()))
655  {
656  // Get corresponding directory name (to handle processors/)
657  const fileName timePath
658  (
659  fileHandler().objectPath
660  (
661  IOobject
662  (
663  "",
664  processorDb.timeName(),
665  processorDb
666  ),
667  word::null
668  )
669  );
670 
671  if (timePath != prevTimePath)
672  {
673  Info<< "Processor " << proci
674  << ": copying " << runTime.timePath() << nl
675  << " to " << timePath << endl;
676  fileHandler().cp(runTime.timePath(), timePath);
677 
678  prevTimePath = timePath;
679  }
680  }
681  }
682  }
683  else
684  {
685  // Decompose the field files
686 
687  // Cached processor meshes and maps. These are only preserved if
688  // running with multiple times.
689  PtrList<Time> processorDbList(mesh.nProcs());
690  PtrList<fvMesh> procMeshList(mesh.nProcs());
691  PtrList<labelIOList> faceProcAddressingList(mesh.nProcs());
692  PtrList<labelIOList> cellProcAddressingList(mesh.nProcs());
693  PtrList<labelIOList> boundaryProcAddressingList(mesh.nProcs());
694  PtrList<fvFieldDecomposer> fieldDecomposerList(mesh.nProcs());
695  PtrList<dimFieldDecomposer> dimFieldDecomposerList(mesh.nProcs());
696  PtrList<labelIOList> pointProcAddressingList(mesh.nProcs());
697  PtrList<pointFieldDecomposer> pointFieldDecomposerList
698  (
699  mesh.nProcs()
700  );
701 
702 
703  // Loop over all times
704  forAll(times, timeI)
705  {
706  runTime.setTime(times[timeI], timeI);
707 
708  Info<< "Time = " << runTime.timeName() << endl;
709 
710  // Search for list of objects for this time
711  IOobjectList objects(mesh, runTime.timeName());
712 
713 
714  // Construct the vol fields
715  // ~~~~~~~~~~~~~~~~~~~~~~~~
716  PtrList<volScalarField> volScalarFields;
717  readFields(mesh, objects, volScalarFields, false);
718  PtrList<volVectorField> volVectorFields;
719  readFields(mesh, objects, volVectorFields, false);
720  PtrList<volSphericalTensorField> volSphericalTensorFields;
721  readFields(mesh, objects, volSphericalTensorFields, false);
722  PtrList<volSymmTensorField> volSymmTensorFields;
723  readFields(mesh, objects, volSymmTensorFields, false);
724  PtrList<volTensorField> volTensorFields;
725  readFields(mesh, objects, volTensorFields, false);
726 
727 
728  // Construct the dimensioned fields
729  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
730  PtrList<DimensionedField<scalar, volMesh>> dimScalarFields;
731  readFields(mesh, objects, dimScalarFields);
732  PtrList<DimensionedField<vector, volMesh>> dimVectorFields;
733  readFields(mesh, objects, dimVectorFields);
734  PtrList<DimensionedField<sphericalTensor, volMesh>>
735  dimSphericalTensorFields;
736  readFields(mesh, objects, dimSphericalTensorFields);
737  PtrList<DimensionedField<symmTensor, volMesh>>
738  dimSymmTensorFields;
739  readFields(mesh, objects, dimSymmTensorFields);
740  PtrList<DimensionedField<tensor, volMesh>> dimTensorFields;
741  readFields(mesh, objects, dimTensorFields);
742 
743 
744  // Construct the surface fields
745  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
746  PtrList<surfaceScalarField> surfaceScalarFields;
747  readFields(mesh, objects, surfaceScalarFields, false);
748  PtrList<surfaceVectorField> surfaceVectorFields;
749  readFields(mesh, objects, surfaceVectorFields, false);
750  PtrList<surfaceSphericalTensorField>
751  surfaceSphericalTensorFields;
752  readFields(mesh, objects, surfaceSphericalTensorFields, false);
753  PtrList<surfaceSymmTensorField> surfaceSymmTensorFields;
754  readFields(mesh, objects, surfaceSymmTensorFields, false);
755  PtrList<surfaceTensorField> surfaceTensorFields;
756  readFields(mesh, objects, surfaceTensorFields, false);
757 
758 
759  // Construct the point fields
760  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
761  const pointMesh& pMesh = pointMesh::New(mesh);
762 
763  PtrList<pointScalarField> pointScalarFields;
764  readFields(pMesh, objects, pointScalarFields, false);
765  PtrList<pointVectorField> pointVectorFields;
766  readFields(pMesh, objects, pointVectorFields, false);
767  PtrList<pointSphericalTensorField> pointSphericalTensorFields;
768  readFields(pMesh, objects, pointSphericalTensorFields, false);
769  PtrList<pointSymmTensorField> pointSymmTensorFields;
770  readFields(pMesh, objects, pointSymmTensorFields, false);
771  PtrList<pointTensorField> pointTensorFields;
772  readFields(pMesh, objects, pointTensorFields, false);
773 
774 
775  // Construct the Lagrangian fields
776  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
777 
778  fileNameList cloudDirs
779  (
781  (
782  runTime.timePath()/cloud::prefix,
783  fileName::DIRECTORY
784  )
785  );
786 
787  // Particles
788  PtrList<Cloud<indexedParticle>> lagrangianPositions
789  (
790  cloudDirs.size()
791  );
792  // Particles per cell
793  PtrList<List<SLList<indexedParticle*>*>> cellParticles
794  (
795  cloudDirs.size()
796  );
797 
798  PtrList<PtrList<labelIOField>> lagrangianLabelFields
799  (
800  cloudDirs.size()
801  );
802  PtrList<PtrList<labelFieldCompactIOField>>
803  lagrangianLabelFieldFields
804  (
805  cloudDirs.size()
806  );
807  PtrList<PtrList<scalarIOField>> lagrangianScalarFields
808  (
809  cloudDirs.size()
810  );
811  PtrList<PtrList<scalarFieldCompactIOField>>
812  lagrangianScalarFieldFields
813  (
814  cloudDirs.size()
815  );
816  PtrList<PtrList<vectorIOField>> lagrangianVectorFields
817  (
818  cloudDirs.size()
819  );
820  PtrList<PtrList<vectorFieldCompactIOField>>
821  lagrangianVectorFieldFields
822  (
823  cloudDirs.size()
824  );
825  PtrList<PtrList<sphericalTensorIOField>>
826  lagrangianSphericalTensorFields
827  (
828  cloudDirs.size()
829  );
830  PtrList<PtrList<sphericalTensorFieldCompactIOField>>
831  lagrangianSphericalTensorFieldFields(cloudDirs.size());
832  PtrList<PtrList<symmTensorIOField>> lagrangianSymmTensorFields
833  (
834  cloudDirs.size()
835  );
836  PtrList<PtrList<symmTensorFieldCompactIOField>>
837  lagrangianSymmTensorFieldFields
838  (
839  cloudDirs.size()
840  );
841  PtrList<PtrList<tensorIOField>> lagrangianTensorFields
842  (
843  cloudDirs.size()
844  );
845  PtrList<PtrList<tensorFieldCompactIOField>>
846  lagrangianTensorFieldFields
847  (
848  cloudDirs.size()
849  );
850 
851  label cloudI = 0;
852 
853  for (const fileName& cloudDir : cloudDirs)
854  {
855  IOobjectList cloudObjects
856  (
857  mesh,
858  runTime.timeName(),
859  cloud::prefix/cloudDir,
860  IOobject::MUST_READ,
861  IOobject::NO_WRITE,
862  false
863  );
864 
865  // Note: look up "positions" for backwards compatibility
866  if
867  (
868  cloudObjects.found("coordinates")
869  || cloudObjects.found("positions")
870  )
871  {
872  // Read lagrangian particles
873  // ~~~~~~~~~~~~~~~~~~~~~~~~~
874 
875  Info<< "Identified lagrangian data set: "
876  << cloudDir << endl;
877 
878  lagrangianPositions.set
879  (
880  cloudI,
881  new Cloud<indexedParticle>
882  (
883  mesh,
884  cloudDir,
885  false
886  )
887  );
888 
889 
890  // Sort particles per cell
891  // ~~~~~~~~~~~~~~~~~~~~~~~
892 
893  cellParticles.set
894  (
895  cloudI,
896  new List<SLList<indexedParticle*>*>
897  (
898  mesh.nCells(),
899  static_cast<SLList<indexedParticle*>*>(nullptr)
900  )
901  );
902 
903  label i = 0;
904 
905  for (indexedParticle& p : lagrangianPositions[cloudI])
906  {
907  p.index() = i++;
908 
909  label celli = p.cell();
910 
911  // Check
912  if (celli < 0 || celli >= mesh.nCells())
913  {
915  << "Illegal cell number " << celli
916  << " for particle with index "
917  << p.index()
918  << " at position "
919  << p.position() << nl
920  << "Cell number should be between 0 and "
921  << mesh.nCells()-1 << nl
922  << "On this mesh the particle should"
923  << " be in cell "
924  << mesh.findCell(p.position())
925  << exit(FatalError);
926  }
927 
928  if (!cellParticles[cloudI][celli])
929  {
930  cellParticles[cloudI][celli] =
931  new SLList<indexedParticle*>();
932  }
933 
934  cellParticles[cloudI][celli]->append(&p);
935  }
936 
937  // Read fields
938  // ~~~~~~~~~~~
939 
940  IOobjectList lagrangianObjects
941  (
942  mesh,
943  runTime.timeName(),
944  cloud::prefix/cloudDirs[cloudI],
945  IOobject::MUST_READ,
946  IOobject::NO_WRITE,
947  false
948  );
949 
951  (
952  cloudI,
953  lagrangianObjects,
954  lagrangianLabelFields
955  );
956 
957  lagrangianFieldDecomposer::readFieldFields
958  (
959  cloudI,
960  lagrangianObjects,
961  lagrangianLabelFieldFields
962  );
963 
965  (
966  cloudI,
967  lagrangianObjects,
968  lagrangianScalarFields
969  );
970 
971  lagrangianFieldDecomposer::readFieldFields
972  (
973  cloudI,
974  lagrangianObjects,
975  lagrangianScalarFieldFields
976  );
977 
979  (
980  cloudI,
981  lagrangianObjects,
982  lagrangianVectorFields
983  );
984 
985  lagrangianFieldDecomposer::readFieldFields
986  (
987  cloudI,
988  lagrangianObjects,
989  lagrangianVectorFieldFields
990  );
991 
993  (
994  cloudI,
995  lagrangianObjects,
996  lagrangianSphericalTensorFields
997  );
998 
999  lagrangianFieldDecomposer::readFieldFields
1000  (
1001  cloudI,
1002  lagrangianObjects,
1003  lagrangianSphericalTensorFieldFields
1004  );
1005 
1007  (
1008  cloudI,
1009  lagrangianObjects,
1010  lagrangianSymmTensorFields
1011  );
1012 
1013  lagrangianFieldDecomposer::readFieldFields
1014  (
1015  cloudI,
1016  lagrangianObjects,
1017  lagrangianSymmTensorFieldFields
1018  );
1019 
1021  (
1022  cloudI,
1023  lagrangianObjects,
1024  lagrangianTensorFields
1025  );
1026 
1027  lagrangianFieldDecomposer::readFieldFields
1028  (
1029  cloudI,
1030  lagrangianObjects,
1031  lagrangianTensorFieldFields
1032  );
1033 
1034  cloudI++;
1035  }
1036  }
1037 
1038  lagrangianPositions.setSize(cloudI);
1039  cellParticles.setSize(cloudI);
1040  lagrangianLabelFields.setSize(cloudI);
1041  lagrangianLabelFieldFields.setSize(cloudI);
1042  lagrangianScalarFields.setSize(cloudI);
1043  lagrangianScalarFieldFields.setSize(cloudI);
1044  lagrangianVectorFields.setSize(cloudI);
1045  lagrangianVectorFieldFields.setSize(cloudI);
1046  lagrangianSphericalTensorFields.setSize(cloudI);
1047  lagrangianSphericalTensorFieldFields.setSize(cloudI);
1048  lagrangianSymmTensorFields.setSize(cloudI);
1049  lagrangianSymmTensorFieldFields.setSize(cloudI);
1050  lagrangianTensorFields.setSize(cloudI);
1051  lagrangianTensorFieldFields.setSize(cloudI);
1052 
1053  Info<< endl;
1054 
1055  // split the fields over processors
1056  for (label proci = 0; proci < mesh.nProcs(); ++proci)
1057  {
1058  Info<< "Processor " << proci << ": field transfer" << endl;
1059 
1060 
1061  // open the database
1062  if (!processorDbList.set(proci))
1063  {
1064  processorDbList.set
1065  (
1066  proci,
1067  new Time
1068  (
1069  Time::controlDictName,
1070  args.rootPath(),
1071  args.caseName()
1072  / ("processor" + Foam::name(proci))
1073  )
1074  );
1075  }
1076  Time& processorDb = processorDbList[proci];
1077 
1078 
1079  processorDb.setTime(runTime);
1080 
1081  // read the mesh
1082  if (!procMeshList.set(proci))
1083  {
1084  procMeshList.set
1085  (
1086  proci,
1087  new fvMesh
1088  (
1089  IOobject
1090  (
1091  regionName,
1092  processorDb.timeName(),
1093  processorDb
1094  )
1095  )
1096  );
1097  }
1098  const fvMesh& procMesh = procMeshList[proci];
1099 
1100  const labelIOList& faceProcAddressing = procAddressing
1101  (
1102  procMeshList,
1103  proci,
1104  "faceProcAddressing",
1105  faceProcAddressingList
1106  );
1107 
1108  const labelIOList& cellProcAddressing = procAddressing
1109  (
1110  procMeshList,
1111  proci,
1112  "cellProcAddressing",
1113  cellProcAddressingList
1114  );
1115 
1116  const labelIOList& boundaryProcAddressing = procAddressing
1117  (
1118  procMeshList,
1119  proci,
1120  "boundaryProcAddressing",
1121  boundaryProcAddressingList
1122  );
1123 
1124 
1125  // FV fields
1126  {
1127  if (!fieldDecomposerList.set(proci))
1128  {
1129  fieldDecomposerList.set
1130  (
1131  proci,
1132  new fvFieldDecomposer
1133  (
1134  mesh,
1135  procMesh,
1137  cellProcAddressing,
1138  boundaryProcAddressing
1139  )
1140  );
1141  }
1142  const fvFieldDecomposer& fieldDecomposer =
1143  fieldDecomposerList[proci];
1144 
1145  fieldDecomposer.decomposeFields(volScalarFields);
1146  fieldDecomposer.decomposeFields(volVectorFields);
1147  fieldDecomposer.decomposeFields
1148  (
1149  volSphericalTensorFields
1150  );
1151  fieldDecomposer.decomposeFields(volSymmTensorFields);
1152  fieldDecomposer.decomposeFields(volTensorFields);
1153 
1154  fieldDecomposer.decomposeFields(surfaceScalarFields);
1155  fieldDecomposer.decomposeFields(surfaceVectorFields);
1156  fieldDecomposer.decomposeFields
1157  (
1158  surfaceSphericalTensorFields
1159  );
1160  fieldDecomposer.decomposeFields
1161  (
1162  surfaceSymmTensorFields
1163  );
1164  fieldDecomposer.decomposeFields(surfaceTensorFields);
1165 
1166  if (times.size() == 1)
1167  {
1168  // Clear cached decomposer
1169  fieldDecomposerList.set(proci, nullptr);
1170  }
1171  }
1172 
1173  // Dimensioned fields
1174  {
1175  if (!dimFieldDecomposerList.set(proci))
1176  {
1177  dimFieldDecomposerList.set
1178  (
1179  proci,
1180  new dimFieldDecomposer
1181  (
1182  mesh,
1183  procMesh,
1185  cellProcAddressing
1186  )
1187  );
1188  }
1189  const dimFieldDecomposer& dimDecomposer =
1190  dimFieldDecomposerList[proci];
1191 
1192  dimDecomposer.decomposeFields(dimScalarFields);
1193  dimDecomposer.decomposeFields(dimVectorFields);
1194  dimDecomposer.decomposeFields(dimSphericalTensorFields);
1195  dimDecomposer.decomposeFields(dimSymmTensorFields);
1196  dimDecomposer.decomposeFields(dimTensorFields);
1197 
1198  if (times.size() == 1)
1199  {
1200  dimFieldDecomposerList.set(proci, nullptr);
1201  }
1202  }
1203 
1204 
1205  // Point fields
1206  if
1207  (
1208  pointScalarFields.size()
1209  || pointVectorFields.size()
1210  || pointSphericalTensorFields.size()
1211  || pointSymmTensorFields.size()
1212  || pointTensorFields.size()
1213  )
1214  {
1215  const labelIOList& pointProcAddressing = procAddressing
1216  (
1217  procMeshList,
1218  proci,
1219  "pointProcAddressing",
1220  pointProcAddressingList
1221  );
1222 
1223  const pointMesh& procPMesh = pointMesh::New(procMesh);
1224 
1225  if (!pointFieldDecomposerList.set(proci))
1226  {
1227  pointFieldDecomposerList.set
1228  (
1229  proci,
1230  new pointFieldDecomposer
1231  (
1232  pMesh,
1233  procPMesh,
1234  pointProcAddressing,
1235  boundaryProcAddressing
1236  )
1237  );
1238  }
1239  const pointFieldDecomposer& pointDecomposer =
1240  pointFieldDecomposerList[proci];
1241 
1242  pointDecomposer.decomposeFields(pointScalarFields);
1243  pointDecomposer.decomposeFields(pointVectorFields);
1244  pointDecomposer.decomposeFields
1245  (
1246  pointSphericalTensorFields
1247  );
1248  pointDecomposer.decomposeFields(pointSymmTensorFields);
1249  pointDecomposer.decomposeFields(pointTensorFields);
1250 
1251 
1252  if (times.size() == 1)
1253  {
1254  pointProcAddressingList.set(proci, nullptr);
1255  pointFieldDecomposerList.set(proci, nullptr);
1256  }
1257  }
1258 
1259 
1260  // If there is lagrangian data write it out
1261  forAll(lagrangianPositions, cloudI)
1262  {
1263  if (lagrangianPositions[cloudI].size())
1264  {
1265  lagrangianFieldDecomposer fieldDecomposer
1266  (
1267  mesh,
1268  procMesh,
1270  cellProcAddressing,
1271  cloudDirs[cloudI],
1272  lagrangianPositions[cloudI],
1273  cellParticles[cloudI]
1274  );
1275 
1276  // Lagrangian fields
1277  {
1278  fieldDecomposer.decomposeFields
1279  (
1280  cloudDirs[cloudI],
1281  lagrangianLabelFields[cloudI]
1282  );
1283  fieldDecomposer.decomposeFieldFields
1284  (
1285  cloudDirs[cloudI],
1286  lagrangianLabelFieldFields[cloudI]
1287  );
1288  fieldDecomposer.decomposeFields
1289  (
1290  cloudDirs[cloudI],
1291  lagrangianScalarFields[cloudI]
1292  );
1293  fieldDecomposer.decomposeFieldFields
1294  (
1295  cloudDirs[cloudI],
1296  lagrangianScalarFieldFields[cloudI]
1297  );
1298  fieldDecomposer.decomposeFields
1299  (
1300  cloudDirs[cloudI],
1301  lagrangianVectorFields[cloudI]
1302  );
1303  fieldDecomposer.decomposeFieldFields
1304  (
1305  cloudDirs[cloudI],
1306  lagrangianVectorFieldFields[cloudI]
1307  );
1308  fieldDecomposer.decomposeFields
1309  (
1310  cloudDirs[cloudI],
1311  lagrangianSphericalTensorFields[cloudI]
1312  );
1313  fieldDecomposer.decomposeFieldFields
1314  (
1315  cloudDirs[cloudI],
1316  lagrangianSphericalTensorFieldFields[cloudI]
1317  );
1318  fieldDecomposer.decomposeFields
1319  (
1320  cloudDirs[cloudI],
1321  lagrangianSymmTensorFields[cloudI]
1322  );
1323  fieldDecomposer.decomposeFieldFields
1324  (
1325  cloudDirs[cloudI],
1326  lagrangianSymmTensorFieldFields[cloudI]
1327  );
1328  fieldDecomposer.decomposeFields
1329  (
1330  cloudDirs[cloudI],
1331  lagrangianTensorFields[cloudI]
1332  );
1333  fieldDecomposer.decomposeFieldFields
1334  (
1335  cloudDirs[cloudI],
1336  lagrangianTensorFieldFields[cloudI]
1337  );
1338  }
1339  }
1340  }
1341 
1342  // Decompose the "uniform" directory in the time region
1343  // directory
1344  decomposeUniform(copyUniform, mesh, processorDb, regionDir);
1345 
1346  // For a multi-region case, also decompose the "uniform"
1347  // directory in the time directory
1348  if (regionNames.size() > 1 && regioni == 0)
1349  {
1350  decomposeUniform(copyUniform, mesh, processorDb);
1351  }
1352 
1353  // We have cached all the constant mesh data for the current
1354  // processor. This is only important if running with
1355  // multiple times, otherwise it is just extra storage.
1356  if (times.size() == 1)
1357  {
1358  boundaryProcAddressingList.set(proci, nullptr);
1359  cellProcAddressingList.set(proci, nullptr);
1360  faceProcAddressingList.set(proci, nullptr);
1361  procMeshList.set(proci, nullptr);
1362  processorDbList.set(proci, nullptr);
1363  }
1364  }
1365 
1366  // Finite area mesh and field decomposition
1367 
1368  IOobject faMeshBoundaryIOobj
1369  (
1370  "faBoundary",
1371  mesh.time().findInstance
1372  (
1373  mesh.dbDir()/polyMesh::meshSubDir,
1374  "boundary"
1375  ),
1376  faMesh::meshSubDir,
1377  mesh,
1378  IOobject::READ_IF_PRESENT,
1379  IOobject::NO_WRITE
1380  );
1381 
1382 
1383  if (faMeshBoundaryIOobj.typeHeaderOk<faBoundaryMesh>(true))
1384  {
1385  Info << "\nFinite area mesh decomposition" << endl;
1386 
1387  faMeshDecomposition aMesh(mesh);
1388 
1389  aMesh.decomposeMesh();
1390 
1391  aMesh.writeDecomposition();
1392 
1393 
1394  // Construct the area fields
1395  // ~~~~~~~~~~~~~~~~~~~~~~~~
1396  PtrList<areaScalarField> areaScalarFields;
1397  readFields(aMesh, objects, areaScalarFields);
1398 
1399  PtrList<areaVectorField> areaVectorFields;
1400  readFields(aMesh, objects, areaVectorFields);
1401 
1402  PtrList<areaSphericalTensorField> areaSphericalTensorFields;
1403  readFields(aMesh, objects, areaSphericalTensorFields);
1404 
1405  PtrList<areaSymmTensorField> areaSymmTensorFields;
1406  readFields(aMesh, objects, areaSymmTensorFields);
1407 
1408  PtrList<areaTensorField> areaTensorFields;
1409  readFields(aMesh, objects, areaTensorFields);
1410 
1411 
1412  // Construct the edge fields
1413  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1414  PtrList<edgeScalarField> edgeScalarFields;
1415  readFields(aMesh, objects, edgeScalarFields);
1416 
1417  Info << endl;
1418 
1419  // Split the fields over processors
1420  for (label procI = 0; procI < mesh.nProcs(); procI++)
1421  {
1422  Info<< "Processor " << procI
1423  << ": finite area field transfer" << endl;
1424 
1425  // open the database
1426  Time processorDb
1427  (
1428  Time::controlDictName,
1429  args.rootPath(),
1430  args.caseName()
1431  / ("processor" + Foam::name(procI))
1432  );
1433 
1434  processorDb.setTime(runTime);
1435 
1436  // Read the mesh
1437  fvMesh procFvMesh
1438  (
1439  IOobject
1440  (
1441  regionName,
1442  processorDb.timeName(),
1443  processorDb
1444  )
1445  );
1446 
1447  faMesh procMesh(procFvMesh);
1448 
1449  // // Does not work. HJ, 15/Aug/2017
1450  // const labelIOList& faceProcAddressing =
1451  // procAddressing
1452  // (
1453  // procMeshList,
1454  // procI,
1455  // "faceProcAddressing",
1456  // faceProcAddressingList
1457  // );
1458 
1459  // const labelIOList& boundaryProcAddressing =
1460  // procAddressing
1461  // (
1462  // procMeshList,
1463  // procI,
1464  // "boundaryProcAddressing",
1465  // boundaryProcAddressingList
1466  // );
1467 
1469  (
1470  IOobject
1471  (
1472  "faceProcAddressing",
1473  "constant",
1474  procMesh.meshSubDir,
1475  procFvMesh,
1476  IOobject::MUST_READ,
1477  IOobject::NO_WRITE
1478  )
1479  );
1480 
1481  labelIOList boundaryProcAddressing
1482  (
1483  IOobject
1484  (
1485  "boundaryProcAddressing",
1486  "constant",
1487  procMesh.meshSubDir,
1488  procFvMesh,
1489  IOobject::MUST_READ,
1490  IOobject::NO_WRITE
1491  )
1492  );
1493 
1494  // FA fields
1495  if
1496  (
1497  areaScalarFields.size()
1498  || areaVectorFields.size()
1499  || areaSphericalTensorFields.size()
1500  || areaSymmTensorFields.size()
1501  || areaTensorFields.size()
1502  || edgeScalarFields.size()
1503  )
1504  {
1505  labelIOList edgeProcAddressing
1506  (
1507  IOobject
1508  (
1509  "edgeProcAddressing",
1510  "constant",
1511  procMesh.meshSubDir,
1512  procFvMesh,
1513  IOobject::MUST_READ,
1514  IOobject::NO_WRITE
1515  )
1516  );
1517 
1518  faFieldDecomposer fieldDecomposer
1519  (
1520  aMesh,
1521  procMesh,
1522  edgeProcAddressing,
1524  boundaryProcAddressing
1525  );
1526 
1527  fieldDecomposer.decomposeFields(areaScalarFields);
1528  fieldDecomposer.decomposeFields(areaVectorFields);
1529  fieldDecomposer.decomposeFields
1530  (
1531  areaSphericalTensorFields
1532  );
1533  fieldDecomposer.decomposeFields
1534  (
1535  areaSymmTensorFields
1536  );
1537  fieldDecomposer.decomposeFields(areaTensorFields);
1538 
1539  fieldDecomposer.decomposeFields(edgeScalarFields);
1540  }
1541  }
1542  }
1543  }
1544  }
1545  }
1546 
1547  Info<< "\nEnd\n" << endl;
1548 
1549  return 0;
1550 }
1551 
1552 
1553 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:130
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:71
runTime
engineTime & runTime
Definition: createEngineTime.H:13
faceProcAddressing
PtrList< labelIOList > & faceProcAddressing
Definition: checkFaceAddressingComp.H:9
OSspecific.H
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::dimless
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
Definition: dimensionSets.H:50
Foam::exists
bool exists(const fileName &name, const bool checkGzip=true, const bool followLink=true)
Does the name exist (as DIRECTORY or FILE) in the file system?
Definition: MSwindows.C:625
faCFD.H
regionProperties.H
Foam::argList::getOrDefault
T getOrDefault(const word &optName, const T &deflt) const
Get a value from the named option if present, or return default.
Definition: argListI.H:280
Foam::read
bool read(const char *buf, int32_t &val)
Same as readInt32.
Definition: int32.H:108
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::fileOperation::flush
virtual void flush() const
Forcibly wait until all output done. Flush any cached data.
Definition: fileOperation.C:1120
Foam::fileOperation::filePath
virtual fileName filePath(const bool checkGlobal, const IOobject &, const word &typeName, const bool search=true) const =0
Search for an object. checkGlobal : also check undecomposed case.
tensorIOField.H
Foam::fileHandler
const fileOperation & fileHandler()
Get current file handler.
Definition: fileOperation.C:1354
IOobjectList.H
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:350
scalarIOField.H
Foam::labelIOList
IOList< label > labelIOList
Label container classes.
Definition: labelIOList.H:44
Foam::instantList
List< instant > instantList
List of instants.
Definition: instantList.H:44
Foam::argList::get
T get(const label index) const
Get a value from the argument at index.
Definition: argListI.H:251
Foam::argList::rootPath
const fileName & rootPath() const
Return root path.
Definition: argListI.H:63
domainDecomposition.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
tensorFieldIOField.H
scalarFieldIOField.H
Foam::wordList
List< word > wordList
A List of words.
Definition: fileName.H:59
aMesh
faMesh aMesh(mesh)
Foam::flush
Ostream & flush(Ostream &os)
Flush stream.
Definition: Ostream.H:342
regionName
Foam::word regionName
Definition: createNamedDynamicFvMesh.H:1
domainDecompositionDryRun.H
faMeshDecomposition.H
vectorFieldIOField.H
Foam::Info
messageStream Info
Information stream (uses stdout - output is on the master only)
dimFieldDecomposer.H
Foam::name
word name(const complex &c)
Return string representation of complex.
Definition: complex.C:76
Foam::Time::timePath
fileName timePath() const
Return current time path.
Definition: Time.H:375
Foam::fileOperation::nProcs
virtual label nProcs(const fileName &dir, const fileName &local="") const
Get number of processor directories/results. Used for e.g.
Definition: fileOperation.C:1059
addRegionOption.H
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:43
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Foam::fileOperation::rmDir
virtual bool rmDir(const fileName &dir, const bool silent=false) const =0
Remove a directory and its contents.
labelFieldIOField.H
Foam::FatalError
error FatalError
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::flatOutput
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
Definition: FlatOutput.H:217
fvFieldDecomposer.H
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::readFields
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
Definition: ReadFieldsTemplates.C:312
Foam::fileOperation::ln
virtual bool ln(const fileName &src, const fileName &dst) const =0
Create a softlink. dst should not exist. Returns true if.
Foam::fileNameList
List< fileName > fileNameList
A List of fileNames.
Definition: fileNameList.H:58
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
vectorIOField.H
Foam::chDir
bool chDir(const fileName &dir)
Change current directory to the one specified and return true on success.
Definition: MSwindows.C:500
setRootCase.H
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:381
Foam::nl
constexpr char nl
Definition: Ostream.H:385
symmTensorFieldIOField.H
Foam::fileOperation::exists
bool exists(IOobject &io) const
Does ioobject exist. Is either a directory (empty name()) or.
Definition: fileOperation.C:522
pointFieldDecomposer.H
Foam::word::null
static const word null
An empty word.
Definition: word.H:77
createTime.H
Foam::List::set
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
Definition: List.H:325
forAllReverse
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Definition: stdFoam.H:309
symmTensorIOField.H
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:275
Foam::fileOperation::cp
virtual bool cp(const fileName &src, const fileName &dst, const bool followLink=true) const =0
Copy, recursively if necessary, the source to the destination.
Foam::cwd
fileName cwd()
The physical or logical current working directory path name.
Definition: MSwindows.C:468
fvCFD.H
decompositionModel.H
Foam::argList::caseName
const fileName & caseName() const
Return case name (parallel run) or global case (serial run)
Definition: argListI.H:69
sphericalTensorIOField.H
args
Foam::argList args(argc, argv)
Foam::mkDir
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
Definition: MSwindows.C:507
sphericalTensorFieldIOField.H
emptyFaPatch.H
Foam::readDir
fileNameList readDir(const fileName &directory, const fileName::Type type=fileName::FILE, const bool filtergz=true, const bool followLink=true)
Read a directory and return the entries as a fileName List.
Definition: MSwindows.C:707
lagrangianFieldDecomposer.H
pointFields.H
labelIOField.H
Foam::isDir
bool isDir(const fileName &name, const bool followLink=true)
Does the name exist as a DIRECTORY in the file system?
Definition: MSwindows.C:643
Foam::argList::found
bool found(const word &optName) const
Return true if the named option is found.
Definition: argListI.H:151
Foam::IOobject::MUST_READ
Definition: IOobject.H:120
faFieldDecomposer.H