FacePostProcessing.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2024 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "FacePostProcessing.H"
27 #include "Pstream.H"
28 #include "ListListOps.H"
29 #include "surfaceWriter.H"
30 #include "globalIndex.H"
31 #include "OSspecific.H"
32 
33 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
34 
35 template<class CloudType>
37 (
38  const word& zoneName,
39  const label zoneI,
40  const label nFaces,
41  const scalar totArea
42 )
43 {
44  // Create the output file if not already created
45  if (log_)
46  {
47  if (debug)
48  {
49  Info<< "Creating output file." << endl;
50  }
51 
52  if (Pstream::master())
53  {
54  // Create directory if does not exist
55  mkDir(this->writeTimeDir());
56 
57  // Open new file at start up
58  outputFilePtr_.set
59  (
60  zoneI,
61  new OFstream
62  (
63  this->writeTimeDir()/(type() + '_' + zoneName + ".dat")
64  )
65  );
66 
67  outputFilePtr_[zoneI]
68  << "# Source : " << type() << nl
69  << "# Face zone : " << zoneName << nl
70  << "# Faces : " << nFaces << nl
71  << "# Area : " << totArea << nl
72  << "# Time" << tab << "mass" << tab << "massFlowRate" << endl;
73  }
74  }
75 }
76 
77 
78 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
79 
80 template<class CloudType>
82 {
83  const fvMesh& mesh = this->owner().mesh();
84  const Time& time = mesh.time();
85  const faceZoneList& mfz = mesh.faceZones();
86  scalar timeNew = time.value();
87  scalar timeElapsed = timeNew - timeOld_;
88 
89  totalTime_ += timeElapsed;
90 
91  const scalar alpha = (totalTime_ - timeElapsed)/totalTime_;
92  const scalar beta = timeElapsed/totalTime_;
93 
94  forAll(faceZoneIndices_, zoneI)
95  {
96  massFlowRate_[zoneI] =
97  alpha*massFlowRate_[zoneI] + beta*mass_[zoneI]/timeElapsed;
98  massTotal_[zoneI] += mass_[zoneI];
99  }
100 
101  const label proci = Pstream::myProcNo();
102 
103  Info<< type() << " output:" << nl;
104 
105  List<scalarField> zoneMassTotal(mass_.size());
106  List<scalarField> zoneMassFlowRate(massFlowRate_.size());
107  forAll(faceZoneIndices_, zoneI)
108  {
109  const word& zoneName = mfz[faceZoneIndices_[zoneI]].name();
110 
111  scalarListList allProcMass(Pstream::nProcs());
112  allProcMass[proci] = massTotal_[zoneI];
113  Pstream::gatherList(allProcMass);
114  zoneMassTotal[zoneI] =
115  ListListOps::combine<scalarList>
116  (
117  allProcMass, accessOp<scalarList>()
118  );
119  const scalar sumMassTotal = sum(zoneMassTotal[zoneI]);
120 
121  scalarListList allProcMassFlowRate(Pstream::nProcs());
122  allProcMassFlowRate[proci] = massFlowRate_[zoneI];
123  Pstream::gatherList(allProcMassFlowRate);
124  zoneMassFlowRate[zoneI] =
125  ListListOps::combine<scalarList>
126  (
127  allProcMassFlowRate, accessOp<scalarList>()
128  );
129  const scalar sumMassFlowRate = sum(zoneMassFlowRate[zoneI]);
130 
131  Info<< " " << zoneName
132  << ": total mass = " << sumMassTotal
133  << "; average mass flow rate = " << sumMassFlowRate
134  << nl;
135 
136  if (outputFilePtr_.set(zoneI))
137  {
138  OFstream& os = outputFilePtr_[zoneI];
139  os << time.name() << token::TAB << sumMassTotal << token::TAB
140  << sumMassFlowRate<< endl;
141  }
142  }
143 
144  Info<< endl;
145 
146 
147  if (surfaceFormat_ != "none")
148  {
149  forAll(faceZoneIndices_, zoneI)
150  {
151  const faceZone& fZone = mfz[faceZoneIndices_[zoneI]];
152 
153  labelList pointToGlobal;
154  labelList uniqueMeshPointLabels;
155  autoPtr<globalIndex> globalPointsPtr =
156  mesh.globalData().mergePoints
157  (
158  fZone.patch().meshPoints(),
159  fZone.patch().meshPointMap(),
160  pointToGlobal,
161  uniqueMeshPointLabels
162  );
163 
164  pointField uniquePoints(mesh.points(), uniqueMeshPointLabels);
165  List<pointField> allProcPoints(Pstream::nProcs());
166  allProcPoints[proci] = uniquePoints;
167  Pstream::gatherList(allProcPoints);
168 
169  faceList faces(fZone.patch().localFaces());
170  forAll(faces, i)
171  {
172  inplaceRenumber(pointToGlobal, faces[i]);
173  }
174  List<faceList> allProcFaces(Pstream::nProcs());
175  allProcFaces[proci] = faces;
176  Pstream::gatherList(allProcFaces);
177 
178  if (Pstream::master())
179  {
180  pointField allPoints
181  (
182  ListListOps::combine<pointField>
183  (
184  allProcPoints, accessOp<pointField>()
185  )
186  );
187 
188  faceList allFaces
189  (
190  ListListOps::combine<faceList>
191  (
192  allProcFaces, accessOp<faceList>()
193  )
194  );
195 
197  (
198  surfaceWriter::New(surfaceFormat_, this->coeffDict())
199  );
200 
201  writer->write
202  (
203  this->writeTimeDir(),
204  fZone.name(),
205  allPoints,
206  allFaces,
207  false,
208  "massTotal",
209  zoneMassTotal[zoneI],
210  "massFlowRate",
211  zoneMassFlowRate[zoneI]
212  );
213  }
214  }
215  }
216 
217 
218  if (resetOnWrite_)
219  {
220  forAll(faceZoneIndices_, zoneI)
221  {
222  massFlowRate_[zoneI] = 0.0;
223  }
224  timeOld_ = timeNew;
225  totalTime_ = 0.0;
226  }
227 
228  forAll(mass_, zoneI)
229  {
230  mass_[zoneI] = 0.0;
231  }
232 
233  // writeProperties();
234 }
235 
236 
237 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
238 
239 template<class CloudType>
241 (
242  const dictionary& dict,
243  CloudType& owner,
244  const word& modelName
245 )
246 :
247  CloudFunctionObject<CloudType>(dict, owner, modelName, typeName),
248  faceZoneIndices_(),
249  surfaceFormat_(this->coeffDict().lookup("surfaceFormat")),
250  resetOnWrite_(this->coeffDict().lookup("resetOnWrite")),
251  totalTime_(0.0),
252  mass_(),
253  massTotal_(),
254  massFlowRate_(),
255  log_(this->coeffDict().lookup("log")),
256  outputFilePtr_(),
257  timeOld_(owner.mesh().time().value())
258 {
259  wordList faceZoneNames(this->coeffDict().lookup("faceZones"));
260  mass_.setSize(faceZoneNames.size());
261  massTotal_.setSize(faceZoneNames.size());
262  massFlowRate_.setSize(faceZoneNames.size());
263 
264  outputFilePtr_.setSize(faceZoneNames.size());
265 
266  DynamicList<label> zoneIDs;
267  const faceZoneList& mfz = owner.mesh().faceZones();
268  const surfaceScalarField& magSf = owner.mesh().magSf();
269  const polyBoundaryMesh& pbm = owner.mesh().boundaryMesh();
270  forAll(faceZoneNames, i)
271  {
272  const word& zoneName = faceZoneNames[i];
273  label zoneI = mfz.findIndex(zoneName);
274  if (zoneI != -1)
275  {
276  zoneIDs.append(zoneI);
277  const faceZone& fz = mfz[zoneI];
278  mass_[i].setSize(fz.size(), 0.0);
279  massTotal_[i].setSize(fz.size(), 0.0);
280  massFlowRate_[i].setSize(fz.size(), 0.0);
281 
282  label nFaces = returnReduce(fz.size(), sumOp<label>());
283  Info<< " " << zoneName << " faces: " << nFaces << nl;
284 
285  scalar totArea = 0.0;
286  forAll(fz, j)
287  {
288  label facei = fz[j];
289  if (facei < owner.mesh().nInternalFaces())
290  {
291  totArea += magSf[fz[j]];
292  }
293  else
294  {
295  label bFacei = facei - owner.mesh().nInternalFaces();
296  label patchi = pbm.patchIndices()[bFacei];
297  const polyPatch& pp = pbm[patchi];
298 
299  if
300  (
301  !magSf.boundaryField()[patchi].coupled()
302  || refCast<const coupledPolyPatch>(pp).owner()
303  )
304  {
305  label localFacei = pp.whichFace(facei);
306  totArea += magSf.boundaryField()[patchi][localFacei];
307  }
308  }
309  }
310  totArea = returnReduce(totArea, sumOp<scalar>());
311 
312  makeLogFile(zoneName, i, nFaces, totArea);
313  }
314  }
315 
316  faceZoneIndices_.transfer(zoneIDs);
317 
318  // readProperties(); AND initialise mass... fields
319 }
320 
321 
322 template<class CloudType>
324 (
326 )
327 :
329  faceZoneIndices_(pff.faceZoneIndices_),
330  surfaceFormat_(pff.surfaceFormat_),
331  resetOnWrite_(pff.resetOnWrite_),
332  totalTime_(pff.totalTime_),
333  mass_(pff.mass_),
334  massTotal_(pff.massTotal_),
335  massFlowRate_(pff.massFlowRate_),
336  log_(pff.log_),
337  outputFilePtr_(),
338  timeOld_(0.0)
339 {}
340 
341 
342 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
343 
344 template<class CloudType>
346 {}
347 
348 
349 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
350 
351 template<class CloudType>
353 {
354  if
355  (
356  this->owner().solution().output()
357  || this->owner().solution().transient()
358  )
359  {
360  const faceZoneList& mfz = this->owner().mesh().faceZones();
361 
362  forAll(faceZoneIndices_, i)
363  {
364  const faceZone& fz = mfz[faceZoneIndices_[i]];
365 
366  label faceId = -1;
367  forAll(fz, j)
368  {
369  if (fz[j] == p.face())
370  {
371  faceId = j;
372  break;
373  }
374  }
375 
376  if (faceId != -1)
377  {
378  mass_[i][faceId] += p.mass()*p.nParticle();
379  }
380  }
381  }
382 }
383 
384 
385 template<class CloudType>
387 (
388  const parcelType& p,
389  const polyPatch& pp
390 )
391 {
392  if (pp.coupled())
393  {
394  preFace(p);
395  }
396 }
397 
398 
399 // ************************************************************************* //
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Templated cloud function object base class.
const CloudType & owner() const
Return const access to the owner cloud.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
const fvMesh & mesh() const
Return references to the mesh.
Definition: DSMCCloudI.H:41
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
Records particle face quantities on used-specified face zone.
virtual void postPatch(const parcelType &p, const polyPatch &pp)
Post-patch hook.
virtual void preFace(const parcelType &p)
Post-face hook.
void write()
Write post-processing info.
virtual ~FacePostProcessing()
Destructor.
FacePostProcessing(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const word & name() const
Return name.
Definition: IOobject.H:310
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void setSize(const label)
Reset size of List.
Definition: List.C:281
Output to file stream.
Definition: OFstream.H:86
const Map< label > & meshPointMap() const
Mesh point map. Given the global point index find its.
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
const List< FaceType > & localFaces() const
Return patch faces addressing into local point list.
label findIndex(const word &key) const
Return the index of the given the key or -1 if not found.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
const MeshType & mesh() const
Return the mesh reference.
Definition: ZoneList.H:109
const word & name() const
Return name.
Definition: Zone.H:184
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
const Type & value() const
Return const reference to value.
const word & name() const
Return const reference to name.
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:59
const primitiveFacePatch & patch() const
Return reference to primitive patch.
Definition: faceZone.C:315
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:418
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
const polyMesh & mesh() const
Return reference to polyMesh.
Definition: fvMesh.H:441
autoPtr< globalIndex > mergePoints(labelList &pointToGlobal, labelList &uniquePoints) const
Helper for merging (collocated!) mesh point data.
Foam::polyBoundaryMesh.
const labelList & patchIndices() const
Boundary face patch indices.
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1515
const faceZoneList & faceZones() const
Return face zones.
Definition: polyMesh.H:446
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:404
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1313
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:360
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
Definition: polyPatch.H:290
label nInternalFaces() const
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:51
const dictionary & coeffDict() const
Return const access to the coefficients dictionary.
Definition: subModelBase.C:128
A class for handling words, derived from string.
Definition: word.H:62
label patchi
label faceId(-1)
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
static const char tab
Definition: Ostream.H:265
messageStream Info
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
void inplaceRenumber(const labelUList &oldToNew, ListType &)
Inplace renumber the values of a list.
static const char nl
Definition: Ostream.H:266
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict
volScalarField & p
mkDir(pdfPath)