Paso 12: Codificación fórmula de Gauss
En nuestro paso anterior encontramos la fórmula:
Π/4=12*arctan(1/18) + 8*arctan(1/57) - 5*arctan(1/239)
podemos fácilmente convertir esto en una función en python si asumimos que hemos definido una función arctan:
Ahora necesitamos crear una función arctan. El método más rápido para implementar sería utilizar la fórmula que vimos en el paso 9:
arctan(x) = x - (x³/3) + (x⁵/5) - (x⁷/7) + (x⁹/9) - (x¹¹/11)...
sin embargo ver como sólo estamos calcular arctans de números en forma de 1 / x tiene más sentido para redefinir la fórmula como:
arctan(1/x) = (1 / x)-(1/3 x ³) + (1/5x⁵) - (1/7x⁷) + (1/9x⁹) - (1/11 x ¹¹)...
Esto nos da la siguiente función en Python:
Esta función es un poco más avanzada que las anteriores funciones no requiere el usuario introducir cuantas iteraciones a ejecutar. En cambio se ve en el valor de getcontext () .prec y busca el valor más pequeño que puede distinguir el programa Python desde cero (por ejemplo si getcontext () .prec = 2 el valor más pequeño que puede distinguir Python desde 0 es 0,1).
Luego compara el valor del último término en la secuencia: si es menor que el valor más pequeño que puede distinguir python desde 0 y no tiene ningún sentido continuar por lo que el programa se detiene.
Si juntamos las dos funciones tenemos el siguiente programa (gauss_pi_method.py):
Si ejecuta esto para encontrar 10.000 dígitos de π debe tomar menos de un minuto. En mi equipo tardó 17 segundos. Es mucho mejor que el anterior mejor programa que tenía, el uno basado en el método del polígono, que tuvo 43 segundos en mi computadora para calcular 100 dígitos de π.
Sin embargo podemos hacer mucho mejor usando la misma función para calcular π si podríamos calcular arctan(1/x) más rápidamente.
Afortunadamente, Euler se subió con una manera de hacer precisamente eso:
arctan(1/x) = (x / (1 + x²)) + ((2 * x) / (3*(1+x²)²)) + ((2 * 4 * x) / (3*5*(1+x²)³)) + ((2 * 4 * 6 * x) / (3*5*7*(1+x²)⁴)) +...
El enésimo término en esta serie está dada por la función f (n):
f (n) = f(n-1) * (2 * n) / ((2*n+1)*(1+x²))
donde el primer término es (x / (1 + x²))
Con el fin de que el código funcione más rápido podemos calcular 1 + x² antes del bucle, así que sólo tenemos que calcular una vez. El código actualizado (gauss_pi_method_accelerated_arctan.py) es:
Este código tuvo poco más de medio segundo para calcular π a 10.000 dígitos, que es aproximadamente 30 veces más rápido que antes!
Hay un truco aseado para hacer el código ejecutar aún más rápido. Hasta ahora hemos estado usando la biblioteca decimal en Python. Sin embargo, si hacemos los cálculos con enteros todo será mucho más rápido. Para ello primero multiplicamos el valor inicial por una gran potencia de 10 y más tarde, cuando queremos utilizar el resultado, dividirlo por ese mismo poder de 10. Aquí es cómo el código se verá (gauss_pi_method_fixed_point.py):
Esta versión del código calcula 10.000 decimales de π en sólo 0,26 segundos, que es casi dos veces tan rápido como el método anterior!
