Go to the documentation of this file.
37 template<
class ParticleType>
43 template<
class ParticleType>
50 "uniform"/cloud::prefix/
name(),
52 IOobject::MUST_READ_IF_MODIFIED,
57 if (dictObj.typeHeaderOk<IOdictionary>(
true))
59 const IOdictionary uniformPropsDict(dictObj);
64 cloud::geometryTypeNames.getOrDefault
68 cloud::geometryType::POSITIONS
71 const word procName(
"processor" +
Foam::name(Pstream::myProcNo()));
73 const dictionary* dictptr = uniformPropsDict.findDict(procName);
77 dictptr->readEntry(
"particleCount", ParticleType::particleCount_);
82 ParticleType::particleCount_ = 0;
87 template<
class ParticleType>
90 IOdictionary uniformPropsDict
96 "uniform"/cloud::prefix/
name(),
105 np[Pstream::myProcNo()] = ParticleType::particleCount_;
107 Pstream::listCombineGather(np, maxEqOp<label>());
108 Pstream::listCombineScatter(np);
113 cloud::geometryTypeNames[geometryType_]
119 uniformPropsDict.add(procName, dictionary());
120 uniformPropsDict.subDict(procName).add(
"particleCount", np[i]);
123 uniformPropsDict.writeObject
125 IOstreamOption(IOstream::ASCII, time().writeCompression()),
131 template<
class ParticleType>
134 readCloudUniformProperties();
136 IOPosition<Cloud<ParticleType>> ioP(*
this, geometryType_);
138 const bool valid = ioP.headerOk();
139 Istream& is = ioP.readStream(checkClass ? typeName :
"", valid);
142 ioP.readData(is, *
this);
148 Pout<<
"Cannot read particle positions file:" <<
nl
149 <<
" " << ioP.objectPath() <<
nl
150 <<
"Assuming the initial cloud contains 0 particles." <<
endl;
154 geometryType_ = cloud::geometryType::COORDINATES;
159 polyMesh_.tetBasePtIs();
165 template<
class ParticleType>
170 const bool checkClass
177 geometryType_(cloud::geometryType::COORDINATES)
181 polyMesh_.tetBasePtIs();
182 polyMesh_.oldCellCentres();
184 initCloud(checkClass);
190 template<
class ParticleType>
193 const word& fieldName,
209 template<
class ParticleType>
210 template<
class DataType>
217 if (
data.size() !=
c.size())
221 <<
" field " <<
data.size()
222 <<
" does not match the number of particles " <<
c.size()
228 template<
class ParticleType>
229 template<
class DataType>
236 if (
data.size() !=
c.size())
240 <<
" field " <<
data.size()
241 <<
" does not match the number of particles " <<
c.size()
247 template<
class ParticleType>
259 return fldNewPtr->store();
266 template<
class ParticleType>
285 if (selectFields.size() && !selectFields.
match(iter()->
name()))
299 auto&
object = *iter();
303 readStoreFile<label>(
object, ioNew)
304 || readStoreFile<scalar>(
object, ioNew)
305 || readStoreFile<vector>(
object, ioNew)
306 || readStoreFile<sphericalTensor>(
object, ioNew)
307 || readStoreFile<symmTensor>(
object, ioNew)
308 || readStoreFile<tensor>(
object, ioNew)
314 <<
"Unhandled field type " << iter()->headerClassName()
321 template<
class ParticleType>
328 template<
class ParticleType>
335 writeCloudUniformProperties();
338 return cloud::writeObject(streamOpt, this->size());
344 template<
class ParticleType>
int debug
Static debugging option.
List< label > labelList
A List of labels.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual void writeFields() const
Write the field data for the cloud of particles Dummy at.
const word cloudName(propsDict.get< word >("cloud"))
A class for handling words, derived from Foam::string.
A primitive field of type <T> with automated input and output.
static constexpr const zero Zero
Global zero (0)
bool readStoreFile(const IOobject &io, const IOobject &ioNew) const
Helper function to store a cloud field on its registry.
void checkFieldIOobject(const Cloud< ParticleType > &c, const IOField< DataType > &data) const
Check lagrangian data field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
A Field of objects of type <T> with automated input and output using a compact storage....
const word & headerClassName() const noexcept
Return name of the class name read from header.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
const word & name() const
IOobject fieldIOobject(const word &fieldName, const IOobject::readOption r) const
Helper to construct IOobject for field and current time.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
bool match(const std::string &text, bool literal=false) const
Smart match as literal or regex, stopping on the first match.
Registry of regIOobjects.
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
Generic templated field type.
The IOstreamOption is a simple container for options an IOstream can normally have.
static word cloudPropertiesName
Name of cloud properties dictionary.
void readFromFiles(objectRegistry &obr, const wordRes &selectFields) const
Read from files into objectRegistry.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
#define forAllIters(container, iter)
Iterate across all elements in the container object.
OBJstream os(runTime.globalPath()/outputName)
errorManip< error > abort(error &err)
List of IOobjects with searching and retrieving facilities.
A cloud is a registry collection of lagrangian particles.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
Base cloud calls templated on particle type.
A List of wordRe with additional matching capabilities.
Cloud(const polyMesh &mesh, const word &cloudName, const IDLList< ParticleType > &particles)
Construct from mesh and a list of particles.
const dimensionedScalar c
Speed of light in a vacuum.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
readOption
Enumeration defining the read options.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Database for solution data, solver performance and other reduced data.
void writeFields(const fvMesh &mesh, const wordHashSet &selectedFields, const bool writeFaceFields)
void checkFieldFieldIOobject(const Cloud< ParticleType > &c, const CompactIOField< Field< DataType >, DataType > &data) const
Check lagrangian data fieldfield.
virtual bool writeObject(IOstreamOption streamOpt, const bool valid) const
Write using stream options.