1 La formule BBP
2 Atteindre le -ième digit d’un nombre sans connaitre les précédents
References

1 La formule BBP

(1)

appelée depuis formule BBP (voir page Plouffe). Nos trois mathématiciens canadiens remarquent que cette expression est très proche d’une décomposition en base 16 de . Dans l’article [1] de 1996 devenu célèbre, ils montrent que l’on peut se servir de 1 pour atteindre le -ième digit de en base (a fortiori en base 16) sans avoir calculé les précédents, ceci en temps presque linéaire ) et en espace !! L’espace requis étant minime, ils appliquent immédiatement ce résultat en calculant le milliardième digit de en base 2.

Arrêtons-nous maintenant sur le principe de ce calcul du -ième digit de en considérant le cas général d’un nombre et d’une formule BBP associée:

2 Atteindre le -ième digit d’un nombre sans connaitre les précédents

Nous suivons ici à quelques adaptations près l’exposé remarquablement clair de Xavier Gourdon et Pascal Sebah sur le sujet [2] qui simplifiait déjà l’article original [1]. X. Gourdon est un grand connaisseur des algorithmes de calcul puisqu’il détient depuis plusieurs années le record de vitesse du calcul des décimales de avec le programme PiFast (plus rapide même que le programme de l’équipe de Kanada).

Pour simplifier le propos, on se restreint à la base , et le -ième bit de désigne le -ième bit de la partie fractionnaire de , c’est-à-dire après la virgule. La partie fractionnaire d’un nombre sera notée .

2.1 Décalage

La première idée repose sur le théorème suivant

Proof. La décomposition en base de s’écrit

(2)

d’où . Le -ième bit de est donc qui est bien le -ième bit de .  __

Example 2.2 Cherchons le 13ème bit de . Comme on a , en particulier, . Il reste donc à calculer le bit de . En passant en binaire, on a

(3)

donc le bit de est , et le bit est .

Chercher le -ième bit de correspond donc à trouver le premier bit de . D’après la formule BBP 1, cela correspond au premier bit de la série

(4)

où et sont les sommes

(5)

Compte tenu de la décroissance rapide de pour , il suffit d’évaluer seulement les premiers termes de (en virgule flottante !), plus exactement suffisamment pour que le reste de la somme soit inférieur à la précision requise (ici bit après la virgule).

Concernant , on peut mettre son terme général sous la forme

(6)

avec et dans . En écrivant la division euclidienne de par , on voit que comme , . Ceci s’écrit encore

(7)

2.2 Exponentiation

Imaginons ainsi que l’on veuille calculer le bit de . Grâce aux principes précédents, cela revient à calculer le bit de et donc le bit de . En vertu de 7, il faut alors trouver tel que donc évaluer .

Le calcul de s’effectue au moyen de la méthode appelée exponentiation binaire (“binary powering” ou “binary exponentiation”), dont l’idée semble avoir été utilisée dès 200 avant J.C [3], et même envisagée par les Egyptiens pour effectuer la multiplication ! Cela revient simplement à utiliser l’arithmétique modulo , autrement dit à se placer dans , en trois étapes :

1. On décompose en puissances de :
2. On calcule de proche en proche chaque puissance de jusqu’à :
3. On utilise les étapes 1 et 2 : .

Pour obtenir plus généralement , on initialise à 1 et on prend la plus grande puissance de inférieure à . Puis on applique l’algorithme suivant équivalent à notre exemple :

1. Si , alors ;  ; FinSi
2. 
3. Si alors ; Aller à l’étape 1. ; FinSi.

Enfin on termine par pour prendre en compte . En raison de la décomposition de en puissance de , on utilise opérations pour ce calcul, ce qui est très faible !

2.3 Etapes finales

Il reste à obtenir et à sommer le tout sur pour obtenir dans 5. Il faut ici remarquer que si le calcul de chaque terme de est effectué en virgule flottante et que, disons, le dernier bit est erroné, alors au pire la propagation des retenues sur toute la somme produira une erreur sur bits. Pour des valeurs de non extravagantes, on pourra donc mener les calculs de en toute sérénité avec une précision arithmétique de bits seulement.

Pour cette raison et en remarquant que pour , on somme termes nécessitant chacun opérations grâce à l’exponentiation, l’algorithme requiert opérations et en espace, où est la complexité de la multiplication de deux entiers de taille bits.

References