Jusqu’à présent, j’avais été généralement peu convaincu des textes francophones à propos des contenus sur Wikipedia, à quelques exceptions près. Très souvent, j’ai constaté que la Wikipedia anglophone est sensiblement plus détaillée que la Wikipedia francophone.
Par exemple, j’ai regardé http://fr.wikipedia.org/wiki/Test_de_primalit%C3%A9_de_Lucas-Lehmer_pour_les_nombres_de_Mersenne
Dans cette page, les explications de l’algorithme sont assez confuses. Cependant, j’ai trouvé un détail très intéressant dans l’article anglophone : http://en.wikipedia.org/wiki/Lucas%E2%80%93Lehmer_primality_test
Voici un pseudo-code :
// Determine if Mp = 2p − 1 is prime
Lucas–Lehmer(p)
var s = 4
var M = 2p − 1
repeat p − 2 times:
s = ((s × s) − 2) mod M
if s = 0 return PRIME else return COMPOSITE
C’est concis, simple, et propre.
A partir de ce pseudo-code, on peut adapter celui-ci à tous les langages informatiques.
Voici le code Perl que j’ai conçu à partir du pseudo-code :
#!/usr/bin/perl
for ($p = 3; $p <= 100; $p++)
{
$s = 4;
$M = (2 ** $p) – 1;
$k = $p – 2;
for ($i = 1; $i <= $k; $i++)
{
$s = (($s * $s) – 2) % $M;
}
if ($s == 0)
{
print « $M = 2 ^ $p – 1 is prime \n »;
}
}
Cela fonctionne, mais rapidement à partir d’un nombre assez grand, c’est « arrondi » à une expression sous forme de puissances en base 10. Et cela fausse le reste des calculs…
Dont acte :
7 = 2 ^ 3 – 1 is prime
31 = 2 ^ 5 – 1 is prime
127 = 2 ^ 7 – 1 is prime
8191 = 2 ^ 13 – 1 is prime
131071 = 2 ^ 17 – 1 is prime
524287 = 2 ^ 19 – 1 is prime
2147483647 = 2 ^ 31 – 1 is prime
7.20575940379279e+16 = 2 ^ 56 – 1 is prime
2.88230376151712e+17 = 2 ^ 58 – 1 is prime
1.15292150460685e+18 = 2 ^ 60 – 1 is prime
4.61168601842739e+18 = 2 ^ 62 – 1 is prime
Ici, au-delà de 2 puissance 31 – 1, tout est faux.
Mais j’ai une arme secrète : j’ai transposé le code Perl en code Python, et Python ici prouve sa toute-puissance : on a enfin des nombres complets, chiffre par chiffre, sans limite pour les nombres entiers !
Voici le code source en langage Python :
#!/usr/bin/python
p=3
while (p <= 100000):
s = 4
M = (2 ** p) – 1
k = p – 2
i=1
while (i <= k):
s = ((s * s) – 2) % M;
i=i+1
if (s == 0):
print « 2 ^ « , p, » – 1 is prime »
p=p+1
Attention cependant, il y a une indentation stricte à respecter dans le script en Python, il peut arriver qu’en collant le code ici l’identation ne soit pas apparente.
Voici les résultats du script Python :
2 ^ 3 – 1 is prime
2 ^ 5 – 1 is prime
2 ^ 7 – 1 is prime
2 ^ 13 – 1 is prime
2 ^ 17 – 1 is prime
2 ^ 19 – 1 is prime
2 ^ 31 – 1 is prime
2 ^ 61 – 1 is prime
2 ^ 89 – 1 is prime
2 ^ 107 – 1 is prime
2 ^ 127 – 1 is prime
2 ^ 521 – 1 is prime
2 ^ 607 – 1 is prime
2 ^ 1279 – 1 is prime
2 ^ 2203 – 1 is prime
2 ^ 2281 – 1 is prime
2 ^ 3217 – 1 is prime
2 ^ 4253 – 1 is prime
2 ^ 4423 – 1 is prime
La série ci-dessus est 100% exacte. 😉
Plus le nombre à tester est grand, plus le temps de calcul est long…
Pour que ça aille beaucoup plus vite :
- il faut un calcul informatique distribué, avec des milliers d’ordinateurs qui fonctionnent en commun
- il faut un ordinateur quantique (le temps de calcul deviendrait alors linéaire et non plus exponentiel)
Je vais bientôt établir le temps de calcul en fonction de l’exposant premier des nombres de Mersenne pour un PC de bureau classique. L’article ici présent sera réédité à cette occasion.
Réédition de l’article le 30 décembre 2014 :
J’ai modifié le programme Python en y ajoutant le module time, afin d’obtenir un outil qui servira à l’horodatage de chaque ligne de résultat.
#!/usr/bin/python
import time
last = time.time()
p=3
while (p <= 100000):
s = 4
M = (2 ** p) – 1
k = p – 2
i=1
while (i <= k):
s = ((s * s) – 2) % M;
i=i+1
if (s == 0):
now = time.time()
chrono = now – last
print « t = « , chrono, » s :: 2 ^ « , p, » – 1 is prime »
p=p+1
Voici le résultat, avec le temps écoulé depuis l’exécution du programme Python :
t = 0.0000138282775879 s :: 2 ^ 3 – 1 is prime
t = 0.00019097328186 s :: 2 ^ 5 – 1 is prime
t = 0.000243902206421 s :: 2 ^ 7 – 1 is prime
t = 0.000355958938599 s :: 2 ^ 13 – 1 is prime
t = 0.000600814819336 s :: 2 ^ 17 – 1 is prime
t = 0.000720024108887 s :: 2 ^ 19 – 1 is prime
t = 0.00141787528992 s :: 2 ^ 31 – 1 is prime
t = 0.00507187843323 s :: 2 ^ 61 – 1 is prime
t = 0.0113019943237 s :: 2 ^ 89 – 1 is prime
t = 0.0172729492188 s :: 2 ^ 107 – 1 is prime
t = 0.0260407924652 s :: 2 ^ 127 – 1 is prime
t = 0.673669815063 s :: 2 ^ 521 – 1 is prime
t = 1.09742283821 s :: 2 ^ 607 – 1 is prime
t = 15.6437869072 s :: 2 ^ 1279 – 1 is prime
t = 121.935223818 s :: 2 ^ 2203 – 1 is prime
t = 138.3559618 s :: 2 ^ 2281 – 1 is prime
t = 673.330420017 s :: 2 ^ 3217 – 1 is prime
t = 1962.14330792 s :: 2 ^ 4253 – 1 is prime
t = 2188.334692 s :: 2 ^ 4423 – 1 is prime
- A partir de ces données chronométriques, j’ai pu estimer que le temps de calcul nécessaire pour découvrir M48 qui est le 48e nombre premier de Mersenne est d’environ 300 années avec un PC de bureau ordinaire (AMD Athlon 64×2 dual core 4600+) si on lance le calcul au-delà de M47 tel que P est supérieur à 43112609 (en testant P jusqu’à 57885161) où M47 = 2^P – 1. Je ne peux donc pas me contenter d’un PC classique pour espérer découvrir de nouveaux nombres premiers de Mersenne encore inconnus à ce jour… Il faut donc au moins des centaines de microprocesseurs de PC pour pouvoir découvrir des nombres premiers de Mersenne dans des délais raisonnables.
En revanche, j’ai apporté une amélioration entre-temps : pour tout nombre premier de Mersenne (M), l’entier P est toujours un nombre premier impair (à l’exception de 2^2 – 1 = 3 où P est premier mais pair). J’ai donc ajouté un test de primalité de P.
#!/usr/bin/python
import time
last = time.time()
p=3
while (p <= 100000):
d=2
lim = int(p ** 0.5)
produit = 1.000
while (d <= lim):
produit = produit * (p % d)
if (d > 2):
d=d+2
if (d == 2):
d=d+1
if (produit != 0):
s = 4
M = (2 ** p) – 1
k = p – 2
i=1
while (i <= k):
s = ((s * s) – 2) % M;
i=i+1
if (s == 0):
now = time.time()
chrono = now – last
print « t = « , chrono, » s :: 2 ^ « , p, » – 1 is prime »
p=p+2
Cette fois, avec la nouvelle amélioration du script Python, le calcul est nettement plus rapide.
Résultats :
t = 6.48498535156e-05 s :: 2 ^ 3 – 1 is prime
t = 0.000236988067627 s :: 2 ^ 5 – 1 is prime
t = 0.000296831130981 s :: 2 ^ 7 – 1 is prime
t = 0.000382900238037 s :: 2 ^ 13 – 1 is prime
t = 0.000496864318848 s :: 2 ^ 17 – 1 is prime
t = 0.000589847564697 s :: 2 ^ 19 – 1 is prime
t = 0.000875949859619 s :: 2 ^ 31 – 1 is prime
t = 0.00200796127319 s :: 2 ^ 61 – 1 is prime
t = 0.0035548210144 s :: 2 ^ 89 – 1 is prime
t = 0.00507688522339 s :: 2 ^ 107 – 1 is prime
t = 0.00654602050781 s :: 2 ^ 127 – 1 is prime
t = 0.151687860489 s :: 2 ^ 521 – 1 is prime
t = 0.218498945236 s :: 2 ^ 607 – 1 is prime
t = 2.25717282295 s :: 2 ^ 1279 – 1 is prime
t = 15.902233839 s :: 2 ^ 2203 – 1 is prime
t = 18.2337968349 s :: 2 ^ 2281 – 1 is prime
t = 62.6072969437 s :: 2 ^ 3217 – 1 is prime
t = 181.48553586 s :: 2 ^ 4253 – 1 is prime
t = 206.73643899 s :: 2 ^ 4423 – 1 is prime
- On n’a plus besoin d’attendre 300 années pour déterminer M48 en testant tous les P supérieurs à 43112609, mais seulement un peu plus de 20 ans environ. Mais 20 ans c’est toujours long comme délai. Le calcul distribué est ainsi toujours une nécessité.
—
© 2014 John Philip C. Manson
jpcmanson 