Euler's constant 0.5772...: Difference between revisions

no edit summary
No edit summary
Line 935:
gamma 0.5772156649015328606065120900824024310421593359399235988057672348848677267776646709369470632917467495 (mpfr_const_euler)
</pre>
 
=={{header|Processing}}==
<lang Processing>
/*********************************************
Subject: Comparing five methods for
computing Euler's constant 0.5772...
// https://rosettacode.org/wiki/Euler%27s_constant_0.5772...
--------------------------------------------*/
double a, b, h, n2, r, u, v;
float floatA, floatB, floatN2;
int k, k2, m, n;
double eps = 1e-6;
 
void setup() {
size(100, 100);
noLoop();
}
 
void draw() {
println("From the definition, err. 3e-10\n");
 
n = 400;
 
h = 1;
 
for (int k = 2; k <= n; k++) {
h += 1.0 / k;
}
//faster convergence: Negoi, 1997
a = log(n +.5 + 1.0 / (24*n));
 
println("Hn ", h);
println("gamma ", h - a);
println("k = ", n);
println("");
 
 
println("Sweeney, 1963, err. idem");
n = 21;
 
double s[] = {0, n};
r = n;
k = 1;
while (r > eps) {
k ++;
r *= (double) n / k;
s[k & 1] = s[k & 1] + r / k;
}
 
// println("gamma %.16f\nk = %d\n\n", s[1] - s[0] - log(n), k);
 
println("Hn ", h);
println("gamma ", s[1] - s[0] - log(n));
println("k = ", k);
println("");
 
println("Bailey, 1988");
n = 5;
floatA = 1;
h = 1;
floatN2 = pow(2, n);
r = 1;
k = 1;
while (abs(floatB - floatA) > eps) {
k += 1;
r *= floatN2 / k;
h += 1.0 / k;
floatB = floatA;
floatA += r * h;
}
floatA *= floatN2 / exp(floatN2);
 
println("gamma ", floatA - n * log(2));
println("k = ", k);
println("");
 
println("Brent-McMillan, 1980");
 
n = 13;
 
floatA = -log(n);
floatB = 1;
u = a;
v = b;
n2 = n * n;
k2 = 0;
k = 0;
 
 
while (abs(floatA) > eps) {
k2 += 2*k + 1;
k += 1;
floatA *= n2 / k;
floatB *= n2 / k2;
floatA = (floatA + floatB) / k;
u += floatA;
v += floatB;
}
println("gamma ", u / v);
println("k ", k);
 
 
 
println("How Euler did it in 1735\n");
//Bernoulli numbers with even indices
 
double[] B2 = new double[11];
 
B2[1] = 1.0;
B2[2] = 1.0/6;
B2[3] = -1.0/30;
B2[4] = 1.0/42;
B2[5] = -1.0/30;
B2[6] = 5.0/66;
B2[7] = -691.0/2730;
B2[8] = 7.0/6;
B2[9] = -3617.0/510;
B2[10]= 43867.0/798;
 
m = 7;
n = 10;
 
//n-th harmonic number
h = 1;
for (k = 2; k <= n; k++) {
h += 1.0 / k;
}
println("Hn ", h);
 
h -= log(n);
println(" -ln ", h);
 
//expansion C = -digamma(1)
a = -1.0 / (2*n);
n2 = n * n;
r = 1;
for (k = 1; k <= m; k++) {
r *= n2;
a += B2[k] / (2*k * r);
}
 
println("err ",a);
println("gamma ", h + a );
println("k = ", n + m);
println("");
println("C = 0.57721566490153286...\n");
}
 
</lang>
 
{{out}}
<pre>
From the definition, err. 3e-10
 
Hn 6.569929751800373
gamma 0.577215823577717
k = 400
 
Sweeney, 1963, err. idem
Hn 6.569929751800373
gamma 0.5772155784070492
k = 68
 
Bailey, 1988
gamma 0.5772176
k = 67
 
Brent-McMillan, 1980
gamma 0.5772157269353247
k 40
How Euler did it in 1735
 
Hn 2.9289682805538177
-ln 0.6263831555843353
err -0.044995839604388355
gamma 0.581387315979947
k = 17
 
C = 0.57721566490153286...
 
</pre>
 
 
 
 
 
=={{header|Raku}}==