Asimov, Leibniz, pi, python, floats y evadirse de la realidad

La matemática tiene cosas sorprendentemente lindas. Y una de las que siempre me gustó de forma especial es la formula de Leibniz para el cálculo de π. También conocida como la serie de Gregory-Leibniz es una simpática manera de expresar el más famoso de los números irracionales:

O para aquellos no gustan de las notación sigma-grande:

Como verán, se trata de una serie sumamente simple y elegante. Y como suele ocurrir con estas cosas, aparecen en la vida de uno en momentos extraños, casi que al rescate.

Anoche me hallé frente a la biblioteca buscando nada en espacial. Necesitaba despejar la cabeza, olvidarme de los asuntos terrenales asociados con las complejas interacciones y relaciones humanas.

Fue en esa circunstancia que me encuentro con De los números y su historia de Isaac Asimov (que puede ser descargado desde aquí). Mientras lo hojeaba vi la formula de Leibniz, promediando el capítulo 6 e inmediatamente atrajo mi atención.

Después de expresar la serie, Asimov explica:

<<Ustedes podrán condenar mi falta de perseverancia, pero los invito a calcular la serie de Leibniz simplemente hasta donde la hemos escrito más arriba, es decir hasta 4/15. Incluso pueden enviarme una postal para darme el resultado. Si al terminar se sienten desilusionados al descubrir que su respuesta no está tan cerca de π como lo está el valor 355/113, no se den por vencidos. Sigan sumando términos. Sumen 4/17 al resultado anterior, luego resten 4/19, después sumen 4/21 y resten 4/23, etcétera. Pueden seguir hasta donde lo deseen, y si alguno de ustedes descubre cuántos términos se requieren para mejorar el valor 355/113, escríbanme unas líneas y no dejen de decírmelo. >>
En efecto, la fracción 355/113 parece ser la forma que mejor balancea precisión y simpleza a la hora de arrimarse a π desde lo números racionales, quedando a solo 2.66764189 × 10-7 de distancia. Mucha gente también utiliza la relación 22/7, aunque con menor precisión.

El desafío propuesto parecía interesante. Resultaba una entretenida oportunidad de jugar con floats en python, que siempre a resultado ser bastante tricky (además de ser un excusa para despejar la cabeza).

Empecé por escribir una versión relativamente obvia de la solución, que me dijo que al llegar al cociente 7497257 (n=3748629 en la sumatoria de Leibniz) el valor de la sumatoria estaría en 3.14159292035, quedando a 2.66764081935 x 10-7 de π. Dado que los decimales (floats, representaciones de punto flotante) en Python no se comportan de forma intuitiva, decidí hacer una segunda versión.

En esta versión traté de mudarme al predecible mundo de los enteros. Y obtuve valores levemente distintos. En el cociente 7497247 (n=3748624) se obtiene la primera aproximación mejor que 355/113, que es 3.141592386825... (a 2.66764162 x 10-7 de π)

Tampoco estoy seguro de que la segunda solución sea la correcta. Me pregunto cual será el promedio de valores en las postales que la viuda Asimov seguramente guarda en una caja de zapatos.

El hecho es que fue entretenido intentarlo y tal vez continúe con la experiencia la próxima vez que quiera desconectarme del mundo. Después de todo, como dijo el mismo Asimov: “Los hombres que se acostumbran a preocuparse por las necesidades de unas máquinas, se vuelven insensibles respecto a las necesidades de los hombres”, y hay veces que volverse insensible se ve, erróneamente, como una propuesta seductora. Hay gente que ahoga penas en alcohol, mucho más sano es ahogarlas matemáticas...

BTW, ya que tenía los datos armé un plot de los primeros 1000 valores de la serie para mostrar de forma gráfica la convergencia.

Ahora es momento de volver a la realidad, que buena o mala, irracional o racional (seguro que real al menos), inexacta o precisa, después de este recreo metal ya no se ve tan mal :-)

UPDATE: Sat, 20 Sep 2008 22:18:18 -0300: La solución estaba realmente cerca (en mi misma ciudad) y se llama Facundo, quien dejó un comentario. Bah.. se llama decimal, pero existe gracias a Facundo :). La biblioteca decimal permite manejar de forma exacta una cantidad arbitraria de decimales. La respuesta correcta es........ (redoblantes de suspenso) con n=(7497259+1)/2! aproxima a Π en 3.141592386825668744985771256 a solo 2.66764125 × 10-7 del número irracional. El nuevo script.. aquí.

Gracias Facundo, gracias decimal, gracias Asimov, gracias Leibniz....

10 thoughts on “Asimov, Leibniz, pi, python, floats y evadirse de la realidad”

  1. Lo único que me falta es aprender Python ;-) Mas que interesante y divertido el tema :-)

    By the way, ya pegué el swirl que me regalaste en DebConf. Ni bien pueda [conseguir alguien que me preste una cámara digital] le saco una foto y hago un blogpost :-D

    ¡Gracias!

  2. Estimado, algunas veces se obtiene mejor precisión en la aproximación al sumar la serie en sentido contrario, es decir, arrancar sumando el n-esimo término mas el n-esimo - 1, ... hasta el primer término. Podrías realizar este test y comprobar el comportamiento del error relativo.
    Por otro lado, la convergencia de esta serie es lentísima, una alternativa para acelarla enormemente es utilizar la transformación de shanks, pero eso es otro tema... no menos interesante.

    Por cierto, muy interesante tu post.
    Saludos.

  3. Buen punto en el de la cuenta regresiva. Los resultados modificando el segundo script:
    n=(7497247+1)/2 -> 3.141592386825242 (la distancia a pi es de 0.000000266764551)

    Así que sí.. algo mal hay.

    El script aquí.

    Gracias por la sugerencia.

    Lo lenta de la convergencia es algo que Asimov quería enfatizar. De ahí el desafío. De cualquier forma.. no pienso pasar mi vida calculando decimales... Me pregunto en que hubiese ocupado su tiempo libre Shanks si hubiese nacido un siglo y medio más tarde :).

  4. Vos sabes que me voy a comprar este libro de Asimov, a ver si en SFE lo consigo hoy mismo...

    Con respecto a lo de Shanks, estaría con Cooley y Tukey, viendo como acelerar el booteo de gnu/linux...
    Saludos.

    Emiliano

  5. No se mucho de programacion, por lo que no revise tus scripts, pero se me ocurre que si en python hay problemas con el punto flotante, pues simplemente bastaria multiplicar por una potencia de 10 segun la cantidad de decimales que quieras probar. Si esto fue lo que hiciste en el segundo intento, entonces no me pongas atencion :P.

  6. Es básicamente la idea del segundo script. De cualquier forma no puedo vivir el en mundo de los naturales totalmente, ya que la serie implica una fracción que no puedo asegurar que dé entera. De cualquier forma, gracias por el comentario :)

Comments are closed.