ChronoMath, une chronologie des MATHÉMATIQUES
à l'usage des professeurs de mathématiques, des étudiants et des élèves des lycées & collèges

Calcul de 100 (ou plus) décimales du nombre e = 2,7181828...

La précision des calculs des micro-ordinateurs est, dans les meilleurs des cas, de l'ordre de 10-16. Sans un programme approprié, même un "gros ordinateur" ne vous fournira qu'une soixantaine de décimales exactes (UNIVAC, CRAY, etc.). Nous allons exaucer un vieux rêve de mathématicien : calculer (au moins...) les 100 premières décimales du nombre e, base des logarithmes népériens.

Pour ce faire, nous devons trouver une formule stable, c'est à dire dont les erreurs d'arrondi ne se cumulent pas et un moyen de s'assurer à chaque calcul 100 décimales exactes.

Concernant la seconde condition, il nous suffira de construire une multiplication par blocs : nous simulons un calcul en base B = 10000 qui accélérera les calculs !

Les systèmes de numération :

Chaque "chiffre" de cette base sera un bloc (tranche de 4 chiffres décimaux) pouvant varier de 0 à 9999.

On pourrait prendre une base plus forte, mais on risque alors des erreurs dans le traitement des nombres entiers sur des ordinateurs ne possédant pas la double précision : si un décimal z est grand, il y aura une perte de décimales, donc de partie entière dans une instruction du type x = B*INT(z)

Nous allons voir que dans les deux cas qui nous intéressent, il suffit de définir la division à 100 décimales (cas de e) et la multiplication à 100 décimales (cas de π) en se limitant au cas où seul le dividende (resp. multiplicande) est un nombre à 100 décimales. Nous utilisons le développement en série obtenu par Euler :

En posant en = 1 + 1/2! + 1/3! + 1/4! + ... + 1/n! , on montre facilement en majorant le reste de la série, l'encadrement :

On déduit de cet encadrement que pour obtenir 100 décimales exactes du nombre e au moyen de en, il suffit de choisir n > 69. Pour ce faire, on peut écrire un petit programme permettant de calculer l'entier N tel que (ln désignant le logarithme népérien) :

Calcul du nombre de termes nécessaire :

Fournissez le nombre de décimales désiré nd, ce programme vous précisera jusqu'à quel rang Nmax vous devrez sommer. Pour 100, 300 et 1000 décimales, il faut Nmax au moins égal à (respectivement) : 70, 170, 450.




<SCRIPT LANGUAGE=JavaScript>
function calculn()
{nd=100 ; n=0
nd=prompt("Combien de decimales :",nd)
if (nd==null) {return} else nd=eval(nd)
s=0
while (s < nd*Math.log(10))
{n=n+10 ; s=0 ;
for (k=2;k<=n;k++){s=s+Math.log(k)}
s=s+Math.log((n+1)/(n+2))
}
alert("Ok pour Nmax="+n)
}
</SCRIPT>

 
Calcul pratique de e :

En se limitant à n = 4, on remarque que e4 = 1 + 1/1! + 1/2! + 1/3! + 1/4! peut s'écrire :

e4 = 1/4! + 1/3! + 1/2! + 1/1! + 1 = (((1/4 + 1) x 1/3 + 1) x 1/2 + 1) x 1/1 + 1

Plus généralement :

en = 1/n! + 1/(n - 1)! + 1/(n - 2)! + 1/(n - 3) + 1/3! + 1/2! + 1/1! + 1
    =
((((((((((1/n + 1) x 1/(n - 1) + 1) x 1/(n - 2) + 1) x 1/(n - 3) + 1) x ... x (1/4 + 1)         x 1/3 + 1) x 1/2 +1) x 1/1 + 1

D'où l'algorithme, très simple :

Condition initiale : e = 1 , n = 70

Programme :  Tant que n > 1
                
diviser e par n (division à 100 décimales)
                
ajouter 1 à e
                
décrémenter n (on soustrait 1 à n)
            
Fin tant que

            Ajouter 1 à e

La division à 100 décimales :

L'ajout de 1 et la décrémentation de n n'engendrent évidemment pas de dégradation de la précision des calculs. Le seul problème est la mise en place d'une division à 100 décimales :

Le principe de la division d'un nombre entier Z par un entier N est fort simple : si B est la base utilisée, Z est de la forme x.Bn+ y.Bn-1 + z.Bn-2 + ...

On effectue d'abord la division de x (première "tranche") par N :

x = Nqn + rn

d'où :

Z = (Nqn + rn)Bn + yBn-1 + ... = Nqn.Bn + (Brn + y)Bn-1 + ...

Le quotient commence donc par qn et on a regroupé rnBn + yBn-1 en (Brn + y)Bn-1 que l'on divise ("on abaisse la tranche suivante") de la même façon. Et ainsi de suite. En définitive, le quotient s'écrira alors :

qnBn + qn-1Bn-1 + ...

Notre dividende sera un nombre de 100 chiffres et le diviseur est l'entier naturel N. Nous procédons par tranches de 4 chiffres (base 10000); il nous faut donc 25 blocs.

Programme de division à 100 décimales :                    

T(0) contiendra la partie entière Q du quotient. le premier reste est R = T(0) - N*Q; nous le décalons dans la tranche suivante T(1) en le multipliant par B = 10000 et nous procédons à sa division par N; ainsi de suite jusqu'à la dernière tranche :

Pour i=0 à 25
  
Faire Q=Ent(Ti/N) , R=Ti - N*Q
       Ti devient Q
       Ti+1 devient Ti+1+B*R
Incrémenter i

Pour utiliser ce programme indépendamment, nombre de décimales au choix :

Cependant, nous voulons nous assurer 100 décimales exactes. Or comme dans toute opération, si l'erreur sur les opérandes est e, il n'est pas assuré que l'erreur sur le résultat soit inférieure ou égale e.

Montrons le sur l'exemple suivant : je veux diviser x = 1/3 par y = 1/7. Le résultat exact est 7/3 , soit 2,3333 avec 4 décimales exactes.

Si nous calculons les opérandes avec 4 décimales exactes, nous avons x = 0,3333 et y = 0,1428 et le quotient à 4 décimales est 2,3340. On constate qu'il y a erreur de l'ordinateur sur la 4ème décimale !

Voyons la théorie : soit x et y les approximations tronquées de deux nombres positifs x et y. Ce faisant, notons e l'erreur absolue sur x et y. Nous avons :

L'erreur absolue sur le quotient est donc inférieure à la valeur absolue de ce résultat et, concernant l'exemple précédent, l'erreur absolue est alors de l'ordre de 9,33.10-4, ce que confirme notre calcul.

L'algorithme utilisé pour le calcul de e est assez réconfortant quant à sa précision. En effet, le diviseur y varie de n à 2 et ce en nombres entiers, donc sans erreur. L'erreur absolue sur chaque quotient est alors , donc au plus égale à e/y, donc au plus égale à e/2.

 Conclusion : le programme de division, donné ci-dessus, ne nécessite pas de corrections pour le calcul de e.

Programmation de la méthode en JavaScript    pour 300 décimales de e
 

On doit avoir ici nmax au moins égal à 170; par prudence, on prend un peu plus... Pour 300 décimales, vu que 300 = 4 75, il nous faut 75 tranches, mais de par le report dans la tranche de rang i+1 dans l'algorithme de division, il nous faut définir une 76è tranche.

Une tranche de 3 chiffres, par exemple 378, signifie en fait 0378 : donc on ajoute des zéros, si nécessaire, en début de tranche convertie préalablement en chaîne.

 fonctions mathématiques usuelles | Traitement des chaînes de caractères



<SCRIPT LANGUAGE=JavaScript>
function go()
{
b=10000;nmax=180;T[0]=1;for(i=1;i<=76;i++) {T[i]=0}
for(n=nmax;n>=2;n--)
{
divise();T[0]++;
}
T[0]++;aff=T[0].toString()+",";
for(i=1;i<=75;i++)
{
t$=T[i].toString()
while (t$.length< 4) {t$="0"+t$}
aff=aff+" "+t$
}
alert(aff)
}

function divise()
{
for(i=0;i<=75;i++)
{
q=Math.floor(T[i]/n);r=T[i]%n
T[i]=q;
T[i+1]=T[i+1]+b*r
}}
</SCRIPT>

Programme équivalent en langage BASIC , Calcul élémentaire de e , Calcul de 1000 décimales de π


© Serge Mehl - www.chronomath.com