populationBalanceSizeDistribution.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) 2017-2023 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 
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace functionObjects
34 {
37  (
41  );
42 }
43 }
44 
45 template<>
46 const char*
48 <
50  6
51 >::names[] =
52 {
53  "numberConcentration",
54  "numberDensity",
55  "volumeConcentration",
56  "volumeDensity",
57  "areaConcentration",
58  "areaDensity"
59 };
60 
61 const Foam::NamedEnum
62 <
64  6
66 
67 template<>
68 const char*
70 <
72  4
73 >::names[] =
74 {
75  "volume",
76  "area",
77  "diameter",
78  "projectedAreaDiameter"
79 };
80 
81 const Foam::NamedEnum
82 <
84  4
87 
88 
89 namespace Foam
90 {
91  template<>
92  const char* NamedEnum
93  <
95  4
96  >::names[] =
97  {
98  "numberConcentration",
99  "volumeConcentration",
100  "areaConcentration",
101  "cellVolume"
102  };
103 }
104 
105 
106 const Foam::NamedEnum
107 <
109  4
110 >
112 
114 
115 
116 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
117 
119 Foam::functionObjects::populationBalanceSizeDistribution::
120 functionTypeSymbolicName()
121 {
122  word functionTypeSymbolicName(word::null);
123 
124  switch (functionType_)
125  {
126  case functionType::numberConcentration:
127  {
128  functionTypeSymbolicName = "N";
129 
130  break;
131  }
132  case functionType::numberDensity:
133  {
134  functionTypeSymbolicName = "n";
135 
136  break;
137  }
138  case functionType::volumeConcentration:
139  {
140  functionTypeSymbolicName = "V";
141 
142  break;
143  }
144  case functionType::volumeDensity:
145  {
146  functionTypeSymbolicName = "v";
147 
148  break;
149  }
150  case functionType::areaConcentration:
151  {
152  functionTypeSymbolicName = "A";
153 
154  break;
155  }
156  case functionType::areaDensity:
157  {
158  functionTypeSymbolicName = "a";
159 
160  break;
161  }
162  }
163 
164  return functionTypeSymbolicName;
165 }
166 
167 
169 Foam::functionObjects::populationBalanceSizeDistribution::
170 coordinateTypeSymbolicName
171 (
172  const coordinateType& cType
173 )
174 {
175  word coordinateTypeSymbolicName(word::null);
176 
177  switch (cType)
178  {
179  case coordinateType::volume:
180  {
181  coordinateTypeSymbolicName = "v";
182 
183  break;
184  }
185  case coordinateType::area:
186  {
187  coordinateTypeSymbolicName = "a";
188 
189  break;
190  }
191  case coordinateType::diameter:
192  {
193  coordinateTypeSymbolicName = "d";
194 
195  break;
196  }
197  case coordinateType::projectedAreaDiameter:
198  {
199  coordinateTypeSymbolicName = "dPa";
200 
201  break;
202  }
203  }
204 
205  return coordinateTypeSymbolicName;
206 }
207 
208 
210 Foam::functionObjects::populationBalanceSizeDistribution::filterField
211 (
212  const scalarField& field
213 ) const
214 {
215  if (all())
216  {
217  return field;
218  }
219  else
220  {
221  return tmp<scalarField>(new scalarField(field, cells()));
222  }
223 }
224 
225 
226 Foam::scalar
227 Foam::functionObjects::populationBalanceSizeDistribution::averageCoordinateValue
228 (
230  const coordinateType& cType
231 )
232 {
233  scalar averageCoordinateValue(Zero);
234 
235  switch (cType)
236  {
237  case coordinateType::volume:
238  {
239  averageCoordinateValue = fi.x().value();
240 
241  break;
242  }
243  case coordinateType::area:
244  {
245  averageCoordinateValue =
246  weightedAverage(fi.a(), fi);
247 
248  break;
249  }
250  case coordinateType::diameter:
251  {
252  averageCoordinateValue =
253  weightedAverage(fi.d(), fi);
254 
255  break;
256  }
257  case coordinateType::projectedAreaDiameter:
258  {
259  averageCoordinateValue =
260  weightedAverage(sqrt(fi.a()/pi), fi);
261 
262  break;
263  }
264  }
265 
266  return averageCoordinateValue;
267 }
268 
269 
270 Foam::scalar
271 Foam::functionObjects::populationBalanceSizeDistribution::weightedAverage
272 (
273  const Foam::scalarField& fld,
275 )
276 {
277  scalar weightedAverage(Zero);
278 
279  switch (weightType_)
280  {
281  case weightType::numberConcentration:
282  {
283  scalarField Ni(filterField(fi*fi.phase()/fi.x().value()));
284 
285  if (gSum(Ni) == 0)
286  {
287  weightedAverage =
288  gSum(filterField(mesh_.V()*fld))/this->V();
289  }
290  else
291  {
292  weightedAverage =
293  gSum(Ni*filterField(fld))/gSum(Ni);
294  }
295 
296  break;
297  }
298  case weightType::volumeConcentration:
299  {
300  scalarField Vi(filterField(fi*fi.phase()));
301 
302  if (gSum(Vi) == 0)
303  {
304  weightedAverage =
305  gSum(filterField(mesh_.V()*fld))/this->V();
306  }
307  else
308  {
309  weightedAverage =
310  gSum(Vi*filterField(fld))/gSum(Vi);
311  }
312 
313  break;
314  }
315  case weightType::areaConcentration:
316  {
317  scalarField Ai(filterField(fi.a().ref()*fi.phase()));
318 
319  if (gSum(Ai) == 0)
320  {
321  weightedAverage =
322  gSum(filterField(mesh_.V()*fld))/this->V();
323  }
324  else
325  {
326  weightedAverage =
327  gSum(Ai*filterField(fld))/gSum(Ai);
328  }
329 
330  break;
331  }
332  case weightType::cellVolume:
333  {
334  weightedAverage =
335  gSum(filterField(mesh_.V()*fld))/this->V();
336 
337  break;
338  }
339  }
340 
341  return weightedAverage;
342 }
343 
344 
345 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
346 
349 (
350  const word& name,
351  const Time& runTime,
352  const dictionary& dict
353 )
354 :
355  fvMeshFunctionObject(name, runTime, dict),
357  file_(obr_, name),
358  mesh_(fvMeshFunctionObject::mesh_),
359  popBal_
360  (
361  obr_.lookupObject<Foam::diameterModels::populationBalanceModel>
362  (
363  dict.lookup("populationBalance")
364  )
365  ),
366  functionType_(functionTypeNames_.read(dict.lookup("functionType"))),
367  coordinateType_(coordinateTypeNames_.read(dict.lookup("coordinateType"))),
368  allCoordinates_
369  (
370  dict.lookupOrDefault<Switch>("allCoordinates", false)
371  ),
372  normalise_(dict.lookupOrDefault<Switch>("normalise", false)),
373  logTransform_
374  (
375  dict.lookupOrDefaultBackwardsCompatible<Switch>
376  (
377  {"logTransform", "geometric"},
378  false
379  )
380  ),
381  weightType_
382  (
383  dict.found("weightType")
384  ? weightTypeNames_.read(dict.lookup("weightType"))
385  : weightType::numberConcentration
386  ),
387  formatterPtr_(nullptr)
388 {
389  read(dict);
390 }
391 
392 
393 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
394 
397 {}
398 
399 
400 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
401 
403 (
404  const dictionary& dict
405 )
406 {
407  Log << type() << " " << name() << ":" << nl;
408 
410 
411  formatterPtr_ = setWriter::New(dict.lookup("setFormat"), dict);
412 
413  return false;
414 }
415 
416 
418 {
419  return true;
420 }
421 
422 
424 {
425  Log << type() << " " << name() << " write:" << nl;
426 
427  const UPtrList<diameterModels::sizeGroup>& sizeGroups =
428  popBal_.sizeGroups();
429 
430  scalarField coordinateValues(sizeGroups.size());
431  scalarField boundaryValues(sizeGroups.size() + 1);
432  scalarField resultValues(sizeGroups.size());
433 
434  forAll(sizeGroups, i)
435  {
436  const diameterModels::sizeGroup& fi = sizeGroups[i];
437 
438  coordinateValues[i] = averageCoordinateValue(fi, coordinateType_);
439  }
440 
441  if
442  (
443  functionType_ == functionType::numberDensity
444  || functionType_ == functionType::volumeDensity
445  || functionType_ == functionType::areaDensity
446  )
447  {
448  boundaryValues.first() = coordinateValues.first();
449  boundaryValues.last() = coordinateValues.last();
450 
451  for (label i = 1; i < boundaryValues.size() - 1; i++)
452  {
453  boundaryValues[i] =
454  0.5*(coordinateValues[i] + coordinateValues[i-1]);
455  }
456 
457  if (logTransform_)
458  {
459  boundaryValues = Foam::log(boundaryValues);
460  }
461  }
462 
463  switch (functionType_)
464  {
465  case functionType::numberConcentration:
466  {
467  forAll(sizeGroups, i)
468  {
469  const diameterModels::sizeGroup& fi = sizeGroups[i];
470 
471  resultValues[i] =
472  gSum(filterField(mesh_.V()*fi*fi.phase()/fi.x()))/this->V();
473  }
474 
475  if (normalise_ && sum(resultValues) != 0)
476  {
477  resultValues /= sum(resultValues);
478  }
479 
480  break;
481  }
482  case functionType::numberDensity:
483  {
484  forAll(sizeGroups, i)
485  {
486  const diameterModels::sizeGroup& fi = sizeGroups[i];
487 
488  resultValues[i] =
489  gSum(filterField(mesh_.V()*fi*fi.phase()/fi.x()))/this->V();
490  }
491 
492  if (normalise_ && sum(resultValues) != 0)
493  {
494  resultValues /= sum(resultValues);
495  }
496 
497  forAll(resultValues, i)
498  {
499  resultValues[i] /= (boundaryValues[i+1] - boundaryValues[i]);
500  }
501 
502  break;
503  }
504  case functionType::volumeConcentration:
505  {
506  forAll(sizeGroups, i)
507  {
508  const diameterModels::sizeGroup& fi = sizeGroups[i];
509 
510  resultValues[i] =
511  gSum(filterField(mesh_.V()*fi*fi.phase()))/this->V();
512  }
513 
514  if (normalise_ && sum(resultValues) != 0)
515  {
516  resultValues /= sum(resultValues);
517  }
518 
519  break;
520  }
521  case functionType::volumeDensity:
522  {
523  forAll(sizeGroups, i)
524  {
525  const diameterModels::sizeGroup& fi = sizeGroups[i];
526 
527  resultValues[i] =
528  gSum(filterField(mesh_.V()*fi*fi.phase()))/this->V();
529  }
530 
531  if (normalise_ && sum(resultValues) != 0)
532  {
533  resultValues /= sum(resultValues);
534  }
535 
536  forAll(resultValues, i)
537  {
538  resultValues[i] /= (boundaryValues[i+1] - boundaryValues[i]);
539  }
540 
541  break;
542  }
543  case functionType::areaConcentration:
544  {
545  forAll(sizeGroups, i)
546  {
547  const diameterModels::sizeGroup& fi = sizeGroups[i];
548 
549  resultValues[i] =
550  gSum
551  (
552  filterField(mesh_.V()*fi.a().ref()*fi*fi.phase()/fi.x())
553  )
554  /this->V();
555  }
556 
557  if (normalise_ && sum(resultValues) != 0)
558  {
559  resultValues /= sum(resultValues);
560  }
561 
562  break;
563  }
564  case functionType::areaDensity:
565  {
566  forAll(sizeGroups, i)
567  {
568  const diameterModels::sizeGroup& fi = sizeGroups[i];
569 
570  resultValues[i] =
571  gSum
572  (
573  filterField(mesh_.V()*fi.a().ref()*fi*fi.phase()/fi.x())
574  )
575  /this->V();
576  }
577 
578  if (normalise_ && sum(resultValues) != 0)
579  {
580  resultValues /= sum(resultValues);
581  }
582 
583  forAll(resultValues, i)
584  {
585  resultValues[i] /= (boundaryValues[i+1] - boundaryValues[i]);
586  }
587 
588  break;
589  }
590  }
591 
592 
593  if (allCoordinates_)
594  {
595  wordList otherCoordinateSymbolicNames(coordinateTypeNames_.size());
596  PtrList<scalarField> otherCoordinateValues(coordinateTypeNames_.size());
597  typedef NamedEnum<coordinateType, 4> namedEnumCoordinateType;
598 
599  forAllConstIter(namedEnumCoordinateType, coordinateTypeNames_, iter)
600  {
601  const coordinateType cType = coordinateTypeNames_[iter.key()];
602 
603  otherCoordinateSymbolicNames[cType] =
604  coordinateTypeSymbolicName(cType);
605 
606  otherCoordinateValues.set
607  (
608  cType,
609  new scalarField(popBal_.sizeGroups().size())
610  );
611 
612  forAll(sizeGroups, i)
613  {
614  const diameterModels::sizeGroup& fi = sizeGroups[i];
615 
616  otherCoordinateValues[cType][i] =
617  averageCoordinateValue(fi, cType);
618  }
619  }
620 
621  if (Pstream::master())
622  {
623  formatterPtr_->write
624  (
625  file_.baseTimeDir(),
626  name(),
627  coordSet
628  (
629  true,
630  coordinateTypeSymbolicName(coordinateType_),
631  coordinateValues
632  ),
633  functionTypeSymbolicName(),
634  resultValues,
635  otherCoordinateSymbolicNames,
636  otherCoordinateValues
637  );
638  }
639  }
640  else
641  {
642  if (Pstream::master())
643  {
644  formatterPtr_->write
645  (
646  file_.baseTimeDir(),
647  name(),
648  coordSet
649  (
650  true,
651  coordinateTypeSymbolicName(coordinateType_),
652  coordinateValues
653  ),
654  functionTypeSymbolicName(),
655  resultValues
656  );
657  }
658  }
659 
660  return true;
661 }
662 
663 
664 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Macros for easy insertion into run-time selection tables.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:54
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
bool set(const label) const
Is element set.
Definition: PtrListI.H:65
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
T & first()
Return the first element of the list.
Definition: UListI.H:114
T & last()
Return the last element of the list.
Definition: UListI.H:128
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: UPtrList.H:66
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
Holds list of sampling positions.
Definition: coordSet.H:51
Single size class fraction field representing a fixed particle volume as defined by the user through ...
Definition: sizeGroup.H:102
const dimensionedScalar & x() const
Return representative volume of the sizeGroup.
Definition: sizeGroupI.H:57
const tmp< volScalarField > a() const
Return representative surface area of the sizeGroup.
Definition: sizeGroup.C:135
const phaseModel & phase() const
Return const-reference to the phase.
Definition: sizeGroupI.H:36
const tmp< volScalarField > d() const
Return representative diameter of the sizeGroup.
Definition: sizeGroup.C:142
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
const Type & value() const
Return const reference to value.
Abstract base-class for Time/database functionObjects.
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Writes out the size distribution determined by a population balance model, either for the entire doma...
static const NamedEnum< coordinateType, 4 > coordinateTypeNames_
Coordinate type names.
populationBalanceSizeDistribution(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
static const NamedEnum< weightType, 4 > weightTypeNames_
Names of the weight types.
static const NamedEnum< functionType, 6 > functionTypeNames_
Function type names.
virtual bool write()
Calculate and write the size distribution.
virtual bool read(const dictionary &)
Read the populationBalanceSizeDistribution data.
virtual bool read(const dictionary &)
Read optional controls.
General run-time selected cell set selection class for fvMesh.
Definition: fvCellSet.H:88
static autoPtr< setWriter > New(const word &writeType, const IOstream::streamFormat writeFormat=IOstream::ASCII, const IOstream::compressionType writeCompression=IOstream::UNCOMPRESSED)
Select given write options.
Definition: setWriter.C:286
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
static const word null
An empty word.
Definition: word.H:77
const cellShapeList & cells
gmvFile<< "tracers "<< particles.size()<< nl;{ pointField positions(particles.size());label particlei=0;forAllConstIter(Cloud< passiveParticle >, particles, iter) { positions[particlei++]=iter().position(mesh);} for(i=0;i< pTraits< point >::nComponents;i++) { forAll(positions, particlei) { gmvFile<< component(positions[particlei], i)<< ' ';} gmvFile<< nl;}}forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
#define Log
Report write to Foam::Info if the local log switch is true.
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
Type gSum(const FieldField< Field, Type > &f)
bool read(const char *, int32_t &)
Definition: int32IO.C:85
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
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar sqrt(const dimensionedScalar &ds)
static const char nl
Definition: Ostream.H:260
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