Fórmula Bailey-Borwein-Plouffe

La fórmula de Bailey-Borwein-Plouffe ( fórmula BBP ) es una fórmula para π . Fue descubierta en 1995 por Simon Plouffe y lleva el nombre de los autores del artículo en el que se publicó, David H. Bailey , Peter Borwein y Plouffe. ^[1] Antes de eso, Plouffe la había publicado en su propio sitio. ^[2] La fórmula es:

${\displaystyle \pi =\sum _{k=0}^{\infty }\left[{\frac {1}{16^{k}}}\left({\frac {4}{8k+1}}-{\frac {2}{8k+4}}-{\frac {1}{8k+5}}-{\frac {1}{8k+6}}\right)\right]}$

La fórmula BBP da lugar a un algoritmo de spigot para calcular el n -ésimo dígito de base 16 (hexadecimal) de π (y, por tanto, también el 4n -ésimo dígito binario de π ) sin calcular los dígitos anteriores. Esto no calcula el n- ésimo dígito decimal de π (es decir, en base 10). ^[3] Pero otra fórmula descubierta por Plouffe en 2022 permite extraer el n- ésimo dígito de π en decimal. ^[4] BBP y algoritmos inspirados en BBP se han utilizado en proyectos como PiHex ^[5] para calcular muchos dígitos de π mediante computación distribuida. La existencia de esta fórmula fue una sorpresa. Se creía ampliamente que calcular el n- ésimo dígito de π es tan difícil como calcular los primeros n dígitos. ^[1]

Desde su descubrimiento, se han utilizado fórmulas de la forma general:

${\displaystyle \alpha =\sum _{k=0}^{\infty }\left[{\frac {1}{b^{k}}}{\frac {p(k)}{q(k)}}\right]}$

Se han descubierto para muchos otros números irracionales , donde y son polinomios con coeficientes enteros y es una base entera . Las fórmulas de esta forma se conocen como fórmulas de tipo BBP . ^[6] Dado un número , no se conoce ningún algoritmo sistemático para encontrar , , y ; dichas fórmulas se descubren experimentalmente . ${\estilo de visualización \alpha}$${\estilo de visualización p(k)}$${\displaystyle q(k)}$${\displaystyle b\geq 2}$${\estilo de visualización \alpha}$${\estilo de visualización p(k)}$${\displaystyle q(k)}$${\estilo de visualización b}$

Una especialización de la fórmula general que ha producido muchos resultados es:

${\displaystyle P(s,b,m,A)=\sum _{k=0}^{\infty }\left[{\frac {1}{b^{k}}}\sum _{j=1}^{m}{\frac {a_{j}}{(mk+j)^{s}}}\right],}$

donde s , b y m son números enteros y es una secuencia de números enteros. La función P conduce a una notación compacta para algunas soluciones. Por ejemplo, la fórmula original de BBP: ${\displaystyle A=(a_{1},a_{2},\puntos ,a_{m})}$

${\displaystyle \pi =\sum _{k=0}^{\infty }\left[{\frac {1}{16^{k}}}\left({\frac {4}{8k+1}}-{\frac {2}{8k+4}}-{\frac {1}{8k+5}}-{\frac {1}{8k+6}}\right)\right]}$

se puede escribir como:

${\displaystyle \pi =P{\bigl (}1,16,8,(4,0,0,-2,-1,-1,0,0){\bigr )}.}$

Algunas de las fórmulas más simples de este tipo que eran bien conocidas antes de BBP y para las cuales la función P conduce a una notación compacta, son:

${\displaystyle {\begin{aligned}\ln {\frac {10}{9}}&={\frac {1}{10}}+{\frac {1}{200}}+{\frac {1}{3\ 000}}+{\frac {1}{40\,000}}+{\frac {1}{500\,000}}+\cdots \\&=\sum _{k=1}^{\infty }{\frac {1}{10^{k}\cdot k}}={\frac {1}{10}}\sum _{k=0}^{\infty }\left[{\frac {1}{10^{k}}}\left({\frac {1}{k+1}}\right)\right]\\&={\frac {1}{10}}P{\bigl (}1,10,1,(1){\bigr )},\end{alineado}}}$

${\displaystyle {\begin{aligned}\ln 2&={\frac {1}{2}}+{\frac {1}{2\cdot 2^{2}}}+{\frac {1}{3\cdot 2^{3}}}+{\frac {1}{4\cdot 2^{4}}}+{\frac {1}{5\cdot 2^{5}}}+\cdots \\&=\sum _{k=1}^{\infty }{\frac {1}{2^{k}\cdot k}}={\frac {1}{2}}\sum _{k=0}^{\infty }\left[{\frac {1}{2^{k}}}\left({\frac {1}{k+1}}\right)\right]\\&={\frac {1}{2}}P{\bigl (}1,2,1,(1){\bigr )}.\end{aligned}}}$

(De hecho, esta identidad es válida para a > 1:

${\displaystyle \ln {\frac {a}{a-1}}=\sum _{k=1}^{\infty }{\frac {1}{a^{k}\cdot k}}}$.)

Plouffe también se inspiró en la serie de potencias arctan de la forma (la notación P también puede generalizarse al caso en que b no es un entero):

${\displaystyle {\begin{aligned}\arctan {\frac {1}{b}}&={\frac {1}{b}}-{\frac {1}{b^{3}3}}+{\frac {1}{b^{5}5}}-{\frac {1}{b^{7}7}}+{\frac {1}{b^{9}9}}+\cdots \\&=\sum _{k=1}^{\infty }\left[{\frac {1}{b^{k}}}{\frac {\sin {\frac {k\pi }{2}}}{k}}\right]={\frac {1}{b}}\sum _{k=0}^{\infty }\left[{\frac {1}{b^{4k}}}\left({\frac {1}{4k+1}}+{\frac {-b^{-2}}{4k+3}}\right)\right]\\&={\frac {1}{b}}P\left(1,b^{4},4,\left(1,0,-b^{-2},0\right)\right).\end{aligned}}}$

Utilizando la función P mencionada anteriormente, la fórmula más simple conocida para π es para s  = 1, pero m  > 1. Muchas fórmulas ahora descubiertas son conocidas para b como un exponente de 2 o 3 y m como un exponente de 2 o algún otro valor rico en factores, pero donde varios de los términos de la secuencia A son cero. El descubrimiento de estas fórmulas implica una búsqueda por computadora para tales combinaciones lineales después de calcular las sumas individuales. El procedimiento de búsqueda consiste en elegir un rango de valores de parámetros para s , b y m , evaluar las sumas hasta muchos dígitos y luego usar un algoritmo de búsqueda de relaciones enteras (típicamente el algoritmo PSLQ de Helaman Ferguson ) para encontrar una secuencia A que sume esas sumas intermedias a una constante bien conocida o quizás a cero.

La fórmula original de suma de π de BBP fue encontrada en 1995 por Plouffe utilizando PSLQ . También se puede representar utilizando la función P :

${\displaystyle {\begin{aligned}\pi &=\sum _{k=0}^{\infty }\left[{\frac {1}{16^{k}}}\left({\frac {4}{8k+1}}-{\frac {2}{8k+4}}-{\frac {1}{8k+5}}-{\frac {1}{8k+6}}\right)\right]\\&=P{\bigl (}1,16,8,(4,0,0,-2,-1,-1,0,0){\bigr )},\end{aligned}}}$

lo que también se reduce a esta relación equivalente de dos polinomios:

${\displaystyle \pi =\sum _{k=0}^{\infty }\left[{\frac {1}{16^{k}}}\left({\frac {120k^{2}+151k+47}{512k^{4}+1024k^{3}+712k^{2}+194k+15}}\right)\right].}$

Se ha demostrado mediante una prueba bastante simple que esta fórmula es igual a π . ^[7]

Nos gustaría definir una fórmula que devuelva el dígito hexadecimal ( )-ésimo (con ) de π . Se requieren algunas manipulaciones para implementar un algoritmo de spigot utilizando esta fórmula. ${\displaystyle n+1}$${\displaystyle n\geq 0}$

Primero debemos reescribir la fórmula como:

${\displaystyle \pi =4\sum _{k=0}^{\infty }{\frac {1}{\left(16^{k}\right)(8k+1)}}-2\sum _{k=0}^{\infty }{\frac {1}{\left(16^{k}\right)(8k+4)}}-\sum _{k=0}^{\infty }{\frac {1}{\left(16^{k}\right)(8k+5)}}-\sum _{k=0}^{\infty }{\frac {1}{\left(16^{k}\right)(8k+6)}}.}$

Ahora, para un valor particular de n y tomando la primera suma, dividimos la suma hasta el infinito entre el término n :

${\displaystyle \sum _{k=0}^{\infty }{\frac {1}{\left(16^{k}\right)(8k+1)}}=\sum _{k=0}^{n}{\frac {1}{\left(16^{k}\right)(8k+1)}}+\sum _{k=n+1}^{\infty }{\frac {1}{\left(16^{k}\right)(8k+1)}}.}$

Ahora multiplicamos por 16 ⁿ , de modo que el punto hexadecimal (la división entre las partes fraccionarias y enteras del número) se desplaza (o permanece, si n = 0 ) a la izquierda del (n+1) -ésimo dígito fraccionario:

${\displaystyle \sum _{k=0}^{\infty }{\frac {16^{n-k}}{8k+1}}=\sum _{k=0}^{n}{\frac {16^{n-k}}{8k+1}}+\sum _{k=n+1}^{\infty }{\frac {16^{n-k}}{8k+1}}.}$

Como solo nos interesa la parte fraccionaria de la suma, observamos nuestros dos términos y nos damos cuenta de que solo la primera suma contiene términos con una parte entera; por el contrario, la segunda suma no contiene términos con una parte entera, ya que el numerador nunca puede ser mayor que el denominador para k  >  n . Por lo tanto, necesitamos un truco para eliminar las partes enteras, que no necesitamos, de los términos de la primera suma, para acelerar y aumentar la precisión de los cálculos. Ese truco es reducir módulo 8 k  + 1. Nuestra primera suma (de cuatro) para calcular la parte fraccionaria se convierte entonces en:

${\displaystyle \sum _{k=0}^{n}{\frac {16^{n-k}{\bmod {(}}8k+1)}{8k+1}}+\sum _{k=n+1}^{\infty }{\frac {16^{n-k}}{8k+1}}.}$

Observe cómo el operador módulo siempre garantiza que solo se conservarán las partes fraccionarias de los términos de la primera suma. Para calcular 16 ^{n − k}  mod (8 k  + 1) de forma rápida y eficiente, el algoritmo de exponenciación modular se realiza en el mismo nivel de bucle, no anidado . Cuando su producto acumulado 16 x se vuelve mayor que uno, se toma el módulo, al igual que para el total acumulado en cada suma.

Ahora, para completar el cálculo, se debe aplicar esto a cada una de las cuatro sumas por turno. Una vez hecho esto, las cuatro sumas se vuelven a colocar en la suma de π :

${\displaystyle 4\Sigma _{1}-2\Sigma _{2}-\Sigma _{3}-\Sigma _{4}.}$

Dado que solo la parte fraccionaria es precisa, para extraer el dígito deseado es necesario eliminar la parte entera de la suma final, multiplicarla por 16 y conservar la parte entera para "quitar" el dígito hexadecimal en la posición deseada (en teoría, los siguientes dígitos hasta la precisión de los cálculos utilizados también serían precisos).

Este proceso es similar a realizar una multiplicación larga , pero solo hay que realizar la suma de algunas columnas del medio. Si bien hay algunos acarreos que no se cuentan, las computadoras generalmente realizan cálculos aritméticos para muchos bits (32 o 64) y redondean, y solo nos interesan los dígitos más significativos. Existe la posibilidad de que un cálculo en particular sea similar a no sumar un número pequeño (por ejemplo, 1) al número 999999999999999, y que el error se propague al dígito más significativo.

Este algoritmo calcula π sin necesidad de tipos de datos personalizados que tengan miles o incluso millones de dígitos. El método calcula el dígito n- ésimo sin calcular los primeros n  − 1 dígitos y puede utilizar tipos de datos pequeños y eficientes. Fabrice Bellard encontró una variante de BBP, la fórmula de Bellard , que es más rápida.

Aunque la fórmula BBP puede calcular directamente el valor de cualquier dígito dado de π con menos esfuerzo computacional que las fórmulas que deben calcular todos los dígitos intermedios, BBP sigue siendo linealítmica ( ), por lo que valores sucesivamente mayores de n requieren cada vez más tiempo para calcularse; es decir, cuanto "más alejado" está un dígito, más tiempo le toma a BBP calcularlo, al igual que los algoritmos estándar de cálculo de π . ^[8]${\displaystyle O(n\log n)}$

DJ Broadhurst proporciona una generalización del algoritmo BBP que puede utilizarse para calcular una serie de otras constantes en tiempo casi lineal y espacio logarítmico. ^[9] Se dan resultados explícitos para la constante de Catalan , , , la constante de Apéry , , (donde es la función zeta de Riemann ), , , , y varios productos de potencias de y . Estos resultados se obtienen principalmente mediante el uso de escaleras de polilogaritmos . ${\displaystyle \pi ^{3}}$${\displaystyle \pi ^{4}}$ ${\displaystyle \zeta (3)}$${\displaystyle \zeta (5)}$${\displaystyle \zeta (x)}$${\displaystyle \log ^{3}2}$${\displaystyle \log ^{4}2}$${\displaystyle \log ^{5}2}$${\displaystyle \pi }$${\displaystyle \log 2}$

• Aproximaciones de π
• Matemáticas experimentales
• La fórmula de Bellard

• DJ Broadhurst, "Escaleras polilogarítmicas, series hipergeométricas y los diezmillonésimos dígitos de ζ(3) y ζ(5)", (1998) arXiv math.CA/9803067

• Richard J. Lipton , "Hacer de un algoritmo un algoritmo — BBP", publicación de blog, 14 de julio de 2010.
• Richard J. Lipton , "La clase de Cook contiene Pi", publicación del blog, 15 de marzo de 2009.
• Bailey, David H. "Un compendio de fórmulas de tipo BBP para constantes matemáticas, actualizado el 15 de agosto de 2017" (PDF) . Consultado el 31 de marzo de 2019 .
• David H. Bailey , "BBP Code Directory", página web con enlaces al código de Bailey que implementa el algoritmo BBP, 8 de septiembre de 2006.