37 template<
class CloudType>
51 if (Pstream::master())
54 mkDir(this->writeTimeDir());
62 this->writeTimeDir()/(
type() +
'_' + zoneName +
".dat")
67 <<
"# Source : " <<
type() <<
nl 68 <<
"# Face zone : " << zoneName <<
nl 69 <<
"# Faces : " << nFaces <<
nl 70 <<
"# Area : " << totArea <<
nl 71 <<
"# Time" <<
tab <<
"mass" <<
tab <<
"massFlowRate" <<
endl;
79 template<
class CloudType>
82 const fvMesh&
mesh = this->owner().mesh();
83 const Time& time =
mesh.time();
85 scalar timeNew = time.value();
86 scalar timeElapsed = timeNew - timeOld_;
88 totalTime_ += timeElapsed;
90 const scalar
alpha = (totalTime_ - timeElapsed)/totalTime_;
91 const scalar
beta = timeElapsed/totalTime_;
93 forAll(faceZoneIDs_, zoneI)
95 massFlowRate_[zoneI] =
96 alpha*massFlowRate_[zoneI] +
beta*mass_[zoneI]/timeElapsed;
97 massTotal_[zoneI] += mass_[zoneI];
100 const label proci = Pstream::myProcNo();
104 List<scalarField> zoneMassTotal(mass_.size());
105 List<scalarField> zoneMassFlowRate(massFlowRate_.size());
106 forAll(faceZoneIDs_, zoneI)
108 const word& zoneName = fzm[faceZoneIDs_[zoneI]].name();
111 allProcMass[proci] = massTotal_[zoneI];
112 Pstream::gatherList(allProcMass);
113 zoneMassTotal[zoneI] =
114 ListListOps::combine<scalarList>
116 allProcMass, accessOp<scalarList>()
118 const scalar sumMassTotal =
sum(zoneMassTotal[zoneI]);
121 allProcMassFlowRate[proci] = massFlowRate_[zoneI];
122 Pstream::gatherList(allProcMassFlowRate);
123 zoneMassFlowRate[zoneI] =
124 ListListOps::combine<scalarList>
126 allProcMassFlowRate, accessOp<scalarList>()
128 const scalar sumMassFlowRate =
sum(zoneMassFlowRate[zoneI]);
130 Log_<<
" " << zoneName
131 <<
": total mass = " << sumMassTotal
132 <<
"; average mass flow rate = " << sumMassFlowRate
135 if (outputFilePtr_.set(zoneI))
137 OFstream&
os = outputFilePtr_[zoneI];
138 os << time.timeName() << token::TAB << sumMassTotal << token::TAB
139 << sumMassFlowRate<<
endl;
146 if (surfaceFormat_ !=
"none")
148 forAll(faceZoneIDs_, zoneI)
150 const faceZone& fZone = fzm[faceZoneIDs_[zoneI]];
154 autoPtr<globalIndex> globalPointsPtr =
155 mesh.globalData().mergePoints
157 fZone().meshPoints(),
158 fZone().meshPointMap(),
160 uniqueMeshPointLabels
164 List<pointField> allProcPoints(Pstream::nProcs());
165 allProcPoints[proci] = uniquePoints;
166 Pstream::gatherList(allProcPoints);
168 faceList faces(fZone().localFaces());
173 List<faceList> allProcFaces(Pstream::nProcs());
174 allProcFaces[proci] = faces;
175 Pstream::gatherList(allProcFaces);
177 if (Pstream::master())
181 ListListOps::combine<pointField>
183 allProcPoints, accessOp<pointField>()
189 ListListOps::combine<faceList>
191 allProcFaces, accessOp<faceList>()
214 (this->writeTimeDir() / fZone.name()),
219 writer->write(
"massTotal", zoneMassTotal[zoneI]);
220 writer->write(
"massFlowRate", zoneMassFlowRate[zoneI]);
228 forAll(faceZoneIDs_, zoneI)
230 massFlowRate_[zoneI] = 0.0;
247 template<
class CloudType>
252 const word& modelName
257 surfaceFormat_(this->coeffDict().getWord(
"surfaceFormat")),
258 resetOnWrite_(this->coeffDict().getBool(
"resetOnWrite")),
259 logToFile_(this->coeffDict().getBool(
"log")),
265 timeOld_(owner.
mesh().time().value())
272 outputFilePtr_.setSize(faceZoneNames.
size());
281 const word& zoneName = faceZoneNames[i];
293 Info<<
" " << zoneName <<
" faces: " << nFaces <<
nl;
295 scalar totArea = 0.0;
296 for (
const label facei : fz)
300 totArea += magSf[facei];
304 const label patchi =
pbm.patchID(facei);
310 || refCast<const coupledPolyPatch>(
pp).
owner()
313 label localFacei =
pp.whichFace(facei);
320 makeLogFile(zoneName, i, nFaces, totArea);
324 faceZoneIDs_.transfer(
zoneIDs);
330 template<
class CloudType>
337 faceZoneIDs_(pff.faceZoneIDs_),
338 surfaceFormat_(pff.surfaceFormat_),
339 resetOnWrite_(pff.resetOnWrite_),
340 logToFile_(pff.logToFile_),
341 totalTime_(pff.totalTime_),
343 massTotal_(pff.massTotal_),
344 massFlowRate_(pff.massFlowRate_),
352 template<
class CloudType>
356 const typename parcelType::trackingData&
td 362 || this->owner().
solution().
transient()
369 const faceZone& fz = fzm[faceZoneIDs_[i]];
374 if (fz[j] ==
p.face())
383 mass_[i][
faceId] +=
p.mass()*
p.nParticle();
const polyBoundaryMesh & pbm
void size(const label n)
Older name for setAddressableSize.
ZoneMesh< faceZone, polyMesh > faceZoneMesh
A ZoneMesh with the type faceZone.
const labelIOList & zoneIDs
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &f1)
void write()
Write post-processing info.
vtk::lineWriter writer(edgeCentres, edgeList::null(), fileName(aMesh.time().globalPath()/"finiteArea-edgesCentres"))
A list of keyword definitions, which are a keyword followed by a number of values (eg...
wallPoints::trackData td(isBlockedFace, regionToBlockSize)
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
FacePostProcessing(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
constexpr char tab
The tab '\t' character(0x09)
Lookup type of boundary radiation properties.
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Perform reduction on a copy, using specified binary operation.
const CloudType & owner() const
Return const access to the owner cloud.
#define forAll(list, i)
Loop across all elements in list.
void inplaceRenumber(const labelUList &oldToNew, IntListType &lists)
Inplace renumber the values (not the indices) of a list of lists.
List< face > faceList
List of faces.
vectorField pointField
pointField is a vectorField.
void setSize(const label n)
Alias for resize()
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
#define Log_
Report write to Foam::Info if the class log switch is true.
label findZoneID(const word &zoneName) const
Find zone index by name, return -1 if not found.
List< scalarList > scalarListList
List of scalarList.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
A class for handling words, derived from Foam::string.
Records particle face quantities on used-specified face zone.
label nInternalFaces() const noexcept
Number of internal faces.
const dictionary & coeffDict() const
Return const access to the coefficients dictionary.
#define DebugInfo
Report an information message using Foam::Info.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO...
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
int debug
Static debugging option.
OBJstream os(runTime.globalPath()/outputName)
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
virtual bool postFace(const parcelType &p, const typename parcelType::trackingData &td)
Post-face hook.
const fvMesh & mesh() const
Return reference to the mesh.
messageStream Info
Information stream (stdout output on master, null elsewhere)
static Ostream & output(Ostream &os, const IntRange< T > &range)
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
Selector class for relaxation factors, solver type and solution.
A subset of mesh faces organised as a primitive patch.
List< label > labelList
A List of labels.
A patch is a list of labels that address the faces in the global face list.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Templated base class for dsmc cloud.
Templated cloud function object base class.
const Boundary & boundaryField() const noexcept
Return const-reference to the boundary field.
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
const dictionary formatOptions(propsDict.subOrEmptyDict("formatOptions", keyType::LITERAL))