Articulo de referencia

Fórmula de 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 recibe su nombre de los autores del artículo en el que s...

La fórmula de Bailey-Borwein-Plouffe ( fórmula BBP ) es una fórmula para π . Fue descubierta en 1995 por Simon Plouffe y recibe su nombre de los autores del artículo en el que se publicó: David H. Bailey , Peter Borwein y Plouffe. [ 1 ] La fórmula es:

π=k=0[116k(48k+128k+418k+518k+6)]{\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 en base 16 (hexadecimal) de π (y, por lo tanto, también el 4n -ésimo dígito binario de π ) sin calcular los dígitos precedentes. Esto no calcula el n -ésimo dígito decimal de π (es decir, en base 10). [ 2 ] En 2022, Plouffe publicó una fórmula que permite extraer el n -ésimo dígito de π en decimal, pero en la práctica requiere calcular los dígitos anteriores, ya que toma el n -ésimo dígito de una aproximación. [ 3 ] Los algoritmos BBP e inspirados en BBP se han utilizado en proyectos como PiHex [ 4 ] para calcular muchos dígitos de π utilizando computación distribuida . La existencia de esta fórmula fue una sorpresa porque 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 formulado fórmulas de la forma general:

α=k=0[1bkpag(k)q(k)]{\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.α{\displaystyle \alpha }, dóndepag(k){\displaystyle p(k)}yq(k){\displaystyle q(k)}son polinomios con coeficientes enteros yb2{\displaystyle b\geq 2}es una base entera . Las fórmulas de esta forma se conocen como fórmulas de tipo BBP . [ 5 ] Dado un númeroα{\displaystyle \alpha }No se conoce ningún algoritmo sistemático para encontrar el adecuado.pag(k){\displaystyle p(k)},q(k){\displaystyle q(k)}, yb{\displaystyle b}Estas fórmulas se descubren experimentalmente .

Especializaciones

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

PAG(s,b,metro,A)=k=0[1bkj=1metroaj(metrok+j)s],{\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, yA=(a1,a2,,ametro){\displaystyle A=(a_{1},a_{2},\dots ,a_{m})}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 BBP original:

π=k=0[116k(48k+128k+418k+518k+6)]{\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:

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

Fórmulas de tipo BBP previamente conocidas

Algunas de las fórmulas más sencillas 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:

ln109=110+1200+13 000+140000+1500000+=k=1110kk=110k=0[110k(1k+1)]=110PAG(1,10,1,(1)),{\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{aligned}}}
ln2=12+1222+1323+1424+1525+=k=112kk=12k=0[12k(1k+1)]=12PAG(1,2,1,(1)).{\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 se cumple para a > 1:

lnaa1=k=11akk{\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 de arcotangente de la forma (la notación P también se puede generalizar al caso en que b no es un número entero):

arctan1b=1b1b33+1b551b77+1b99+=k=1[1bkpecadokπ2k]=1bk=0[1b4k(14k+1+b24k+3)]=1bPAG(1,b4,4,(1,0,b2,0)).{\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}}}

La búsqueda de nuevas igualdades

Utilizando la función P mencionada anteriormente, la fórmula más simple conocida para π es para s  =  1, pero m  >  1. Se conocen muchas fórmulas ahora descubiertas para b como exponente de 2 o 3 y m como 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 computacional de 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 BBP para π

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

π=k=0[116k(48k+128k+418k+518k+6)]=PAG(1,16,8,(4,0,0,2,1,1,0,0)),{\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}}}

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

π=k=0[116k(120k2+151k+47512k4+1024k3+712k2+194k+15)].{\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 sencilla que esta fórmula es igual a π . [ 6 ]

Algoritmo de extracción de dígitos BBP para π

Nos gustaría definir una fórmula que devuelva el (norte+1{\displaystyle n+1})-o (connorte0{\displaystyle n\geq 0}) dígito hexadecimal de π . Se requieren algunas manipulaciones para implementar un algoritmo de espiga usando esta fórmula.

Primero debemos reescribir la fórmula como:

π=4k=01(16k)(8k+1)2k=01(16k)(8k+4)k=01(16k)(8k+5)k=01(16k)(8k+6).{\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 a través del n- ésimo término:

k=01(16k)(8k+1)=k=0norte1(16k)(8k+1)+k=norte+11(16k)(8k+1).{\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 n , de modo que el punto hexadecimal (la división entre la parte fraccionaria y la parte entera del número) se desplaza (o permanece, si n = 0 ) a la izquierda del (n+1) -ésimo dígito fraccionario:

k=016nortek8k+1=k=0norte16nortek8k+1+k=norte+116nortek8k+1.{\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, con el fin de acelerar y aumentar la precisión de los cálculos. Ese truco consiste en reducir módulo 8k + 1. Nuestra primera suma (de cuatro) para calcular la parte fraccionaria queda entonces:    

k=0norte16nortekmod(8k+1)8k+1+k=norte+116nortek8k+1.{\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 nk  mod  (8 k  +  1) de forma rápida y eficiente, el algoritmo de exponenciación modular se realiza en el mismo nivel de bucle, sin anidamiento . Cuando su producto 16 x acumulado es mayor que uno, se toma el módulo, al igual que para el total acumulado en cada suma.

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

4Σ12Σ2Σ3Σ4.{\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 "extraer" 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 requiere sumar algunas columnas centrales. Si bien existen algunos acarreos que no se contabilizan, las computadoras suelen realizar operaciones aritméticas con muchos bits (32 o 64) y redondean, y solo nos interesa el dígito o dígitos más significativos. Existe la posibilidad de que un cálculo 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.

BBP en comparación con otros métodos de cálculo de π

Este algoritmo calcula π sin necesidad de utilizar tipos de datos personalizados con miles o incluso millones de dígitos. El método calcula el enésimo dígito sin calcular los primeros n  1 dígitos y puede utilizar tipos de datos pequeños y eficientes. Fabrice Bellard descubrió 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 (O(norteregistronorte){\displaystyle O(n\log n)}), 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 tarda BBP en calcularlo, como ocurre con los algoritmos estándar de cálculo de π . [ 7 ]

Generalizaciones

DJ Broadhurst proporciona una generalización del algoritmo BBP que puede utilizarse para calcular varias otras constantes en tiempo casi lineal y espacio logarítmico. [ 8 ] Se dan resultados explícitos para la constante de Catalan ,π3{\displaystyle \pi ^{3}},π4{\displaystyle \pi ^{4}}, la constante de Apéryζ(3){\displaystyle \zeta (3)},ζ(5){\displaystyle \zeta (5)}, (dóndeζ(incógnita){\displaystyle \zeta (x)}es la función zeta de Riemann ),registro32{\displaystyle \log ^{3}2},registro42{\displaystyle \log ^{4}2},registro52{\displaystyle \log ^{5}2}y diversos productos de potencias deπ{\displaystyle \pi }yregistro2{\displaystyle \log 2}Estos resultados se obtienen principalmente mediante el uso de escaleras polilogarítmicas .

Véase también

Referencias

  1. 1 2 Bailey, David H.; Borwein, Peter B.; Plouffe, Simon (1997). "Sobre el cálculo rápido de varias constantes polilogarítmicas" . Matemáticas de la computación . 66 (218): 903– 913. doi : 10.1090/S0025-5718-97-00856-9 . hdl : 2060/19970009337 . MR 1415794 . 
  2. Gourdon, Xavier (12 de febrero de 2003). "Cálculo del enésimo dígito" (PDF) . Consultado el 4 de noviembre de 2020 .
  3. Plouffe, Simon (2022). "Una fórmula para el dígito decimal o binario $n^ { \rm th } $ de $π$ y potencias de $π$". arXiv : 2201.12601 [ math.NT ].
  4. "Créditos PiHex" . Centro de Matemáticas Experimentales y Constructivas . Universidad Simon Fraser. 21 de marzo de 1999. Archivado del original el 10 de junio de 2017. Consultado el 30 de marzo de 2018 .
  5. ^ Weisstein, Eric W. "Fórmula BBP" . MundoMatemático .
  6. Bailey, David H.; Borwein, Jonathan M.; Borwein, Peter B.; Plouffe, Simon (1997). "La búsqueda de pi". Mathematical Intelligencer . 19 (1): 50– 57. doi : 10.1007/BF03024340 . MR 1439159. S2CID 14318695 .  
  7. Bailey, David H. (8 de septiembre de 2006). "El algoritmo BBP para Pi" (PDF) . Recuperado el 17 de enero de 2013. Los tiempos de ejecución del algoritmo BBP... aumentan aproximadamente de forma lineal con la posición d .
  8. DJ Broadhurst, "Escaleras polilogarítmicas, series hipergeométricas y los dígitos diez millones de ζ(3) y ζ(5)" , (1998) arXiv math.CA/9803067

Lecturas adicionales

  • DJ Broadhurst, "Escaleras polilogarítmicas, series hipergeométricas y los dígitos diez millones de ζ(3) y ζ(5)" , (1998) arXiv math.CA/9803067
  • Richard J. Lipton , " Cómo hacer que un algoritmo sea un algoritmo — BBP ", entrada de blog, 14 de julio de 2010.
  • Richard J. Lipton , " La clase de Cook contiene Pi ", entrada de blog, 15 de marzo de 2009.
  • Bailey, David H. "Un compendio de fórmulas tipo BBP para constantes matemáticas, actualizado el 15 de agosto de 2017" (PDF) . Consultado el 31 de marzo de 2019 .
  • David H. Bailey , " Directorio de código BBP ", página web con enlaces al código de Bailey que implementa el algoritmo BBP, 8 de septiembre de 2006.