RosinRammler.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) 2011-2020 OpenFOAM Foundation
9 Copyright (C) 2021 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 "RosinRammler.H"
30#include "MathFunctions.H"
32
33// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34
35namespace Foam
36{
37namespace distributionModels
38{
41}
42}
43
44// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45
47(
48 const dictionary& dict,
50)
51:
52 distributionModel(typeName, dict, rndGen),
53 lambda_(distributionModelDict_.getCompat<scalar>("lambda", {{"d", 2106}})),
54 n_(distributionModelDict_.get<scalar>("n"))
55{
56 const word parcelBasisType =
57 dict.getOrDefault<word>("parcelBasisType", "none");
58
59 if (parcelBasisType == "mass")
60 {
62 << "Selected parcel basis type: " << parcelBasisType << nl
63 << " Please consider to use massRosinRammler distribution model"
64 << endl;
65 }
66
67 if (lambda_ < VSMALL || n_ < VSMALL)
68 {
70 << "Scale/Shape parameter cannot be equal to or less than zero:"
71 << " lambda = " << lambda_
72 << " n = " << n_
73 << exit(FatalError);
74 }
75
76 check();
77}
78
79
81:
83 lambda_(p.lambda_),
84 n_(p.n_)
85{}
86
87
88// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
89
91{
92 const scalar u = rndGen_.sample01<scalar>();
93 const scalar qMin = pow(minValue_/lambda_, n_);
94 const scalar qMax = pow(maxValue_/lambda_, n_);
95 const scalar r = scalar(1) - exp(-qMax + qMin);
96 return lambda_*pow(qMin - log(scalar(1) - u*r), scalar(1)/n_);
97}
98
99
101{
102 // (C:Eq. 5)
103 const scalar a = scalar(1)/lambda_ + scalar(1);
104 const scalar qMax = pow(maxValue_/n_, lambda_);
105 const scalar qMin = pow(minValue_/n_, lambda_);
106 const scalar gMax = Math::incGamma_P(a, qMax);
107 const scalar gMin = Math::incGamma_P(a, qMin);
108
109 return n_/(exp(-qMin) - exp(-qMax))*(gMax - gMin);
110}
111
112
113// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Random number generator.
Definition: Random.H:60
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A library of runtime-selectable doubly-truncated probability distribution models. Returns random samp...
Particle-size distribution model wherein random samples are drawn from the doubly-truncated two-param...
Definition: RosinRammler.H:222
virtual scalar meanValue() const
Return the theoretical mean of the distribution.
Definition: RosinRammler.C:100
virtual scalar sample() const
Sample the distribution.
Definition: RosinRammler.C:90
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
volScalarField & p
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
#define WarningInFunction
Report a warning using Foam::Warning.
scalar incGamma_P(const scalar a, const scalar x)
Lower incomplete gamma function.
Definition: incGamma.C:450
Namespace for OpenFOAM.
dimensionedScalar exp(const dimensionedScalar &ds)
static void check(const int retVal, const char *what)
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
error FatalError
Type gMin(const FieldField< Field, Type > &f)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Type gMax(const FieldField< Field, Type > &f)
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
Random rndGen
Definition: createFields.H:23