populationBalanceMoments.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) 2022-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 
46 namespace Foam
47 {
48  template<>
49  const char* NamedEnum
50  <
52  4
53  >::names[] = {"integerMoment", "mean", "variance", "stdDev"};
54 }
55 
56 
57 const Foam::NamedEnum
58 <
60  4
61 >
63 
64 
65 namespace Foam
66 {
67  template<>
68  const char* NamedEnum
69  <
71  3
72  >::names[] = {"volume", "area", "diameter"};
73 }
74 
75 
76 const Foam::NamedEnum
77 <
79  3
80 >
82 
83 
84 namespace Foam
85 {
86  template<>
87  const char* NamedEnum
88  <
90  3
91  >::names[] =
92  {
93  "numberConcentration",
94  "volumeConcentration",
95  "areaConcentration"
96  };
97 }
98 
99 
100 const Foam::NamedEnum
101 <
103  3
104 >
106 
107 
108 namespace Foam
109 {
110  template<>
111  const char* NamedEnum
112  <
114  3
115  >::names[] = {"arithmetic", "geometric", "notApplicable"};
116 }
117 
118 
119 const Foam::NamedEnum
120 <
122  3
123 >
125 
126 
127 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
128 
130 Foam::functionObjects::populationBalanceMoments::coordinateTypeSymbolicName()
131 {
132  word coordinateTypeSymbolicName(word::null);
133 
134  switch (coordinateType_)
135  {
137  {
138  coordinateTypeSymbolicName = "v";
139 
140  break;
141  }
143  {
144  coordinateTypeSymbolicName = "a";
145 
146  break;
147  }
149  {
150  coordinateTypeSymbolicName = "d";
151 
152  break;
153  }
154  }
155 
156  return coordinateTypeSymbolicName;
157 }
158 
159 
161 Foam::functionObjects::populationBalanceMoments::weightTypeSymbolicName()
162 {
163  word weightTypeSymbolicName(word::null);
164 
165  switch (weightType_)
166  {
167  case weightType::numberConcentration:
168  {
169  weightTypeSymbolicName = "N";
170 
171  break;
172  }
173  case weightType::volumeConcentration:
174  {
175  weightTypeSymbolicName = "V";
176 
177  break;
178  }
179  case weightType::areaConcentration:
180  {
181  weightTypeSymbolicName = "A";
182 
183  break;
184  }
185  }
186 
187  return weightTypeSymbolicName;
188 }
189 
190 
191 Foam::word Foam::functionObjects::populationBalanceMoments::defaultFldName()
192 {
193  word meanName
194  (
195  meanType_ == meanType::geometric
196  ? word(meanTypeNames_[meanType_]).capitalise()
197  : word("")
198  );
199 
200  return
201  word
202  (
204  (
205  "weighted"
206  + meanName
207  + word(momentTypeNames_[momentType_]).capitalise()
208  + "("
209  + weightTypeSymbolicName()
210  + ","
211  + coordinateTypeSymbolicName()
212  + ")",
213  popBal_.name()
214  )
215  );
216 }
217 
218 
220 Foam::functionObjects::populationBalanceMoments::integerMomentFldName()
221 {
222  return
223  word
224  (
226  (
227  word(momentTypeNames_[momentType_])
228  + Foam::name(order_)
229  + "("
230  + weightTypeSymbolicName()
231  + ","
232  + coordinateTypeSymbolicName()
233  + ")",
234  popBal_.name()
235  )
236  );
237 }
238 
239 
240 void Foam::functionObjects::populationBalanceMoments::setDimensions
241 (
243  momentType momType
244 )
245 {
246  switch (momType)
247  {
248  case momentType::integerMoment:
249  {
250  switch (coordinateType_)
251  {
252  case coordinateType::volume:
253  {
254  fld.dimensions().reset
255  (
256  pow(dimVolume, order_)/dimVolume
257  );
258 
259  break;
260  }
261  case coordinateType::area:
262  {
263  fld.dimensions().reset
264  (
265  pow(dimArea, order_)/dimVolume
266  );
267 
268  break;
269  }
270  case coordinateType::diameter:
271  {
272  fld.dimensions().reset
273  (
274  pow(dimLength, order_)/dimVolume
275  );
276 
277  break;
278  }
279  }
280 
281  switch (weightType_)
282  {
283  case weightType::volumeConcentration:
284  {
285  fld.dimensions().reset(fld.dimensions()*dimVolume);
286 
287  break;
288  }
289  case weightType::areaConcentration:
290  {
291  fld.dimensions().reset(fld.dimensions()*dimArea);
292 
293  break;
294  }
295  default:
296  {
297  break;
298  }
299  }
300 
301  break;
302  }
303  case momentType::mean:
304  {
305  switch (coordinateType_)
306  {
307  case coordinateType::volume:
308  {
309  fld.dimensions().reset(dimVolume);
310 
311  break;
312  }
313  case coordinateType::area:
314  {
315  fld.dimensions().reset(dimArea);
316 
317  break;
318  }
319  case coordinateType::diameter:
320  {
321  fld.dimensions().reset(dimLength);
322 
323  break;
324  }
325  }
326 
327  break;
328  }
329  case momentType::variance:
330  {
331  switch (coordinateType_)
332  {
333  case coordinateType::volume:
334  {
335  fld.dimensions().reset(sqr(dimVolume));
336 
337  break;
338  }
339  case coordinateType::area:
340  {
341  fld.dimensions().reset(sqr(dimArea));
342 
343  break;
344  }
345  case coordinateType::diameter:
346  {
347  fld.dimensions().reset(sqr(dimLength));
348 
349  break;
350  }
351  }
352 
353  if (meanType_ == meanType::geometric)
354  {
355  fld.dimensions().reset(dimless);
356  }
357 
358  break;
359  }
360  case momentType::stdDev:
361  {
362  switch (coordinateType_)
363  {
364  case coordinateType::volume:
365  {
366  fld.dimensions().reset(dimVolume);
367 
368  break;
369  }
370  case coordinateType::area:
371  {
372  fld.dimensions().reset(dimArea);
373 
374  break;
375  }
376  case coordinateType::diameter:
377  {
378  fld.dimensions().reset(dimLength);
379 
380  break;
381  }
382  }
383 
384  if (meanType_ == meanType::geometric)
385  {
386  fld.dimensions().reset(dimless);
387  }
388 
389  break;
390  }
391  }
392 }
393 
394 
396 Foam::functionObjects::populationBalanceMoments::totalConcentration()
397 {
398  tmp<volScalarField> tTotalConcentration
399  (
401  (
402  "totalConcentration",
403  mesh_,
405  )
406  );
407 
408  volScalarField& totalConcentration = tTotalConcentration.ref();
409 
410  switch (weightType_)
411  {
412  case weightType::volumeConcentration:
413  {
414  totalConcentration.dimensions().reset
415  (
416  totalConcentration.dimensions()*dimVolume
417  );
418 
419  break;
420  }
421  case weightType::areaConcentration:
422  {
423  totalConcentration.dimensions().reset
424  (
425  totalConcentration.dimensions()*dimArea
426  );
427 
428  break;
429  }
430  default:
431  {
432  break;
433  }
434  }
435 
436  forAll(popBal_.sizeGroups(), i)
437  {
438  const Foam::diameterModels::sizeGroup& fi = popBal_.sizeGroups()[i];
439 
440  switch (weightType_)
441  {
442  case weightType::numberConcentration:
443  {
444  totalConcentration += fi*fi.phase()/fi.x();
445 
446  break;
447  }
448  case weightType::volumeConcentration:
449  {
450  totalConcentration += fi*fi.phase();
451 
452  break;
453  }
454  case weightType::areaConcentration:
455  {
456  totalConcentration += fi.a()*fi*fi.phase()/fi.x();
457 
458  break;
459  }
460  }
461  }
462 
463  return tTotalConcentration;
464 }
465 
466 
468 Foam::functionObjects::populationBalanceMoments::mean()
469 {
470  tmp<volScalarField> tMean
471  (
473  (
474  "mean",
475  mesh_,
477  )
478  );
479 
480  volScalarField& mean = tMean.ref();
481 
482  setDimensions(mean, momentType::mean);
483 
484  volScalarField totalConcentration(this->totalConcentration());
485 
486  forAll(popBal_.sizeGroups(), i)
487  {
488  const Foam::diameterModels::sizeGroup& fi = popBal_.sizeGroups()[i];
489 
490  volScalarField concentration(fi*fi.phase()/fi.x());
491 
492  switch (weightType_)
493  {
494  case weightType::volumeConcentration:
495  {
496  concentration *= fi.x();
497 
498  break;
499  }
500  case weightType::areaConcentration:
501  {
502  concentration *= fi.a();
503 
504  break;
505  }
506  default:
507  {
508  break;
509  }
510  }
511 
512  switch (meanType_)
513  {
514  case meanType::geometric:
515  {
516  mean.dimensions().reset(dimless);
517 
518  switch (coordinateType_)
519  {
520  case coordinateType::volume:
521  {
522  dimensionedScalar unitVolume(dimVolume, 1);
523 
524  mean +=
525  Foam::log(fi.x()/unitVolume)
526  *concentration/totalConcentration;
527 
528  break;
529  }
530  case coordinateType::area:
531  {
532  dimensionedScalar unitArea(dimArea, 1);
533 
534  mean +=
535  Foam::log(fi.a()/unitArea)
536  *concentration/totalConcentration;
537 
538  break;
539  }
540  case coordinateType::diameter:
541  {
542  dimensionedScalar unitLength(dimLength, 1);
543 
544  mean +=
545  Foam::log(fi.d()/unitLength)
546  *concentration/totalConcentration;
547 
548  break;
549  }
550  }
551 
552  break;
553  }
554  default:
555  {
556  switch (coordinateType_)
557  {
558  case coordinateType::volume:
559  {
560  mean += fi.x()*concentration/totalConcentration;
561 
562  break;
563  }
564  case coordinateType::area:
565  {
566  mean += fi.a()*concentration/totalConcentration;
567 
568  break;
569  }
570  case coordinateType::diameter:
571  {
572  mean += fi.d()*concentration/totalConcentration;
573 
574  break;
575  }
576  }
577 
578  break;
579  }
580  }
581  }
582 
583  if (meanType_ == meanType::geometric)
584  {
585  mean = exp(mean);
586 
587  setDimensions(mean, momentType::mean);
588  }
589 
590  return tMean;
591 }
592 
593 
595 Foam::functionObjects::populationBalanceMoments::variance()
596 {
597  tmp<volScalarField> tVariance
598  (
600  (
601  "variance",
602  mesh_,
604  )
605  );
606 
607  volScalarField& variance = tVariance.ref();
608 
609  setDimensions(variance, momentType::variance);
610 
611  volScalarField totalConcentration(this->totalConcentration());
612  volScalarField mean(this->mean());
613 
614  forAll(popBal_.sizeGroups(), i)
615  {
616  const Foam::diameterModels::sizeGroup& fi = popBal_.sizeGroups()[i];
617 
618  volScalarField concentration(fi*fi.phase()/fi.x());
619 
620  switch (weightType_)
621  {
622  case weightType::volumeConcentration:
623  {
624  concentration *= fi.x();
625 
626  break;
627  }
628  case weightType::areaConcentration:
629  {
630  concentration *= fi.a();
631 
632  break;
633  }
634  default:
635  {
636  break;
637  }
638  }
639 
640  switch (meanType_)
641  {
642  case meanType::geometric:
643  {
644  switch (coordinateType_)
645  {
646  case coordinateType::volume:
647  {
648  variance +=
649  sqr(Foam::log(fi.x()/mean))
650  *concentration/totalConcentration;
651 
652  break;
653  }
654  case coordinateType::area:
655  {
656  variance +=
657  sqr(Foam::log(fi.a()/mean))
658  *concentration/totalConcentration;
659 
660  break;
661  }
662  case coordinateType::diameter:
663  {
664  variance +=
665  sqr(Foam::log(fi.d()/mean))
666  *concentration/totalConcentration;
667 
668  break;
669  }
670  }
671 
672  break;
673  }
674  default:
675  {
676  switch (coordinateType_)
677  {
678  case coordinateType::volume:
679  {
680  variance +=
681  sqr(fi.x() - mean)*concentration/totalConcentration;
682 
683  break;
684  }
685  case coordinateType::area:
686  {
687  variance +=
688  sqr(fi.a() - mean)*concentration/totalConcentration;
689 
690  break;
691  }
692  case coordinateType::diameter:
693  {
694  variance +=
695  sqr(fi.d() - mean)*concentration/totalConcentration;
696 
697  break;
698  }
699  }
700 
701  break;
702  }
703  }
704  }
705 
706  return tVariance;
707 }
708 
709 
711 Foam::functionObjects::populationBalanceMoments::stdDev()
712 {
713  switch (meanType_)
714  {
715  case meanType::geometric:
716  {
717  return exp(sqrt(this->variance()));
718  }
719  default:
720  {
721  return sqrt(this->variance());
722  }
723  }
724 }
725 
726 
727 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
728 
730 (
731  const word& name,
732  const Time& runTime,
733  const dictionary& dict
734 )
735 :
736  fvMeshFunctionObject(name, runTime, dict),
737  popBal_
738  (
739  obr_.lookupObject<Foam::diameterModels::populationBalanceModel>
740  (
741  dict.lookup("populationBalance")
742  )
743  ),
744  momentType_(momentTypeNames_.read(dict.lookup("momentType"))),
745  coordinateType_(coordinateTypeNames_.read(dict.lookup("coordinateType"))),
746  weightType_
747  (
748  dict.found("weightType")
749  ? weightTypeNames_.read(dict.lookup("weightType"))
750  : weightType::numberConcentration
751  ),
752  meanType_(meanType::notApplicable),
753  order_(-1),
754  fldPtr_(nullptr)
755 {
756  read(dict);
757 }
758 
759 
760 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
761 
763 {}
764 
765 
766 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
767 
768 bool
770 {
772 
773  switch (momentType_)
774  {
775  case momentType::integerMoment:
776  {
777  order_ = dict.lookup<int>("order");
778 
779  break;
780  }
781  default:
782  {
783  meanType_ =
784  dict.found("meanType")
785  ? meanTypeNames_.read(dict.lookup("meanType"))
786  : meanType::arithmetic;
787 
788  break;
789  }
790  }
791 
792  switch (momentType_)
793  {
794  case momentType::integerMoment:
795  {
796  fldPtr_.set
797  (
798  new volScalarField
799  (
800  IOobject
801  (
802  this->integerMomentFldName(),
803  mesh_.time().name(),
804  mesh_,
807  ),
808  mesh_,
810  )
811  );
812 
813  volScalarField& integerMoment = fldPtr_();
814 
815  setDimensions(integerMoment, momentType::integerMoment);
816 
817  break;
818  }
819  case momentType::mean:
820  {
821  fldPtr_.set
822  (
823  new volScalarField
824  (
825  IOobject
826  (
827  this->defaultFldName(),
828  mesh_.time().name(),
829  mesh_,
832  ),
833  this->mean()
834  )
835  );
836 
837  break;
838  }
839  case momentType::variance:
840  {
841  fldPtr_.set
842  (
843  new volScalarField
844  (
845  IOobject
846  (
847  this->defaultFldName(),
848  mesh_.time().name(),
849  mesh_,
852  ),
853  this->variance()
854  )
855  );
856 
857  break;
858  }
859  case momentType::stdDev:
860  {
861  fldPtr_.set
862  (
863  new volScalarField
864  (
865  IOobject
866  (
867  this->defaultFldName(),
868  mesh_.time().name(),
869  mesh_,
872  ),
873  this->stdDev()
874  )
875  );
876 
877  break;
878  }
879  }
880 
881  return true;
882 }
883 
884 
886 {
887  switch (momentType_)
888  {
889  case momentType::integerMoment:
890  {
891  volScalarField& integerMoment = fldPtr_();
892 
893  integerMoment = Zero;
894 
895  forAll(popBal_.sizeGroups(), i)
896  {
898  popBal_.sizeGroups()[i];
899 
900  volScalarField concentration(fi*fi.phase()/fi.x());
901 
902  switch (weightType_)
903  {
904  case weightType::volumeConcentration:
905  {
906  concentration *= fi.x();
907 
908  break;
909  }
910  case weightType::areaConcentration:
911  {
912  concentration *= fi.a();
913 
914  break;
915  }
916  default:
917  {
918  break;
919  }
920  }
921 
922  switch (coordinateType_)
923  {
924  case coordinateType::volume:
925  {
926  integerMoment +=
927  pow(fi.x(), order_)*concentration;
928 
929  break;
930  }
931  case coordinateType::area:
932  {
933  integerMoment +=
934  pow(fi.a(), order_)*concentration;
935 
936  break;
937  }
938  case coordinateType::diameter:
939  {
940  integerMoment +=
941  pow(fi.d(), order_)*concentration;
942 
943  break;
944  }
945  }
946  }
947 
948  break;
949  }
950  case momentType::mean:
951  {
952  fldPtr_() = this->mean();
953 
954  break;
955  }
956  case momentType::variance:
957  {
958  fldPtr_() = this->variance();
959 
960  break;
961  }
962  case momentType::stdDev:
963  {
964  fldPtr_() = sqrt(this->variance());
965 
966  break;
967  }
968  }
969 
970  return true;
971 }
972 
973 
975 {
976  writeObject(fldPtr_->name());
977 
978  return true;
979 }
980 
981 
982 // ************************************************************************* //
bool found
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
static word groupName(Name name, const word &group)
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:54
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
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:58
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:37
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:162
Abstract base-class for Time/database functionObjects.
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Calculates and writes out integral (integer moments) or mean properties (mean, variance,...
static const NamedEnum< coordinateType, 3 > coordinateTypeNames_
Names of the coordinate types.
static const NamedEnum< meanType, 3 > meanTypeNames_
Names of the mean types.
static const NamedEnum< momentType, 4 > momentTypeNames_
Names of the moment types.
static const NamedEnum< weightType, 3 > weightTypeNames_
Names of the weight types.
coordinateType
Enumeration for the coordinate types.
populationBalanceMoments(const word &name, const Time &runTime, const dictionary &)
Construct from Time and dictionary.
virtual bool execute()
Calculate the moment fields.
virtual bool write()
Write the moment fields.
virtual bool read(const dictionary &)
Read the data.
virtual bool read(const dictionary &)
Read optional controls.
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
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))
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
dimensionedScalar exp(const dimensionedScalar &ds)
bool read(const char *, int32_t &)
Definition: int32IO.C:85
dimensionedSymmTensor sqr(const dimensionedVector &dv)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
const dimensionSet dimLength
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimVolume
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
const dimensionSet dimArea
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dictionary dict