Apéry's constant: Difference between revisions
m
→{{header|Phix}}: use pygments
(Added Algol 68) |
m (→{{header|Phix}}: use pygments) |
||
(3 intermediate revisions by one other user not shown) | |||
Line 256:
</syntaxhighlight>
{{out}}
<div style=background-color:#f8f9fa;padding:1em;white-space:pre-wrap;font-family:monospace;line-height:1.2em;>Apéry's constant via Mathematica's Zeta:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
Apéry's constant via reciprocal cubes (accurate to 6 decimal places):
1.202056<span style=color:red;>4036593442854830714115115999903483212709031775135036540966118572571921400836130084123260473112</span>
Apéry's constant via Markov's summation:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
Apéry's constant via Wedeniwski's summation:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
</
=={{header|Nim}}==
Line 346 ⟶ 349:
First 20 terms of Wedeniwski representation truncated to 100 decimal places:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
</pre>
=={{header|PARI/GP}}==
<syntaxhighlight lang="PARI/GP">
\\ Set the precision to, say, 110 digits
default(realprecision, 110);
\\ Function to display Apéry's constant
my_show_apery_constant(expr, caption) = printf("%s %.101g\n\n", caption, expr);
\\ Apéry's constant via PARI/GP's zeta function
my_show_apery_constant(zeta(3), "Apéry's constant via PARI/GP's Zeta:\n");
\\ Apéry's constant via reciprocal cubes
my_show_apery_constant(sum(k = 1, 1000, 1/k^3), "Apéry's constant via reciprocal cubes:\n");
\\ Apéry's constant via Markov's summation
my_show_apery_constant((5/2) * sum(k = 1, 158, (-1)^(k - 1) * (k!)^2 / ((2*k)! * k^3)), "Apéry's constant via Markov's summation:\n");
\\ Apéry's constant via Wedeniwski's summation
my_show_apery_constant(1/24 * sum(k = 0, 19, (-1)^k * ((2*k + 1)!)^3 * ((2*k)!)^3 * (k!)^3 * (126392*k^5 + 412708*k^4 + 531578*k^3 + 336367*k^2 + 104000*k + 12463) / (((3*k + 2)!) * ((4*k + 3)!)^3)), "Apéry's constant via Wedeniwski's summation:\n");
</syntaxhighlight>
{{out}}
<pre>
Apéry's constant via PARI/GP's Zeta:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
Apéry's constant via reciprocal cubes:
1.2020564036593442854830714115115999903483212709031775135036540966118572571921400836130084123260473112
Apéry's constant via Markov's summation:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
Apéry's constant via Wedeniwski's summation:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
</pre>
Line 390 ⟶ 430:
Ugh. If you ran this on the James Webb, you might just be able to pick out a faint small print outline of the word "elegant".<br>
Still, at least it is not like you do this sort of thing every day... and I got to fix a couple of bugs in my mpfr.js code.
<!--
<syntaxhighlight lang="phix">
with javascript_semantics
requires("1.0.2") -- (missing mpfr_ui_pow_ui() and bug in mpfr_mul_d(), both in mpfr.js)
include mpfr.e
mpfr_set_default_precision(-100)
mpfr {d,a,w,t} = mpfr_inits(4)
mpz {z,pk} = mpz_inits(2)
for k=1 to 1000 do
mpfr_ui_pow_ui(t,k,3)
mpfr_si_div(t,1,t)
mpfr_add(d,d,t)
end for
for k=1 to 158 do
mpz_fac_ui(z,k)
mpz_mul(z,z,z)
mpfr_set_z(t,z)
mpz_fac_ui(z,2*k)
mpfr_div_z(t,t,z)
mpz_ui_pow_ui(z,k,3)
mpfr_div_z(t,t,z)
if even(k) then
mpfr_sub(a,a,t)
else
mpfr_add(a,a,t)
end if
end for
mpfr_mul_d(a,a,5/2)
for k=0 to 19 do
mpz_ui_pow_ui(z,k,5)
mpz_mul_si(z,z,126392)
mpz_ui_pow_ui(pk,k,4)
mpz_mul_si(pk,pk,412708)
mpz_add(z,z,pk)
mpz_ui_pow_ui(pk,k,3)
mpz_mul_si(pk,pk,531578)
mpz_add(z,z,pk)
mpz_add_si(z,z,k*k*336367)
mpz_add_si(z,z,k*104000)
mpz_add_si(z,z,12463)
mpfr_set_z(t,z)
mpz_fac_ui(z,2*k+1)
mpz_pow_ui(z,z,3)
mpfr_mul_z(t,t,z)
mpz_fac_ui(z,2*k)
mpz_pow_ui(z,z,3)
mpfr_mul_z(t,t,z)
mpz_fac_ui(z,k)
mpz_pow_ui(z,z,3)
mpfr_mul_z(t,t,z)
mpz_fac_ui(z,3*k+2)
mpfr_div_z(t,t,z)
mpz_fac_ui(z,4*k+3)
mpz_pow_ui(z,z,3)
mpfr_div_z(t,t,z)
if odd(k) then
mpfr_sub(w,w,t)
else
mpfr_add(w,w,t)
end if
end for
mpfr_div_si(w,w,24)
constant fmt = """
Actual value to 100 decimal places:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
First 1000 terms of zeta(3) truncated to 100 decimal places. (accurate to 6 decimal places):
%s
First 158 terms of Markov / Apery representation truncated to 100 decimal places:
%s
First 20 terms of Wedeniwski representation truncated to 100 decimal places:
%s
"""
string direct = mpfr_get_fixed(d,100),
mapery = mpfr_get_fixed(a,100),
wdnski = mpfr_get_fixed(w,100)
printf(1,fmt,{direct,mapery,wdnski})
</syntaxhighlight>
{{out}}
<pre>
Line 489 ⟶ 530:
<small>Last digit of the 1000 terms line is 2 under pwa/p2js...</small><br>
As per Wren, you can verify or completely replace all this with mpfr_zeta_ui(w,3) <small>[on desktop/Phix only, not supported under pwa/p2js]</small>
=={{header|Python}}==
<syntaxhighlight lang="Python">
from sympy import zeta, factorial
from decimal import Decimal, getcontext
# Set the desired precision
getcontext().prec = 120
def my_sympy_format_to_decimal(sympy_result):
return Decimal(str(sympy_result.evalf(getcontext().prec)))
def print_apery_constant(description, value):
print(f"{description}:\n{str(value)[:102]}")
# Apéry's constant via SymPy's zeta function
zeta_3_str = str(zeta(3).evalf(getcontext().prec))
zeta_3_decimal = Decimal(zeta_3_str)
print_apery_constant("Apéry's constant via SymPy's zeta", zeta_3_decimal)
# Apéry's constant via Riemann summation of 1/(k cubed)
def apery_r(nterms=1_000):
total = sum(Decimal('1') / Decimal(k) ** 3 for k in range(1, nterms + 1))
return total
print_apery_constant("Apéry's constant via reciprocal cubes", apery_r())
# Apéry's constant via Markov's summation
def apery_m(nterms=158):
total = Decimal(2.5) * sum(
(Decimal(1) if k % 2 != 0 else Decimal(-1)) *
my_sympy_format_to_decimal(factorial(k) ** 2) /
my_sympy_format_to_decimal(factorial(2*k) * (k ** 3) )
for k in range(1, nterms + 1)
)
return total
print_apery_constant("Apéry's constant via Markov's summation", apery_m())
# Apéry's constant via Wedeniwski's summation
def apery_w(nterms=20):
total = Decimal('1') / Decimal('24') * sum(
(Decimal('1') if k % 2 == 0 else Decimal('-1')) *
my_sympy_format_to_decimal(factorial(2 * k + 1)) ** 3 *
my_sympy_format_to_decimal(factorial(2 * k)) ** 3 *
my_sympy_format_to_decimal(factorial(k)) ** 3 *
(Decimal('126392') * Decimal(k) ** 5 +
Decimal('412708') * Decimal(k) ** 4 +
Decimal('531578') * Decimal(k) ** 3 +
Decimal('336367') * Decimal(k) ** 2 +
Decimal('104000') * Decimal(k) +
Decimal('12463')) /
(my_sympy_format_to_decimal(factorial(3 * k + 2)) * my_sympy_format_to_decimal(factorial(4 * k + 3)) ** 3)
for k in range(0, nterms + 1)
)
return total
print_apery_constant("Apéry's constant via Wedeniwski's summation", apery_w())
</syntaxhighlight>
{{out}}
<pre>
Apéry's constant via SymPy's zeta:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
Apéry's constant via reciprocal cubes:
1.2020564036593442854830714115115999903483212709031775135036540966118572571921400836130084123260473111
Apéry's constant via Markov's summation:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
Apéry's constant via Wedeniwski's summation:
1.2020569031595942853997381615114499907649862923404988817922715553418382057863130901864558736093352581
</pre>
=={{header|Raku}}==
|