Go to the documentation of this file.
36 template<
class CloudType>
50 if (Pstream::master())
53 mkDir(this->writeTimeDir());
61 this->writeTimeDir()/(
type() +
'_' + zoneName +
".dat")
66 <<
"# Source : " <<
type() <<
nl
67 <<
"# Face zone : " << zoneName <<
nl
68 <<
"# Faces : " << nFaces <<
nl
69 <<
"# Area : " << totArea <<
nl
70 <<
"# Time" <<
tab <<
"nParcels" <<
tab <<
"mass" <<
endl;
78 template<
class CloudType>
81 const typename parcelType::trackingData& td
84 Info<< this->modelName() <<
" output:" <<
nl;
91 const word& zoneName = fzm[faceZoneIDs_[i]].name();
96 Info<<
" faceZone " << zoneName
97 <<
": removed " << zoneNParcels
98 <<
" parcels with mass " << zoneMass
106 template<
class CloudType>
114 List<label> allZoneNParcels(faceZoneIDs_.size(), 0);
121 if (outputFilePtr_.set(i))
125 << allZoneNParcels[i] << token::TAB
126 << allZoneMass[i] <<
endl;
141 this->setModelProperty(
"mass", allZoneMass);
142 this->setModelProperty(
"nParcels", allZoneNParcels);
148 template<
class CloudType>
153 const word& modelName
160 typeId_(this->coeffDict().
template getOrDefault<label>(
"parcelType", -1)),
161 log_(this->coeffDict().getBool(
"log")),
162 resetOnWrite_(this->coeffDict().getBool(
"resetOnWrite")),
163 resetOnStart_(this->coeffDict().getBool(
"resetOnStart")),
166 const wordList faceZoneNames(this->coeffDict().
lookup(
"faceZones"));
168 nParcels_.setSize(faceZoneNames.size(), 0);
169 mass_.setSize(faceZoneNames.size(), 0.0);
171 if (!resetOnStart_ && Pstream::master())
173 this->getModelProperty(
"mass", mass_);
174 this->getModelProperty(
"nParcels", nParcels_);
177 outputFilePtr_.setSize(faceZoneNames.size());
186 const word& zoneName = faceZoneNames[i];
195 Info<<
" " << zoneName <<
" faces: " << nFaces <<
nl;
197 scalar totArea = 0.0;
198 for (label facei : fz)
202 totArea += magSf[facei];
207 label patchi = pbm.
patchID()[bFacei];
213 || refCast<const coupledPolyPatch>(pp).owner()
223 makeLogFile(zoneName, i, nFaces, totArea);
227 faceZoneIDs_.transfer(
zoneIDs);
231 template<
class CloudType>
238 faceZoneIDs_(rpf.faceZoneIDs_),
239 nParcels_(rpf.nParcels_),
241 typeId_(rpf.typeId_),
243 resetOnWrite_(rpf.resetOnWrite_),
244 resetOnStart_(rpf.resetOnStart_),
251 template<
class CloudType>
258 if ((typeId_ >= 0) && (
p.typeId() != typeId_))
266 this->owner().solution().output()
267 || this->owner().solution().transient()
274 const faceZone& fz = fzm[faceZoneIDs_[i]];
275 if (fz.found(
p.face()))
278 mass_[i] +=
p.mass()*
p.nParticle();
279 keepParticle =
false;
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
void write()
Write post-processing info.
label nInternalFaces() const
Number of internal faces.
virtual void postFace(const parcelType &p, bool &keepParticle)
Post-face hook.
static word timeName(const scalar t, const int precision=precision_)
Removes parcels that hit user-specified face zone faces.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const labelIOList & zoneIDs
#define forAll(list, i)
Loop across all elements in list.
const faceZoneMesh & faceZones() const
Return face zone mesh.
RemoveParcels(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
const labelList & patchID() const
Per boundary face label the patch index.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
A subset of mesh faces organised as a primitive patch.
const fvMesh & mesh() const
Return reference to the mesh.
messageStream Info
Information stream (uses stdout - output is on the master only)
A patch is a list of labels that address the faces in the global face list.
Lookup type of boundary radiation properties.
Templated base class for dsmc cloud.
virtual void postEvolve(const typename parcelType::trackingData &td)
Post-evolve hook.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh data needed to do the Finite Volume discretisation.
label findZoneID(const word &zoneName) const
Find zone index by name, return -1 if not found.
Output to file stream, using an OSstream.
#define DebugInfo
Report an information message using Foam::Info.
label whichFace(const label l) const
Return label of face in patch from global face label.
Templated cloud function object base class.
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
const Boundary & boundaryField() const
Return const-reference to the boundary field.