Go to the documentation of this file.
42 namespace functionObjects
58 static void writeList(Ostream&
os,
const string& header,
const UList<T>&
L)
61 os << header.c_str() <<
nl;
85 const fileName& commsDir,
86 const word& regionGroupName,
87 const wordRe& groupName
102 void Foam::functionObjects::externalCoupled::readColumns
105 const label nColumns,
106 autoPtr<IFstream>& masterFilePtr,
107 List<scalarField>& data
111 const globalIndex globalFaces(nRows);
123 List<scalarField>
values(nColumns);
126 const label procNRows = globalFaces.localSize(proci);
130 values[columni].setSize(procNRows);
133 for (label rowi = 0; rowi < procNRows; ++rowi)
138 if (!masterFilePtr().good())
141 <<
"Trying to read data for processor " << proci
143 <<
". Does your file have as many rows as there are"
144 <<
" patch faces (" << globalFaces.size()
148 masterFilePtr().getLine(line);
150 while (line.empty() || line[0] ==
'#');
152 IStringStream lineStr(line);
154 for (label columni = 0; columni < nColumns; ++columni)
156 lineStr >>
values[columni][rowi];
161 UOPstream str(proci, pBufs);
165 pBufs.finishedSends();
173 void Foam::functionObjects::externalCoupled::readLines
176 autoPtr<IFstream>& masterFilePtr,
181 const globalIndex globalFaces(nRows);
194 const label procNRows = globalFaces.localSize(proci);
196 UOPstream toProc(proci, pBufs);
198 for (label rowi = 0; rowi < procNRows; ++rowi)
203 if (!masterFilePtr().good())
206 <<
"Trying to read data for processor " << proci
208 <<
". Does your file have as many rows as there are"
209 <<
" patch faces (" << globalFaces.size()
213 masterFilePtr().getLine(line);
215 while (line.empty() || line[0] ==
'#');
224 pBufs.finishedSends();
228 for (label rowi = 0; rowi < nRows; ++rowi)
231 lines << line.c_str() <<
nl;
262 osPointsPtr() <<
"// Group: " << groupName <<
endl;
263 osFacesPtr() <<
"// Group: " << groupName <<
endl;
265 Info<< typeName <<
": writing geometry to " << dir <<
endl;
279 mesh.boundaryMesh().patchSet
285 for (
const label patchi : patchIDs)
289 mesh.globalData().mergePoints
300 collectedPoints[proci] =
pointField(
mesh.points(), uniquePointIDs);
304 faceList& patchFaces = collectedFaces[proci];
305 patchFaces =
p.localFaces();
319 allPoints.append(collectedPoints[proci]);
320 allFaces.
append(collectedFaces[proci]);
323 Info<< typeName <<
": mesh " <<
mesh.name()
324 <<
", patch " <<
p.name()
325 <<
": writing " <<
allPoints.size() <<
" points to "
326 << osPointsPtr().name() <<
nl
327 << typeName <<
": mesh " <<
mesh.name()
328 <<
", patch " <<
p.name()
329 <<
": writing " << allFaces.size() <<
" faces to "
330 << osFacesPtr().name() <<
endl;
333 const string entryHeader =
334 patchKey +
' ' +
mesh.name() +
' ' +
p.name();
337 writeList(osFacesPtr(), entryHeader, allFaces);
382 void Foam::functionObjects::externalCoupled::checkOrder
391 <<
"regionNames " <<
regionNames <<
" not in alphabetical order :"
397 void Foam::functionObjects::externalCoupled::initCoupling()
399 if (initialisedCoupling_)
405 forAll(regionGroupNames_, regioni)
407 const word& compName = regionGroupNames_[regioni];
417 const labelList& groups = regionToGroups_[compName];
419 for (
const label groupi : groups)
421 const wordRe& groupName = groupNames_[groupi];
423 bool geomExists =
false;
426 fileName dir(groupDir(commDirectory(), compName, groupName));
430 ||
isFile(dir/
"patchFaces");
451 initialisedCoupling_ =
true;
455 void Foam::functionObjects::externalCoupled::performCoupling()
468 const auto action = waitForSlave();
480 lastTrigger_ = time_.timeIndex();
485 action != time_.stopAt()
486 && action != Time::stopAtControls::saUnknown
489 Info<<
type() <<
": slave requested action "
492 time_.stopAt(action);
510 initialisedCoupling_(
false)
528 !initialisedCoupling_
529 || (time_.timeIndex() >= lastTrigger_ + calcFrequency_)
574 for (
const entry& dEntry : allRegionsDict)
576 if (!dEntry.isDict())
579 <<
"Regions must be specified in dictionary format"
583 const wordRe regionGroupName(dEntry.keyword());
590 regionGroupNames_.append(compositeName(
regionNames));
593 for (
const entry& dEntry : regionDict)
595 if (!dEntry.isDict())
598 <<
"Regions must be specified in dictionary format"
602 const wordRe groupName(dEntry.keyword());
605 const label nGroups = groupNames_.size();
609 auto fnd = regionToGroups_.find(regionGroupNames_.last());
612 fnd().append(nGroups);
616 regionToGroups_.insert
618 regionGroupNames_.last(),
622 groupNames_.append(groupName);
629 Info<<
type() <<
": Communicating with regions:" <<
endl;
630 for (
const word& compName : regionGroupNames_)
633 const labelList& groups = regionToGroups_[compName];
634 for (
const label groupi : groups)
636 const wordRe& groupName = groupNames_[groupi];
638 Info<<
indent <<
"patchGroup: " << groupName <<
"\t"
641 <<
indent <<
"Reading fields: "
642 << groupReadFields_[groupi]
644 <<
indent <<
"Writing fields: "
645 << groupWriteFields_[groupi]
658 for (
const word& compName : regionGroupNames_)
660 const labelList& groups = regionToGroups_[compName];
661 for (
const label groupi : groups)
663 const wordRe& groupName = groupNames_[groupi];
665 fileName dir(groupDir(commDirectory(), compName, groupName));
669 Log <<
type() <<
": creating communications directory "
683 forAll(regionGroupNames_, regioni)
685 const word& compName = regionGroupNames_[regioni];
695 const labelList& groups = regionToGroups_[compName];
697 for (
const label groupi : groups)
699 const wordRe& groupName = groupNames_[groupi];
706 readData<scalar>(
meshes, groupName, fieldName)
707 || readData<vector>(
meshes, groupName, fieldName)
708 || readData<sphericalTensor>(
meshes, groupName, fieldName)
709 || readData<symmTensor>(
meshes, groupName, fieldName)
710 || readData<tensor>(
meshes, groupName, fieldName)
716 <<
"Field " << fieldName <<
" in regions " << compName
717 <<
" was not found." <<
endl;
727 forAll(regionGroupNames_, regioni)
729 const word& compName = regionGroupNames_[regioni];
739 const labelList& groups = regionToGroups_[compName];
741 for (
const label groupi : groups)
743 const wordRe& groupName = groupNames_[groupi];
750 writeData<scalar>(
meshes, groupName, fieldName)
751 || writeData<vector>(
meshes, groupName, fieldName)
752 || writeData<sphericalTensor>(
meshes, groupName, fieldName)
753 || writeData<symmTensor>(
meshes, groupName, fieldName)
754 || writeData<tensor>(
meshes, groupName, fieldName)
760 <<
"Field " << fieldName <<
" in regions " << compName
761 <<
" was not found." <<
endl;
776 Log <<
type() <<
": removing data files written by master" <<
nl;
778 for (
const word& compName : regionGroupNames_)
780 const labelList& groups = regionToGroups_[compName];
781 for (
const label groupi : groups)
783 const wordRe& groupName = groupNames_[groupi];
790 groupDir(commDirectory(), compName, groupName)
806 Log <<
type() <<
": removing data files written by slave" <<
nl;
808 for (
const word& compName : regionGroupNames_)
810 const labelList& groups = regionToGroups_[compName];
811 for (
const label groupi : groups)
813 const wordRe& groupName = groupNames_[groupi];
820 groupDir(commDirectory(), compName, groupName)
A keyword and a list of tokens is an 'entry'.
List< label > labelList
A List of labels.
virtual void readDataMaster()
Read data files (all regions, all fields) on master (OpenFOAM)
vectorField pointField
pointField is a vectorField.
const vector L(dict.get< vector >("L"))
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
static constexpr int masterNo() noexcept
Process index of the master (always 0)
virtual bool read(const dictionary &dict)
Read and set the function object if its data have changed.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
A class for handling file names.
static word defaultRegion
Return the default region name.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
bool read(const char *buf, int32_t &val)
Same as readInt32.
Virtual base class for function objects with a reference to Time.
static word validate(const std::string &s, const bool prefix=false)
Construct validated word (no invalid characters).
A class representing the concept of 1 (one) that can be used to avoid manipulating objects known to b...
static bool master(const label communicator=worldComm)
Am I the master process.
static void scatter(const List< commsStruct > &comms, T &Value, const int tag, const label comm)
Scatter data. Distribute without modification. Reverse of gather.
bool isFile(const fileName &name, const bool checkGzip=true, const bool followLink=true)
Does the name exist as a FILE in the file system?
bool rm(const fileName &file)
Remove a file (or its gz equivalent), returning true if successful.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const T * set(const label i) const
Return const pointer to element (can be nullptr),.
A class for handling character strings derived from std::string.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
void clear() noexcept
Clear the addressed list, i.e. set the size to zero.
A wordRe is a Foam::word, but can contain a regular expression for matching words or strings.
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
#define forAll(list, i)
Loop across all elements in list.
List< word > wordList
A List of words.
static const Enum< stopAtControls > stopAtControlNames
Names for stopAtControls.
messageStream Info
Information stream (stdout output on master, null elsewhere)
A patch is a list of labels that address the faces in the global face list.
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
Foam::PtrList< Foam::fvMesh > meshes(regionNames.size())
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
void writeList(vtk::formatter &fmt, const UList< uint8_t > &values)
Write a list of uint8_t values.
virtual bool read(const dictionary &dict)
Read and set the function object if its data have changed.
virtual bool end()
Called when Time::run() determines that the time-loop exits.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const wordRes fieldNames(propsDict.getOrDefault< wordRes >("fields", wordRes()))
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
virtual void removeDataMaster() const
Remove data files written by master (OpenFOAM)
Mesh data needed to do the Finite Volume discretisation.
virtual void writeDataMaster() const
Write data files (all regions, all fields) from master (OpenFOAM)
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
labelList findStrings(const regExp &matcher, const UList< StringType > &input, const bool invert=false)
Return list indices for strings matching the regular expression.
errorManip< error > abort(error &err)
Ostream & indent(Ostream &os)
Indent stream.
List< Key > sortedToc() const
The table of contents (the keys) in sorted order.
static string patchKey
Name of patch key, e.g. '// Patch:' when looking for start of patch data.
errorManipArg< error, int > exit(error &err, const int errNo=1)
virtual bool end()
Called when Time::run() determines that the time-loop exits.
Output to file stream, using an OSstream.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
void reset(autoPtr< T > &&other) noexcept
Delete managed object and set to new given pointer.
virtual void removeDataSlave() const
Remove data files written by slave (external code)
virtual bool write()
Write, currently a no-op.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Encapsulates the logic for coordinating between OpenFOAM and an external application.
virtual bool execute()
Called at each ++ or += of the time-loop.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
static rangeType allProcs(const label communicator=worldComm)
Range of process indices for all processes.
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
static MinMax< T > ge(const T &minVal)
A semi-infinite range from minVal to the type max.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
A List of wordRe with additional matching capabilities.
defineTypeNameAndDebug(ObukhovLength, 0)
static const word null
An empty word.
Reads fields from the time directories and adds them to the mesh database for further post-processing...
externalCoupled(const word &name, const Time &runTime, const dictionary &dict)
Construct given time and dictionary.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Operations on lists of strings.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
bool readDict(const dictionary &dict)
Read communication settings from dictionary.
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
static constexpr direction size() noexcept
The number of elements in the VectorSpace = Ncmpts.
static void writeGeometry(const UPtrList< const fvMesh > &meshes, const fileName &commsDir, const wordRe &groupName)
Write geometry for the group as region/patch.
Begin list [isseparator].
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
static word compositeName(const wordList &)
#define WarningInFunction
Report a warning using Foam::Warning.
static label nProcs(const label communicator=worldComm)
Number of processes in parallel run, and 1 for serial run.
void writeFields(const fvMesh &mesh, const wordHashSet &selectedFields, const bool writeFaceFields)
bool isDir(const fileName &name, const bool followLink=true)
Does the name exist as a DIRECTORY in the file system?