
You don't need to compare prev
and new
during each iteration.
The difference between the new and the previous sum is simply the current term : $$\frac{8}{(2i+1)^2}$$
If you want this term to be smaller than error
, you can solve:
$$\mathrm{error} > \frac{8}{(2i+1)^2}\\
\iff (2i+1)^2 > \frac{8}{error}\\
\iff 2i+1 > \sqrt{\frac{8}{error}}\\
\iff i > \frac{\sqrt{\frac{8}{error}} - 1}{2}\\
$$
Now that you know how many terms your series should have, you can return the result directly:
def approx_pi_squared(error):
n = int(((8 / error)**0.5 - 1) / 2) + 1
return sum(8 / (2 * i + 1)**2 for i in range(n))
Better formulas
Adding a delta
Note that error
represents how small the terms are, not how close approx_pi_squared
is from π²:
>>> import math
>>> approx_pi_squared(1e-10)
9.869590258918535
>>> math.pi**2 - approx_pi_squared(1e-7)
0.0004472271895057389
>>> math.pi**2 - approx_pi_squared(1e-10)
1.414217082285063e-05
Even with more than 140 000 terms, the series only gives the 3 first digits of π². This formula is very simple but converges too slowly.
What's very interesting, though, is that the difference between math.pi**2
and approx_pi_squared(error)
seems very close to \$\sqrt{2\mathrm{error}}\$. It seems to hold true for any error
, so we could update the function:
def approx_pi_squared(error):
n = int(((8 / error)**0.5 - 1) / 2) + 1
delta = (2 * error)**0.5
return sum(8 / (2 * i + 1)**2 for i in range(n)) + delta
approx_pi_squared(1e-10)
now returns 10 correct digits for π².
This new formula is unproven, so use at your own risk!
BBP-Type Formula
There are many π² formulas, so feel free to pick another one. For example:
def approx_pi_squared(error):
n = int(((12 / error)**0.5 - 1) / 2) + 1
return 12 * sum((-1)**i / (i + 1)**2 for i in range(n))
error
seems to have the same order of magnitude as math.pi**2 - approx_pi_squared(error)
now:
>>> math.pi**2 - approx_pi_squared(1e-9)
2.0001476030984122e-09
>>> math.pi**2 - approx_pi_squared(1e-10)
-1.9977974829998857e-10
delta
looks like (-1)**n * 2 * error
.
With Sympy
You can delegate the job to sympy
and be sure you'll get a correct result with arbitrary precision:
>>> from sympy import Sum, N, pi, oo, init_printing
>>> from sympy.abc import i, n
>>> init_printing()
>>> Sum(8/(2*i+1)**2, (i, 0, oo))
∞
____
╲
╲ 8
╲ ──────────
╱ 2
╱ (2⋅i + 1)
╱
‾‾‾‾
i = 0
>>> N(Sum(8/(2*i+1)**2, (i, 0, oo)), 100)
9.869604401089358618834490999876151135313699407240790626413349376220044822419205243001773403718552232
>>> N(pi**2, 100)
9.869604401089358618834490999876151135313699407240790626413349376220044822419205243001773403718552232
return prev
(8 / (n * n) <= error < 8 / ((n - 1) * (n - 1))
) \$\endgroup\$