diff --git a/src/OpenFOAM/Make/files b/src/OpenFOAM/Make/files
index bcd75424a5..5e9e8d99b7 100644
--- a/src/OpenFOAM/Make/files
+++ b/src/OpenFOAM/Make/files
@@ -37,6 +37,7 @@ $(ints)/lists/labelListIOList.C
primitives/Scalar/doubleScalar/doubleScalar.C
primitives/Scalar/floatScalar/floatScalar.C
primitives/Scalar/scalar/scalar.C
+primitives/Scalar/scalar/invIncGamma.C
primitives/Scalar/lists/scalarList.C
primitives/Scalar/lists/scalarIOList.C
primitives/Scalar/lists/scalarListIOList.C
diff --git a/src/OpenFOAM/global/constants/mathematical/mathematicalConstants.H b/src/OpenFOAM/global/constants/mathematical/mathematicalConstants.H
index d286493919..fa4850fb5a 100644
--- a/src/OpenFOAM/global/constants/mathematical/mathematicalConstants.H
+++ b/src/OpenFOAM/global/constants/mathematical/mathematicalConstants.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -52,6 +52,9 @@ namespace mathematical
const scalar twoPi(2*pi);
const scalar piByTwo(0.5*pi);
+ //- Euler's constant
+ const scalar Eu(0.57721566490153286060651209);
+
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace mathematical
diff --git a/src/OpenFOAM/primitives/Scalar/Scalar.H b/src/OpenFOAM/primitives/Scalar/Scalar.H
index d3cf2beca4..34a7051535 100644
--- a/src/OpenFOAM/primitives/Scalar/Scalar.H
+++ b/src/OpenFOAM/primitives/Scalar/Scalar.H
@@ -134,6 +134,7 @@ transFunc(atanh)
transFunc(erf)
transFunc(erfc)
transFunc(lgamma)
+transFunc(tgamma)
transFunc(j0)
transFunc(j1)
diff --git a/src/OpenFOAM/primitives/Scalar/scalar/invIncGamma.C b/src/OpenFOAM/primitives/Scalar/scalar/invIncGamma.C
new file mode 100644
index 0000000000..aa44c195df
--- /dev/null
+++ b/src/OpenFOAM/primitives/Scalar/scalar/invIncGamma.C
@@ -0,0 +1,316 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+Global
+ Foam::invIncGamma
+
+Description
+ Function to calculate the inverse of the
+ normalized incomplete gamma function.
+
+ The algorithm is described in detain in reference:
+ \verbatim
+ DiDonato, A. R., & Morris Jr, A. H. (1986).
+ Computation of the incomplete gamma function ratios and their inverse.
+ ACM Transactions on Mathematical Software (TOMS), 12(4), 377-393.
+ \endverbatim
+
+ All equation numbers in the following code refer to the above text and
+ the algorithm in the function 'invIncGamma' is described in section 4.
+
+\*---------------------------------------------------------------------------*/
+
+#include "mathematicalConstants.H"
+using namespace Foam::constant::mathematical;
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+static scalar minimaxs(const scalar P)
+{
+ // Eqn. 32
+
+ static const scalar a[] =
+ {
+ 3.31125922108741,
+ 11.6616720288968,
+ 4.28342155967104,
+ 0.213623493715853
+ };
+
+ static const scalar b[] =
+ {
+ 6.61053765625462,
+ 6.40691597760039,
+ 1.27364489782223,
+ 0.03611708101884203
+ };
+
+ const scalar t = P < 0.5 ? sqrt(-2*log(P)) : sqrt(-2*log(1 - P));
+
+ const scalar s =
+ t
+ - (a[0] + t*(a[1] + t*(a[2] + t*a[3])))
+ /(1 + t*(b[0] + t*(b[1] + t*(b[2] + t*b[3]))));
+
+ return P < 0.5 ? -s : s;
+}
+
+
+static scalar Sn(const scalar a, const scalar x)
+{
+ //- Eqn. 34
+
+ scalar Sn = 1;
+ scalar Si = 1;
+
+ for (int i=1; i<100; i++)
+ {
+ Si *= x/(a + i);
+ Sn += Si;
+
+ if (Si < 1e-4) break;
+ }
+
+ return Sn;
+}
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+//- Inverse normalized incomplete gamma function
+Foam::scalar Foam::invIncGamma(const scalar a, const scalar P)
+{
+ const scalar Q = 1 - P;
+
+ if (a == 1)
+ {
+ return -log(Q);
+ }
+ else if (a < 1)
+ {
+ const scalar Ga = tgamma(a);
+ const scalar B = Q*Ga;
+
+ if (B > 0.6 || (B >= 0.45 && a >= 0.3))
+ {
+ // Eqn. 21
+ const scalar u =
+ (B*Q > 1e-8) ? pow(P*Ga*a, 1/a) : exp((-Q/a) - Eu);
+
+ return u/(1 - (u/(a + 1)));
+ }
+ else if (a < 0.3 && B >= 0.35)
+ {
+ // Eqn. 22
+ const scalar t = exp(-Eu - B);
+ const scalar u = t*exp(t);
+
+ return t*exp(u);
+ }
+ else if (B > 0.15 || a >= 0.3)
+ {
+ // Eqn. 23
+ const scalar y = -log(B);
+ const scalar u = y - (1 - a)*log(y);
+
+ return y - (1 - a)*log(u) - log(1 + (1 - a)/(1 + u));
+ }
+ else if (B > 0.1)
+ {
+ // Eqn. 24
+ const scalar y = -log(B);
+ const scalar u = y - (1 - a)*log(y);
+
+ return y
+ - (1 - a)*log(u)
+ - log
+ (
+ (sqr(u) + 2*(3 - a)*u + (2 - a)*(3 - a))
+ /(sqr(u) + (5 - a)*u + 2)
+ );
+ }
+ else
+ {
+ // Eqn. 25:
+ const scalar y = -log(B);
+ const scalar c1 = (a - 1)*log(y);
+ const scalar c12 = c1*c1;
+ const scalar c13 = c12*c1;
+ const scalar c14 = c12*c12;
+ const scalar a2 = a*a;
+ const scalar a3 = a2*a;
+ const scalar c2 = (a - 1)*(1 + c1);
+ const scalar c3 = (a - 1)*(-(c12/2) + (a - 2)*c1 + (3*a - 5)/2);
+ const scalar c4 =
+ (a - 1)
+ *(
+ (c13/3)
+ - (3*a - 5)*c12/2
+ + (a2 - 6*a + 7)*c1
+ + (11*a2 - 46*a + 47)/6
+ );
+ const scalar c5 =
+ (a - 1)*(-(c14/4)
+ + (11*a - 17)*c13/6
+ + (-3*a2 + 13*a - 13)*c12
+ + (2*a3 - 25*a2 + 72*a - 61)*c1/2
+ + (25*a3 - 195*a2 + 477*a - 379)/12);
+ const scalar y2 = y*y;
+ const scalar y3 = y2*y;
+ const scalar y4 = y2*y2;
+
+ return y + c1 + (c2/y) + (c3/y2) + (c4/y3) + (c5/y4);
+ }
+ }
+ else
+ {
+ // Eqn. 31:
+ scalar s = minimaxs(P);
+
+ const scalar s2 = sqr(s);
+ const scalar s3 = s*s2;
+ const scalar s4 = s2*s2;
+ const scalar s5 = s*s4;
+ const scalar sqrta = sqrt(a);
+
+ const scalar w =
+ a + s*sqrta + (s2 - 1)/3
+ + (s3 - 7*s)/(36*sqrta)
+ - (3*s4 + 7*s2 - 16)/(810*a)
+ + (9*s5 + 256*s3 - 433*s)/(38880*a*sqrta);
+
+ if (a >= 500 && mag(1 - w/a) < 1e-6)
+ {
+ return w;
+ }
+ else if (P > 0.5)
+ {
+ if (w < 3*a)
+ {
+ return w;
+ }
+ else
+ {
+ const scalar D = max(scalar(2), scalar(a*(a - 1)));
+ const scalar lnGa = lgamma(a);
+ const scalar lnB = log(Q) + lnGa;
+
+ if (lnB < -2.3*D)
+ {
+ // Eqn. 25:
+ const scalar y = -lnB;
+ const scalar c1 = (a - 1)*log(y);
+ const scalar c12 = c1*c1;
+ const scalar c13 = c12*c1;
+ const scalar c14 = c12*c12;
+ const scalar a2 = a*a;
+ const scalar a3 = a2*a;
+
+ const scalar c2 = (a - 1)*(1 + c1);
+ const scalar c3 =
+ (a - 1)
+ *(
+ - (c12/2)
+ + (a - 2)*c1
+ + (3*a - 5)/2
+ );
+ const scalar c4 =
+ (a - 1)
+ *(
+ (c13/3)
+ - (3*a - 5)*c12/2
+ + (a2 - 6*a + 7)*c1
+ + (11*a2 - 46*a + 47)/6
+ );
+ const scalar c5 =
+ (a - 1)
+ *(
+ - (c14/4)
+ + (11*a - 17)*c13/6
+ + (-3*a2 + 13*a - 13)*c12
+ + (2*a3 - 25*a2 + 72*a - 61)*c1/2
+ + (25*a3 - 195*a2 + 477*a - 379)/12
+ );
+
+ const scalar y2 = y*y;
+ const scalar y3 = y2*y;
+ const scalar y4 = y2*y2;
+
+ return y + c1 + (c2/y) + (c3/y2) + (c4/y3) + (c5/y4);
+ }
+ else
+ {
+ // Eqn. 33
+ const scalar u =
+ -lnB + (a - 1)*log(w) - log(1 + (1 - a)/(1 + w));
+
+ return -lnB + (a - 1)*log(u) - log(1 + (1 - a)/(1 + u));
+ }
+ }
+ }
+ else
+ {
+ scalar z = w;
+ const scalar ap1 = a + 1;
+
+ if (w < 0.15*ap1)
+ {
+ // Eqn. 35
+ const scalar ap2 = a + 2;
+ const scalar v = log(P) + lgamma(ap1);
+ z = exp((v + w)/a);
+ s = log1p(z/ap1*(1 + z/ap2));
+ z = exp((v + z - s)/a);
+ s = log1p(z/ap1*(1 + z/ap2));
+ z = exp((v + z - s)/a);
+ s = log1p(z/ap1*(1 + z/ap2*(1 + z/(a + 3))));
+ z = exp((v + z - s)/a);
+ }
+
+ if (z <= 0.01*ap1 || z > 0.7*ap1)
+ {
+ return z;
+ }
+ else
+ {
+ // Eqn. 36
+ const scalar lnSn = log(Sn(a, z));
+ const scalar v = log(P) + lgamma(ap1);
+ z = exp((v + z - lnSn)/a);
+
+ return z*(1 - (a*log(z) - z - v + lnSn)/(a - z));
+ }
+ }
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/primitives/Scalar/scalar/scalar.H b/src/OpenFOAM/primitives/Scalar/scalar/scalar.H
index f4a4a0f63b..7c7c00cefc 100644
--- a/src/OpenFOAM/primitives/Scalar/scalar/scalar.H
+++ b/src/OpenFOAM/primitives/Scalar/scalar/scalar.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -81,6 +81,14 @@ namespace Foam
#endif
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// Additional transcendental functions
+
+namespace Foam
+{
+ //- Inverse normalized incomplete gamma function
+ scalar invIncGamma(const scalar a, const scalar P);
+}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/distributionModels/Make/files b/src/lagrangian/distributionModels/Make/files
index c7b621cf58..aee429b06c 100644
--- a/src/lagrangian/distributionModels/Make/files
+++ b/src/lagrangian/distributionModels/Make/files
@@ -7,6 +7,7 @@ general/general.C
multiNormal/multiNormal.C
normal/normal.C
RosinRammler/RosinRammler.C
+massRosinRammler/massRosinRammler.C
uniform/uniform.C
LIB = $(FOAM_LIBBIN)/libdistributionModels
diff --git a/src/lagrangian/distributionModels/RosinRammler/RosinRammler.C b/src/lagrangian/distributionModels/RosinRammler/RosinRammler.C
index f2e5ac0965..2e5a045619 100644
--- a/src/lagrangian/distributionModels/RosinRammler/RosinRammler.C
+++ b/src/lagrangian/distributionModels/RosinRammler/RosinRammler.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -30,11 +30,11 @@ License
namespace Foam
{
- namespace distributionModels
- {
- defineTypeNameAndDebug(RosinRammler, 0);
- addToRunTimeSelectionTable(distributionModel, RosinRammler, dictionary);
- }
+namespace distributionModels
+{
+ defineTypeNameAndDebug(RosinRammler, 0);
+ addToRunTimeSelectionTable(distributionModel, RosinRammler, dictionary);
+}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/distributionModels/RosinRammler/RosinRammler.H b/src/lagrangian/distributionModels/RosinRammler/RosinRammler.H
index 64cc13ca3f..972a48d047 100644
--- a/src/lagrangian/distributionModels/RosinRammler/RosinRammler.H
+++ b/src/lagrangian/distributionModels/RosinRammler/RosinRammler.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -100,7 +100,7 @@ public:
// Member Functions
- //- Sample the distributionModel
+ //- Sample the distributionModel
virtual scalar sample() const;
//- Return the minimum value
diff --git a/src/lagrangian/distributionModels/distributionModel/distributionModel.H b/src/lagrangian/distributionModels/distributionModel/distributionModel.H
index fe5e876544..5e712e10ec 100644
--- a/src/lagrangian/distributionModels/distributionModel/distributionModel.H
+++ b/src/lagrangian/distributionModels/distributionModel/distributionModel.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -36,6 +36,7 @@ Description
- multi-normal
- normal
- Rosin-Rammler
+ - Mass-based Rosin-Rammler
- uniform
The distributionModel is tabulated in equidistant nPoints, in an interval.
diff --git a/src/lagrangian/distributionModels/exponential/exponential.C b/src/lagrangian/distributionModels/exponential/exponential.C
index 01298de33c..013e5c39bb 100644
--- a/src/lagrangian/distributionModels/exponential/exponential.C
+++ b/src/lagrangian/distributionModels/exponential/exponential.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -30,11 +30,11 @@ License
namespace Foam
{
- namespace distributionModels
- {
- defineTypeNameAndDebug(exponential, 0);
- addToRunTimeSelectionTable(distributionModel, exponential, dictionary);
- }
+namespace distributionModels
+{
+ defineTypeNameAndDebug(exponential, 0);
+ addToRunTimeSelectionTable(distributionModel, exponential, dictionary);
+}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/distributionModels/exponential/exponential.H b/src/lagrangian/distributionModels/exponential/exponential.H
index 0bcadfc0b8..ff2473e2b9 100644
--- a/src/lagrangian/distributionModels/exponential/exponential.H
+++ b/src/lagrangian/distributionModels/exponential/exponential.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -93,7 +93,7 @@ public:
// Member Functions
- //- Sample the distributionModel
+ //- Sample the distributionModel
virtual scalar sample() const;
//- Return the minimum value
diff --git a/src/lagrangian/distributionModels/fixedValue/fixedValue.C b/src/lagrangian/distributionModels/fixedValue/fixedValue.C
index 048a9482f3..f8f0c234a9 100644
--- a/src/lagrangian/distributionModels/fixedValue/fixedValue.C
+++ b/src/lagrangian/distributionModels/fixedValue/fixedValue.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -30,11 +30,11 @@ License
namespace Foam
{
- namespace distributionModels
- {
- defineTypeNameAndDebug(fixedValue, 0);
- addToRunTimeSelectionTable(distributionModel, fixedValue, dictionary);
- }
+namespace distributionModels
+{
+ defineTypeNameAndDebug(fixedValue, 0);
+ addToRunTimeSelectionTable(distributionModel, fixedValue, dictionary);
+}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/distributionModels/fixedValue/fixedValue.H b/src/lagrangian/distributionModels/fixedValue/fixedValue.H
index 3de10f43b6..32fe917b56 100644
--- a/src/lagrangian/distributionModels/fixedValue/fixedValue.H
+++ b/src/lagrangian/distributionModels/fixedValue/fixedValue.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -84,7 +84,7 @@ public:
// Member Functions
- //- Sample the distributionModel
+ //- Sample the distributionModel
virtual scalar sample() const;
//- Return the minimum value
diff --git a/src/lagrangian/distributionModels/general/general.C b/src/lagrangian/distributionModels/general/general.C
index 5199ff8074..78663997a1 100644
--- a/src/lagrangian/distributionModels/general/general.C
+++ b/src/lagrangian/distributionModels/general/general.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -30,11 +30,11 @@ License
namespace Foam
{
- namespace distributionModels
- {
- defineTypeNameAndDebug(general, 0);
- addToRunTimeSelectionTable(distributionModel, general, dictionary);
- }
+namespace distributionModels
+{
+ defineTypeNameAndDebug(general, 0);
+ addToRunTimeSelectionTable(distributionModel, general, dictionary);
+}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/distributionModels/massRosinRammler/massRosinRammler.C b/src/lagrangian/distributionModels/massRosinRammler/massRosinRammler.C
new file mode 100644
index 0000000000..015343b711
--- /dev/null
+++ b/src/lagrangian/distributionModels/massRosinRammler/massRosinRammler.C
@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+\*---------------------------------------------------------------------------*/
+
+#include "massRosinRammler.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace distributionModels
+{
+ defineTypeNameAndDebug(massRosinRammler, 0);
+ addToRunTimeSelectionTable(distributionModel, massRosinRammler, dictionary);
+}
+}
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::distributionModels::massRosinRammler::massRosinRammler
+(
+ const dictionary& dict,
+ cachedRandom& rndGen
+)
+:
+ distributionModel(typeName, dict, rndGen),
+ minValue_(readScalar(distributionModelDict_.lookup("minValue"))),
+ maxValue_(readScalar(distributionModelDict_.lookup("maxValue"))),
+ d_(readScalar(distributionModelDict_.lookup("d"))),
+ n_(readScalar(distributionModelDict_.lookup("n")))
+{
+ check();
+}
+
+
+Foam::distributionModels::massRosinRammler::massRosinRammler
+(
+ const massRosinRammler& p
+)
+:
+ distributionModel(p),
+ minValue_(p.minValue_),
+ maxValue_(p.maxValue_),
+ d_(p.d_),
+ n_(p.n_)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::distributionModels::massRosinRammler::~massRosinRammler()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+Foam::scalar Foam::distributionModels::massRosinRammler::sample() const
+{
+ scalar d;
+
+ // Re-sample if the calculated d is out of the physical range
+ do
+ {
+ const scalar a = 3/n_ + 1;
+ const scalar P = rndGen_.sample01();
+ const scalar x = invIncGamma(a, P);
+ d = d_*pow(x, 1/n_);
+ } while (d < minValue_ || d > maxValue_);
+
+ return d;
+}
+
+
+Foam::scalar Foam::distributionModels::massRosinRammler::minValue() const
+{
+ return minValue_;
+}
+
+
+Foam::scalar Foam::distributionModels::massRosinRammler::maxValue() const
+{
+ return maxValue_;
+}
+
+
+Foam::scalar Foam::distributionModels::massRosinRammler::meanValue() const
+{
+ return d_;
+}
+
+
+// ************************************************************************* //
diff --git a/src/lagrangian/distributionModels/massRosinRammler/massRosinRammler.H b/src/lagrangian/distributionModels/massRosinRammler/massRosinRammler.H
new file mode 100644
index 0000000000..f1745326b6
--- /dev/null
+++ b/src/lagrangian/distributionModels/massRosinRammler/massRosinRammler.H
@@ -0,0 +1,137 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+Class
+ Foam::massRosinRammler
+
+Description
+ Mass-based Rosin-Rammler distributionModel.
+
+ Corrected form of the Rosin-Rammler distribution taking into account the
+ varying number of particels per parces for for fixed-mass parcels. This
+ distribution should be used when
+ \verbatim
+ parcelBasisType mass;
+ \endverbatim
+
+ See equation 10 in reference:
+ \verbatim
+ Yoon, S. S., Hewson, J. C., DesJardin, P. E., Glaze, D. J.,
+ Black, A. R., & Skaggs, R. R. (2004).
+ Numerical modeling and experimental measurements of a high speed
+ solid-cone water spray for use in fire suppression applications.
+ International Journal of Multiphase Flow, 30(11), 1369-1388.
+ \endverbatim
+
+SourceFiles
+ massRosinRammler.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef massRosinRammler_H
+#define massRosinRammler_H
+
+#include "distributionModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace distributionModels
+{
+
+/*---------------------------------------------------------------------------*\
+ Class massRosinRammler Declaration
+\*---------------------------------------------------------------------------*/
+
+class massRosinRammler
+:
+ public distributionModel
+{
+ // Private data
+
+ //- Distribution minimum
+ scalar minValue_;
+
+ //- Distribution maximum
+ scalar maxValue_;
+
+ //- Characteristic droplet size
+ scalar d_;
+
+ //- Empirical dimensionless constant to specify the distribution width,
+ // sometimes referred to as the dispersion coefficient
+ scalar n_;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("massRosinRammler");
+
+
+ // Constructors
+
+ //- Construct from components
+ massRosinRammler(const dictionary& dict, cachedRandom& rndGen);
+
+ //- Construct copy
+ massRosinRammler(const massRosinRammler& p);
+
+ //- Construct and return a clone
+ virtual autoPtr clone() const
+ {
+ return autoPtr(new massRosinRammler(*this));
+ }
+
+
+ //- Destructor
+ virtual ~massRosinRammler();
+
+
+ // Member Functions
+
+ //- Sample the distributionModel
+ virtual scalar sample() const;
+
+ //- Return the minimum value
+ virtual scalar minValue() const;
+
+ //- Return the maximum value
+ virtual scalar maxValue() const;
+
+ //- Return the mean value
+ virtual scalar meanValue() const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace distributionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/lagrangian/distributionModels/multiNormal/multiNormal.C b/src/lagrangian/distributionModels/multiNormal/multiNormal.C
index 9a8bf423ec..1cbb5a5856 100644
--- a/src/lagrangian/distributionModels/multiNormal/multiNormal.C
+++ b/src/lagrangian/distributionModels/multiNormal/multiNormal.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -30,11 +30,11 @@ License
namespace Foam
{
- namespace distributionModels
- {
- defineTypeNameAndDebug(multiNormal, 0);
- addToRunTimeSelectionTable(distributionModel, multiNormal, dictionary);
- }
+namespace distributionModels
+{
+ defineTypeNameAndDebug(multiNormal, 0);
+ addToRunTimeSelectionTable(distributionModel, multiNormal, dictionary);
+}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/distributionModels/normal/normal.C b/src/lagrangian/distributionModels/normal/normal.C
index a585ee59a1..1fc16376a9 100644
--- a/src/lagrangian/distributionModels/normal/normal.C
+++ b/src/lagrangian/distributionModels/normal/normal.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -31,11 +31,11 @@ License
namespace Foam
{
- namespace distributionModels
- {
- defineTypeNameAndDebug(normal, 0);
- addToRunTimeSelectionTable(distributionModel, normal, dictionary);
- }
+namespace distributionModels
+{
+ defineTypeNameAndDebug(normal, 0);
+ addToRunTimeSelectionTable(distributionModel, normal, dictionary);
+}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/distributionModels/uniform/uniform.C b/src/lagrangian/distributionModels/uniform/uniform.C
index 05f7a21bef..72a0d2fc3d 100644
--- a/src/lagrangian/distributionModels/uniform/uniform.C
+++ b/src/lagrangian/distributionModels/uniform/uniform.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -30,11 +30,11 @@ License
namespace Foam
{
- namespace distributionModels
- {
- defineTypeNameAndDebug(uniform, 0);
- addToRunTimeSelectionTable(distributionModel, uniform, dictionary);
- }
+namespace distributionModels
+{
+ defineTypeNameAndDebug(uniform, 0);
+ addToRunTimeSelectionTable(distributionModel, uniform, dictionary);
+}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/distributionModels/uniform/uniform.H b/src/lagrangian/distributionModels/uniform/uniform.H
index 21d0e1c14b..20eb4a9c53 100644
--- a/src/lagrangian/distributionModels/uniform/uniform.H
+++ b/src/lagrangian/distributionModels/uniform/uniform.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -88,7 +88,7 @@ public:
// Member Functions
- //- Sample the distributionModel
+ //- Sample the distributionModel
virtual scalar sample() const;
//- Return the minimum value