56 const HashSet<word>& masterTimeDirSet,
63 if (!masterTimeDirSet.found(timeDirs[timei].name()))
73 int main(
int argc,
char *argv[])
77 "Reconstruct fields of a parallel case" 90 "specify a list of fields to be reconstructed. Eg, '(U T p)' - " 91 "regular expressions not currently supported" 96 "skip reconstructing fields" 102 "specify a list of lagrangian fields to be reconstructed. Eg, '(U d)' -" 103 "regular expressions not currently supported, " 104 "positions always included." 109 "skip reconstructing lagrangian positions and fields" 114 "skip reconstructing cellSets, faceSets, pointSets" 119 "only reconstruct new times (i.e. that do not exist already)" 125 HashSet<word> selectedFields;
135 Info<<
"Skipping reconstructing fields" 143 Info<<
"Skipping reconstructing lagrangian positions and fields" 150 if (noReconstructSets)
152 Info<<
"Skipping reconstructing cellSets, faceSets and pointSets" 157 HashSet<word> selectedLagrangianFields;
163 <<
"Cannot specify noLagrangian and lagrangianFields " 164 <<
"options together." 180 <<
"No processor* directories found" 185 const_cast<fileOperation&
>(
fileHandler()).setNProcs(nProcs);
188 PtrList<Time> databases(nProcs);
208 databases[0].times(),
218 if (timeDirs.empty())
232 HashSet<word> masterTimeDirSet(2*masterTimeDirs.size());
235 masterTimeDirSet.insert(masterTimeDirs[i].
name());
241 && regionNames.size() == 1
243 && haveAllTimes(masterTimeDirSet, timeDirs)
246 Info<<
"All times already reconstructed.\n\nEnd\n" <<
endl;
253 databases[proci].setTime(
runTime);
256 forAll(regionNames, regioni)
258 const word& regionName = regionNames[regioni];
261 Info<<
"\n\nReconstructing fields for mesh " << regionName <<
nl 277 processorMeshes procMeshes(databases, regionName);
287 if (newTimes && masterTimeDirSet.found(timeDirs[timei].name()))
289 Info<<
"Skipping time " << timeDirs[timei].name()
303 databases[proci].setTime(timeDirs[timei], timei);
315 procMeshes.reconstructPoints(mesh);
317 else if (meshStat != procStat)
320 <<
"readUpdate for the reconstructed mesh:" 322 <<
"readUpdate for the processor meshes :" 324 <<
"These should be equal or your addressing" 325 <<
" might be incorrect." 326 <<
" Please check your time directories for any " 327 <<
"mesh directories." <<
endl;
334 procMeshes.meshes()[0],
335 databases[0].timeName()
341 Info<<
"Reconstructing FV fields" <<
nl <<
endl;
343 fvFieldReconstructor fvReconstructor
347 procMeshes.faceProcAddressing(),
348 procMeshes.cellProcAddressing(),
349 procMeshes.boundaryProcAddressing()
352 fvReconstructor.reconstructFvVolumeInternalFields<scalar>
357 fvReconstructor.reconstructFvVolumeInternalFields<
vector>
362 fvReconstructor.reconstructFvVolumeInternalFields
368 fvReconstructor.reconstructFvVolumeInternalFields<
symmTensor>
373 fvReconstructor.reconstructFvVolumeInternalFields<
tensor>
379 fvReconstructor.reconstructFvVolumeFields<scalar>
384 fvReconstructor.reconstructFvVolumeFields<
vector>
394 fvReconstructor.reconstructFvVolumeFields<
symmTensor>
399 fvReconstructor.reconstructFvVolumeFields<
tensor>
405 fvReconstructor.reconstructFvSurfaceFields<scalar>
410 fvReconstructor.reconstructFvSurfaceFields<
vector>
420 fvReconstructor.reconstructFvSurfaceFields<
symmTensor>
425 fvReconstructor.reconstructFvSurfaceFields<
tensor>
431 if (fvReconstructor.nReconstructed() == 0)
439 Info<<
"Reconstructing point fields" <<
nl <<
endl;
442 PtrList<pointMesh> pMeshes(procMeshes.meshes().size());
449 new pointMesh(procMeshes.meshes()[proci])
453 pointFieldReconstructor pointReconstructor
457 procMeshes.pointProcAddressing(),
458 procMeshes.boundaryProcAddressing()
461 pointReconstructor.reconstructFields<scalar>
466 pointReconstructor.reconstructFields<
vector>
476 pointReconstructor.reconstructFields<
symmTensor>
481 pointReconstructor.reconstructFields<
tensor>
487 if (pointReconstructor.nReconstructed() == 0)
503 HashTable<IOobjectList> cloudObjects;
507 fileName lagrangianDir
511 databases[proci].timePath()
518 if (!lagrangianDir.empty())
532 cloudObjects.find(cloudDirs[i]);
534 if (iter == cloudObjects.end())
537 IOobjectList sprayObjs
539 procMeshes.meshes()[proci],
540 databases[proci].timeName(),
544 IOobject* positionsPtr =
545 sprayObjs.lookup(word(
"positions"));
549 cloudObjects.insert(cloudDirs[i], sprayObjs);
556 if (cloudObjects.size())
561 const word cloudName = string::validate<word>
567 const IOobjectList& sprayObjs = iter();
569 Info<<
"Reconstructing lagrangian fields for cloud " 570 << cloudName <<
nl <<
endl;
577 procMeshes.faceProcAddressing(),
578 procMeshes.cellProcAddressing()
580 reconstructLagrangianFields<label>
586 selectedLagrangianFields
588 reconstructLagrangianFieldFields<label>
594 selectedLagrangianFields
596 reconstructLagrangianFields<scalar>
602 selectedLagrangianFields
604 reconstructLagrangianFieldFields<scalar>
610 selectedLagrangianFields
612 reconstructLagrangianFields<vector>
618 selectedLagrangianFields
620 reconstructLagrangianFieldFields<vector>
626 selectedLagrangianFields
628 reconstructLagrangianFields<sphericalTensor>
634 selectedLagrangianFields
636 reconstructLagrangianFieldFields<sphericalTensor>
642 selectedLagrangianFields
644 reconstructLagrangianFields<symmTensor>
650 selectedLagrangianFields
652 reconstructLagrangianFieldFields<symmTensor>
658 selectedLagrangianFields
660 reconstructLagrangianFields<tensor>
666 selectedLagrangianFields
668 reconstructLagrangianFieldFields<tensor>
674 selectedLagrangianFields
685 if (!noReconstructSets)
688 HashTable<label> cSetNames;
689 HashTable<label> fSetNames;
690 HashTable<label> pSetNames;
692 forAll(procMeshes.meshes(), proci)
694 const fvMesh& procMesh = procMeshes.meshes()[proci];
707 IOobjectList cSets(objects.lookupClass(cellSet::typeName));
710 cSetNames.insert(iter.key(), cSetNames.size());
713 IOobjectList fSets(objects.lookupClass(faceSet::typeName));
716 fSetNames.insert(iter.key(), fSetNames.size());
718 IOobjectList pSets(objects.lookupClass(pointSet::typeName));
721 pSetNames.insert(iter.key(), pSetNames.size());
725 if (cSetNames.size() || fSetNames.size() || pSetNames.size())
728 PtrList<cellSet> cellSets(cSetNames.size());
729 PtrList<faceSet> faceSets(fSetNames.size());
730 PtrList<pointSet> pointSets(pSetNames.size());
732 Info<<
"Reconstructing sets:" <<
endl;
733 if (cSetNames.size())
736 << cSetNames.sortedToc() <<
endl;
738 if (fSetNames.size())
741 << fSetNames.sortedToc() <<
endl;
743 if (pSetNames.size())
746 << pSetNames.sortedToc() <<
endl;
750 forAll(procMeshes.meshes(), proci)
752 const fvMesh& procMesh = procMeshes.meshes()[proci];
763 procMeshes.cellProcAddressing()[proci];
767 objects.lookupClass(cellSet::typeName)
773 const cellSet procSet(*iter());
774 label setI = cSetNames[iter.key()];
775 if (!cellSets.set(setI))
788 cellSet& cSet = cellSets[setI];
793 cSet.insert(cellMap[iter.key()]);
799 procMeshes.faceProcAddressing()[proci];
803 objects.lookupClass(faceSet::typeName)
809 const faceSet procSet(*iter());
810 label setI = fSetNames[iter.key()];
811 if (!faceSets.set(setI))
824 faceSet& fSet = faceSets[setI];
829 fSet.insert(
mag(faceMap[iter.key()])-1);
834 procMeshes.pointProcAddressing()[proci];
838 objects.lookupClass(pointSet::typeName)
843 const pointSet propSet(*iter());
844 label setI = pSetNames[iter.key()];
845 if (!pointSets.set(setI))
858 pointSet& pSet = pointSets[setI];
863 pSet.insert(pointMap[iter.key()]);
879 pointSets[i].write();
887 PtrList<hexRef8Data> procData(procMeshes.meshes().size());
889 forAll(procMeshes.meshes(), procI)
891 const fvMesh& procMesh = procMeshes.meshes()[procI];
901 procMesh.time().timeName(),
914 const PtrList<labelIOList>& cellAddr =
915 procMeshes.cellProcAddressing();
917 UPtrList<const labelList> cellMaps(cellAddr.size());
920 cellMaps.set(i, &cellAddr[i]);
923 const PtrList<labelIOList>& pointAddr =
924 procMeshes.pointProcAddressing();
926 UPtrList<const labelList> pointMaps(pointAddr.size());
929 pointMaps.set(i, &pointAddr[i]);
932 UPtrList<const hexRef8Data> procRefs(procData.size());
935 procRefs.set(i, &procData[i]);
943 mesh.time().timeName(),
963 databases[0].timePath()/regionDir/
"uniform" 981 databases[0].timePath()/
"uniform"
List< instant > instantList
List of instants.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
fileName timePath() const
Return current time path.
const fileName & rootPath() const
Return root path.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
static word defaultRegion
Return the default region name.
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool optionFound(const word &opt) const
Return true if the named option is found.
static void noParallel()
Remove the parallel options.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Vector< scalar > vector
A scalar version of the templated Vector.
void reconstructLagrangianPositions(const polyMesh &mesh, const word &cloudName, const PtrList< fvMesh > &meshes, const PtrList< labelIOList > &faceProcAddressing, const PtrList< labelIOList > &cellProcAddressing)
static const pointMesh & New(const polyMesh &mesh)
virtual fileNameList readDir(const fileName &, const fileType=fileType::file, const bool filterVariants=true, const bool followLink=true) const =0
Read a directory and return the entries as a string list.
const word & regionDir(const word ®ionName)
SymmTensor< scalar > symmTensor
SymmTensor of scalars.
friend class const_iterator
Declare friendship with the const_iterator.
static void addOption(const word &opt, const string ¶m="", const string &usage="")
Add to an option to validOptions with usage information.
static const word null
An empty word.
List< label > labelList
A List of labels.
virtual void setTime(const Time &)
Reset the time and time-index to those of the given time.
const fileOperation & fileHandler()
Get current file handler.
static word controlDictName
The default control dictionary name (normally "controlDict")
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
virtual bool cp(const fileName &src, const fileName &dst, const bool followLink=true) const =0
Copy, recursively if necessary, the source to the destination.
instantList select(const instantList &) const
Select a list of Time values that are within the ranges.
const fileName & caseName() const
Return case name (parallel run) or global case (serial run)
static const word prefix
The prefix to local: lagrangian.
word name(const complex &)
Return a string representation of a complex.
fileName path() const
Return the path to the caseName.
List< word > wordList
A List of words.
virtual label nProcs(const fileName &dir, const fileName &local="") const
Get number of processor directories/results. Used for e.g.
#define WarningInFunction
Report a warning using Foam::Warning.
virtual bool isDir(const fileName &, const bool followLink=true) const =0
Does the name exist as a directory in the file system?
instantList times() const
Search the case for valid time directories.
dimensioned< scalar > mag(const dimensioned< Type > &)
const word cloudName(propsDict.lookup("cloudName"))
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
static void addNote(const string &)
Add extra notes for the usage information.
SphericalTensor< scalar > sphericalTensor
SphericalTensor of scalars.
readUpdateState
Enumeration defining the state of the mesh after a read update.
Foam::argList args(argc, argv)
List< fileName > fileNameList
A List of fileNames.
Tensor< scalar > tensor
Tensor of scalars.
wordList selectRegionNames(const argList &args, const Time &runTime)
static void addOptions(const bool constant=true, const bool withZero=false)
Add the options handled by timeSelector to argList::validOptions.
IStringStream optionLookup(const word &opt) const
Return an IStringStream from the named option.