sizeDistribution.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017-2019 OpenFOAM Foundation
9 Copyright (C) 2019-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "sizeDistribution.H"
30#include "sizeGroup.H"
32
33// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34
35namespace Foam
36{
37namespace functionObjects
38{
41}
42}
43
44
45const Foam::Enum
46<
48>
50({
51 {selectionModeTypes::rtCellZone, "cellZone"},
52 {selectionModeTypes::rtAll, "all"},
53});
54
55
56const Foam::Enum
57<
59>
61({
62 {functionTypes::ftNdf, "numberDensity"},
63 {functionTypes::ftVdf, "volumeDensity"},
64 {functionTypes::ftNc, "numberConcentration"},
65 {functionTypes::ftMom, "moments"},
66});
67
68
69const Foam::Enum
70<
72>
74({
75
76 {abszissaTypes::atDiameter, "diameter"},
77 {abszissaTypes::atVolume, "volume"},
78});
79
80
81// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
82
84(
85 const dictionary& dict
86)
87{
88 switch (functionType_)
89 {
90 case ftNdf:
91 {
92 break;
93 }
94
95 case ftVdf:
96 {
97 break;
98 }
99
100 case ftNc:
101 {
102 break;
103 }
104
105 case ftMom:
106 {
107 break;
108 }
109
110 default:
111 {
113 (
114 dict,
115 "functionType",
118 ) << exit(FatalIOError);
119 }
120 }
121
122 switch (abszissaType_)
123 {
124 case atDiameter:
125 {
126 break;
127 }
128
129 case atVolume:
130 {
131 break;
132 }
133
134 default:
135 {
137 (
138 dict,
139 "abszissaType",
142 ) << exit(FatalIOError);
143 }
144 }
145
147
148 if (nCells_ == 0)
149 {
151 << type() << " " << name() << ": "
153 << '(' << selectionModeTypeName_ << "):" << nl
154 << " Selection has no cells" << exit(FatalIOError);
155 }
156
157 volume_ = volume();
158
159 Info<< type() << " " << name() << ":"
161 << '(' << selectionModeTypeName_ << "):" << nl
162 << " total cells = " << nCells_ << nl
163 << " total volume = " << volume_
164 << nl << endl;
165}
166
167
169{
170 switch (selectionModeType_)
171 {
172 case rtCellZone:
173 {
174 dict().readEntry("cellZone", selectionModeTypeName_);
175
176 label zoneId =
177 mesh().cellZones().findZoneID(selectionModeTypeName_);
178
179 if (zoneId < 0)
180 {
182 << "Unknown cellZone name: " << selectionModeTypeName_
183 << ". Valid cellZone names are: "
184 << mesh().cellZones().names()
185 << nl << exit(FatalIOError);
186 }
187
188 cellId_ = mesh().cellZones()[zoneId];
189 nCells_ = returnReduce(cellId_.size(), sumOp<label>());
190 break;
191 }
192
193 case rtAll:
194 {
195 cellId_ = identity(mesh().nCells());
196 nCells_ = returnReduce(cellId_.size(), sumOp<label>());
197 break;
198 }
199
200 default:
201 {
203 (
204 dict_,
205 "selectionMode",
207 selectionModeTypeNames_
208 ) << exit(FatalIOError);
209 }
210 }
211}
212
213
215{
216 return gSum(filterField(mesh().V()));
217}
218
219
221{
223
224 allValues[Pstream::myProcNo()] = field;
225
226 Pstream::gatherList(allValues);
227
228 if (Pstream::master())
229 {
230 field =
231 ListListOps::combine<scalarField>
232 (
233 allValues,
235 );
236 }
237}
238
239
242(
243 const scalarField& field
244) const
245{
246 return tmp<scalarField>(new scalarField(field, cellId_));
247}
248
249
251(
252 const label i
253)
254{
255 OFstream& file = this->file();
256
257 switch (functionType_)
258 {
259 case ftNdf:
260 {
261 writeHeader(file, "Number density function");
262 break;
263 }
264
265 case ftVdf:
266 {
267 writeHeader(file, "Volume density function");
268 break;
269 }
270
271 case ftNc:
272 {
273 writeHeader(file, "Number concentration");
274 break;
275 }
276
277 case ftMom:
278 {
279 writeHeader(file, "Moments");
280 break;
281 }
282 }
283
284 switch (abszissaType_)
285 {
286 case atVolume:
287 {
288 writeCommented(file, "Time/volume");
289 break;
290 }
291
292 case atDiameter:
293 {
294 writeCommented(file, "Time/diameter");
295 break;
296 }
297 }
298
299 switch (functionType_)
300 {
301 case ftMom:
302 {
303 for (label i = 0; i <= momentOrder_; i++)
304 {
305 file() << tab << i;
306 }
307
308 break;
309 }
310 default:
311 {
312 forAll(popBal_.sizeGroups(), sizeGroupi)
313 {
314 const diameterModels::sizeGroup& fi =
315 popBal_.sizeGroups()[sizeGroupi];
316
317 switch (abszissaType_)
318 {
319 case atDiameter:
320 {
321 file() << tab << fi.d().value();
322
323 break;
324 }
325
326 case atVolume:
327 {
328 file() << tab << fi.x().value();
329
330 break;
331 }
332 }
333 }
334
335 break;
336 }
337 }
338
339 file << endl;
340}
341
342
343// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
344
346(
347 const word& name,
348 const Time& runTime,
349 const dictionary& dict
350)
351:
353 writeFile(obr_, name),
354 dict_(dict),
355 selectionModeType_
356 (
357 selectionModeTypeNames_.get("selectionMode", dict)
358 ),
359 selectionModeTypeName_(word::null),
360 functionType_(functionTypeNames_.get("functionType", dict)),
361 abszissaType_(abszissaTypeNames_.get("abszissaType", dict)),
362 nCells_(0),
363 cellId_(),
364 volume_(0.0),
365 writeVolume_(dict.getOrDefault("writeVolume", false)),
366 popBal_
367 (
368 obr_.lookupObject<Foam::diameterModels::populationBalanceModel>
369 (
370 dict.get<word>("populationBalance")
371 )
372 ),
373 N_(popBal_.sizeGroups().size()),
374 momentOrder_(dict.getOrDefault<label>("momentOrder", 0)),
375 normalize_(dict.getOrDefault("normalize", false)),
376 sumN_(0.0),
377 sumV_(0.0)
378{
379 read(dict);
382}
383
384
385// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
386
388{}
389
390
391// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
392
394{
395 if (dict != dict_)
396 {
397 dict_ = dict;
398 }
399
402
403 initialise(dict);
404
405 return true;
406}
407
408
410{
411 return true;
412}
413
414
416{
417 writeFileHeader();
418 writeCurrentTime(file());
419
420 Log << type() << " " << name() << " write" << nl;
421
422 scalarField V(filterField(mesh().V()));
423 combineFields(V);
424
425 sumN_ = 0;
426 sumV_ = 0;
427
428 forAll(N_, i)
429 {
430 const Foam::diameterModels::sizeGroup& fi = popBal_.sizeGroups()[i];
431
432 const volScalarField& alpha = fi.VelocityGroup().phase();
433
434 scalarField Ni(fi*alpha/fi.x());
435 scalarField values(filterField(Ni));
436 scalarField V(filterField(mesh().V()));
437
438 // Combine onto master
439 combineFields(values);
440 combineFields(V);
441
442 if (Pstream::master())
443 {
444 // Calculate volume-averaged number concentration
445 N_[i] = sum(V*values)/sum(V);
446 }
447
448 sumN_ += N_[i];
449
450 sumV_ += N_[i]*fi.x().value();
451 }
452
453 if (Pstream::master())
454 {
455 switch (functionType_)
456 {
457 case ftMom:
458 {
459 for (label m = 0; m <= momentOrder_; m++)
460 {
461 scalar result(0.0);
462
463 forAll(N_, i)
464 {
466 popBal_.sizeGroups()[i];
467
468 switch (abszissaType_)
469 {
470 case atVolume:
471 {
472 result += pow(fi.x().value(), m)*N_[i];
473
474 break;
475 }
476
477 case atDiameter:
478 {
479 result += pow(fi.d().value(), m)*N_[i];
480
481 break;
482 }
483 }
484 }
485
486 file() << tab << result;
487 }
488
489 break;
490 }
491
492 default:
493 {
494 forAll(popBal_.sizeGroups(), i)
495 {
497 popBal_.sizeGroups()[i];
498
499 scalar result(0.0);
500 scalar delta(0.0);
501
502 switch (abszissaType_)
503 {
504 case atVolume:
505 {
506 delta = popBal_.v()[i+1].value()
507 - popBal_.v()[i].value();
508
509 break;
510 }
511
512 case atDiameter:
513 {
514 const scalar& formFactor =
516
517 delta =
518 pow
519 (
520 popBal_.v()[i+1].value()
521 /formFactor,
522 1.0/3.0
523 )
524 - pow
525 (
526 popBal_.v()[i].value()
527 /formFactor,
528 1.0/3.0
529 );
530
531 break;
532 }
533 }
534
535 switch (functionType_)
536 {
537 case ftNdf:
538 {
539 if (normalize_ == true)
540 {
541 result = N_[i]/delta/sumN_;
542 }
543 else
544 {
545 result = N_[i]/delta;
546 }
547
548 break;
549 }
550
551 case ftVdf:
552 {
553 if (normalize_ == true)
554 {
555 result = N_[i]*fi.x().value()/delta/sumV_;
556 }
557 else
558 {
559 result = N_[i]*fi.x().value()/delta;
560 }
561
562 break;
563 }
564
565 case ftNc:
566 {
567 if (normalize_ == true)
568 {
569 result = N_[i]/sumN_;
570 }
571 else
572 {
573 result = N_[i];
574 }
575
576 break;
577 }
578
579 default:
580 {
581 break;
582 }
583 }
584
585 file()<< tab << result;
586 }
587 }
588 }
589 }
590 {
591 file()<< endl;
592 }
593
594 Log << endl;
595
596 return true;
597}
598
599
600// ************************************************************************* //
scalar delta
#define Log
Definition: PDRblock.C:35
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
Output to file stream, using an OSstream.
Definition: OFstream.H:57
label nProcs() const noexcept
Number of ranks associated with PstreamBuffers.
static void gatherList(const List< commsStruct > &comms, List< T > &values, const int tag, const label comm)
virtual bool read()
Re-read model coefficients if they have changed.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
const phaseModel & phase() const
Return the phase.
This class represents a single sizeGroup belonging to a velocityGroup. The main property of a sizeGro...
Definition: sizeGroup.H:99
const dimensionedScalar & d() const
Return representative diameter of the sizeGroup.
Definition: sizeGroupI.H:52
const dimensionedScalar & x() const
Return representative volume of the sizeGroup.
Definition: sizeGroupI.H:59
const velocityGroup & VelocityGroup() const
Return const-reference to the velocityGroup.
Definition: sizeGroupI.H:45
const dimensionedScalar & formFactor() const
Return the form factor.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
const Type & value() const
Return const reference to value.
Abstract base-class for Time/database function objects.
const word & name() const noexcept
Return the name of this functionObject.
virtual const word & type() const =0
Runtime type information.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Computes the power of an input volScalarField.
Definition: pow.H:217
This function object calculates and outputs information about the size distribution of the dispersed ...
label nCells_
Global number of cells.
word selectionModeTypeName_
Name of selection.
selectionModeTypes
Selection mode type enumeration.
selectionModeTypes selectionModeType_
Selection mode type.
void setCellZoneCells()
Set cells to evaluate based on a cell zone.
scalar volume() const
Calculate and return volume of the evaluated cell zone.
abszissaTypes abszissaType_
Abszissa type.
virtual bool read(const dictionary &dict)
Read from dictionary.
functionTypes functionType_
Function type.
static const Enum< selectionModeTypes > selectionModeTypeNames_
Selection mode type names.
const dictionary & dict() const
Return the reference to the construction dictionary.
void writeFileHeader(const label i=0)
Output file header information.
static const Enum< abszissaTypes > abszissaTypeNames_
Abszissa type names.
static const Enum< functionTypes > functionTypeNames_
Function type names.
void combineFields(scalarField &field)
Combine fields from all processor domains into single field.
tmp< scalarField > filterField(const scalarField &field) const
Filter field according to cellIds.
scalar volume_
Total volume of the evaluated selection.
functionTypes
Function type enumeration.
Base class for writing single files from the function objects.
Definition: writeFile.H:120
virtual autoPtr< OFstream > createFile(const word &name, scalar timeValue) const
Return autoPtr to a new file for a given time.
Definition: writeFile.C:80
virtual void resetFile(const word &name)
Reset internal file pointer to new file with new name.
Definition: writeFile.C:134
int myProcNo() const noexcept
Return processor number.
splitCell * master() const
Definition: splitCell.H:113
A class for managing temporary objects.
Definition: tmp.H:65
void initialise()
Initialise integral-scale box properties.
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
rDeltaTY field()
dynamicFvMesh & mesh
engineTime & runTime
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
Definition: error.H:478
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
Namespace for OpenFOAM.
Type gSum(const FieldField< Field, Type > &f)
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
Definition: labelList.C:38
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
static void writeHeader(Ostream &os, const word &fieldName)
messageStream Info
Information stream (stdout output on master, null elsewhere)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
IOerror FatalIOError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
constexpr char tab
The tab '\t' character(0x09)
Definition: Ostream.H:52
volScalarField & alpha
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333