Go to the documentation of this file.
59 static scalar
calcQE11(
const scalar a,
const scalar
x,
const int e = 30)
68 for (
n = 1; (2*
n) <=
e;
n++)
70 const scalar a_2nm1 = a_2np1;
71 const scalar b_2nm1 = b_2np1;
73 a_2n = a_2nm1 + (
n - a)*a_2n;
74 b_2n = b_2nm1 + (
n - a)*b_2n;
76 a_2np1 =
x*a_2n +
n*a_2nm1;
77 b_2np1 =
x*b_2n +
n*b_2nm1;
92 static scalar
calcPE15(
const scalar a,
const scalar
x,
const int nmax = 20)
97 for (
int n = 1;
n <= nmax;
n++)
103 const scalar
R = (
exp(-
x)*
pow(
x, a))/tgamma(a);
105 return R/a*(1 +
sum);
110 static scalar
calcQE16(
const scalar a,
const scalar
x,
const int N = 20)
115 for (
int n = 1;
n <= (
N - 1);
n++)
121 const scalar
R = (
exp(-
x)*
pow(
x, a))/tgamma(a);
123 return R/
x*(1 +
sum);
138 static const scalar D0[] =
140 -0.333333333333333E-00,
141 0.833333333333333E-01,
142 -0.148148148148148E-01,
143 0.115740740740741E-02,
144 0.352733686067019E-03,
145 -0.178755144032922E-03,
146 0.391926317852244E-04,
147 -0.218544851067999E-05,
148 -0.185406221071516E-05,
149 0.829671134095309E-06,
150 -0.176659527368261E-06,
151 0.670785354340150E-08,
152 0.102618097842403E-07,
153 -0.438203601845335E-08
156 static const scalar D1[] =
158 -0.185185185185185E-02,
159 -0.347222222222222E-02,
160 0.264550264550265E-02,
161 -0.990226337448560E-03,
162 0.205761316872428E-03,
163 -0.401877572016461E-06,
164 -0.180985503344900E-04,
165 0.764916091608111E-05,
166 -0.161209008945634E-05,
167 0.464712780280743E-08,
168 0.137863344691572E-06,
169 -0.575254560351770E-07,
170 0.119516285997781E-07
173 static const scalar D2[] =
175 0.413359788359788E-02,
176 -0.268132716049383E-02,
177 0.771604938271605E-03,
178 0.200938786008230E-05,
179 -0.107366532263652E-03,
180 0.529234488291201E-04,
181 -0.127606351886187E-04,
182 0.342357873409614E-07,
183 0.137219573090629E-05,
184 -0.629899213838006E-06,
185 0.142806142060642E-06
188 const scalar u = 1/a;
200 + D0[3]*
pow3(z) + D0[2]*
sqr(z) + D0[1]*z + D0[0];
203 D1[4]*
pow4(z) + D1[3]*
pow3(z) + D1[2]*
sqr(z) + D1[1]*z + D1[0];
205 const scalar C2 = D2[1]*z + D2[0];
207 return C2*
sqr(u) + C1*u + C0;
211 const scalar C0 = D0[2]*
sqr(z) + D0[1]*z + D0[0];
212 const scalar C1 = D1[1]*z + D1[0];
213 const scalar C2 = D2[1]*z + D2[0];
215 return C2*
sqr(u) + C1*u + C0;
228 const scalar BIG = 14;
229 const scalar x0 = 17;
230 const scalar e0 = 0.025;
258 for (
int n = 1;
n <= 10;
n++)
263 const scalar J = -a*
sum;
268 return 1 - (
pow(
x, a)*(1 - J))/tgamma(a + 1);
274 const scalar
H = 1/(tgamma(a + 1)) - 1;
276 return (
pow(
x, a)*J -
L)/tgamma(a + 1) -
H;
282 const scalar
R = (
exp(-
x)*
pow(
x, a))/tgamma(a);
289 const scalar
sigma = fabs(1 -
x/a);
296 const scalar
y = a*
phi;
298 const scalar E = 0.5 - (1 -
y/3)*
sqrt(
y/
pi);
325 const scalar
y = a*
phi;
353 const scalar
R = (
exp(-
x)*
pow(
x, a))/tgamma(a);
367 if (a >
x ||
x >= x0)
377 const scalar
R = (
exp(-
x)*
pow(
x, a))/tgamma(a);
389 if (floor(2*a) == 2*a)
396 for (
int n = 0;
n <= (a - 1);
n++)
409 for (
int n = 1;
n <= i;
n++)
418 else if (
x <=
max(a,
log(10.0)))
426 const scalar
R = (
exp(-
x)*
pow(
x, a))/tgamma(a);
const vector L(dict.get< vector >("L"))
scalar incGamma_Q(const scalar a, const scalar x)
Upper incomplete gamma function.
static scalar calcTE18(const scalar a, const scalar e0, const scalar x, const scalar lambda, const scalar sigma, const scalar phi)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
dimensionedScalar exp(const dimensionedScalar &ds)
static scalar calcQE16(const scalar a, const scalar x, const int N=20)
dimensionedScalar erf(const dimensionedScalar &ds)
scalar incGammaRatio_Q(const scalar a, const scalar x)
Normalized upper incomplete gamma function.
dimensionedScalar pow4(const dimensionedScalar &ds)
#define R(A, B, C, D, E, F, K, M)
dimensionedScalar pow6(const dimensionedScalar &ds)
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensionedScalar erfc(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
scalar incGamma_P(const scalar a, const scalar x)
Lower incomplete gamma function.
dimensionedScalar pow5(const dimensionedScalar &ds)
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
dimensionedScalar log(const dimensionedScalar &ds)
static scalar calcPE15(const scalar a, const scalar x, const int nmax=20)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr scalar pi(M_PI)
scalar incGammaRatio_P(const scalar a, const scalar x)
Normalized lower incomplete gamma function.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
label factorial(label n)
Evaluate n! : 0 < n <= 12.
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
const Vector< label > N(dict.get< Vector< label >>("N"))
static scalar calcQE11(const scalar a, const scalar x, const int e=30)