exponential.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-2016 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 "exponential.H"
31
32// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33
34namespace Foam
35{
36namespace distributionModels
37{
40}
41}
42
43// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44
46(
47 const dictionary& dict,
49)
50:
51 distributionModel(typeName, dict, rndGen),
52 lambda_(distributionModelDict_.get<scalar>("lambda"))
53{
54 if (lambda_ < VSMALL)
55 {
57 << "Rate parameter cannot be equal to or less than zero:" << nl
58 << " lambda = " << lambda_
59 << exit(FatalError);
60 }
61
62 check();
63}
64
65
67:
69 lambda_(p.lambda_)
70{}
71
72
73// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
74
76{
77 const scalar u = rndGen_.sample01<scalar>();
78 const scalar qMin = exp(-lambda_*minValue_);
79 const scalar qMax = exp(-lambda_*maxValue_);
80 return -(scalar(1)/lambda_)*log(qMin + u*(qMax - qMin));
81}
82
83
85{
86 return scalar(1)/lambda_;
87}
88
89
90// ************************************************************************* //
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...
virtual void check() const
Check that the distribution model is valid.
Particle-size distribution model wherein random samples are drawn from the doubly-truncated exponenti...
Definition: exponential.H:192
virtual scalar meanValue() const
Return the theoretical mean of the distribution.
Definition: exponential.C:84
virtual scalar sample() const
Sample the distribution.
Definition: exponential.C:75
#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
Namespace for OpenFOAM.
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
error FatalError
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
dictionary dict
Random rndGen
Definition: createFields.H:23