Nombres premiers: crible des multiplications
Le principe de l'algorithme consiste à parcourir la
table de multiplication pour neutraliser dans le crible les éléments qui correspondent à un produit.
Avantage de l'algorithme: temps d'exécution rapide et programmation simple (à part les dernières optimisations).
Désavantage: utilise un tableau (bool) de tous les nombres candidats. (voir: Optimiasation C)
Cet article correspond à un essai personnel (et peut-être original ?).
Connaitriez-vous d'autres articles similaires ? Ai-je réinventé la roue ?
Proposez-vous un autre nom que "crible des multiplications" ou "crible des produits" ?
C'est par hasard que j'ai remarqué dernièrement que l'étape
Mult4.js correspond à peu près à l'algorithme connu sous le nom de
Crible de Sundaram (fichier Summary.rtf).
Cet algorithme sera spécifiquement traité sous
Mult5.js (classique) et
Mult6.js (optimisé).
Une version actualisée de l'introduction aux nombres premiers, complétée de mesures des temps d'exécution, peut être obtenue ici
Nombres premiers.
Autres liens sur les nombres premiers:
de division (Contient la définition)
Essais de division avec liste
Crible d'Eratosthène
Crible d'Eratosthène (optimisation C)
Crible d'Atkin
Mult0.js
Voici pour commencer une interprétation très naïve de l'algorithme basé sur la table de multiplication:
// i 2 3 4 5 6 7 8 9 ...
// j
// 2 4 6 8 10 12 14 16 18
// 3 6 9 12 15 18 21 24 27
// 4 8 12 16 20 24 28 32 36
// 5 10 15 20 25 30 35 40 45
// 6 12 18 24 30 36 42 48 54
// 7 14 21 28 35 42 49 56 63
// 8 16 24 32 40 48 56 64 72
// 9 18 27 36 45 54 63 72 81
// ...
function Mult(mx) { // Mult0.js: naïf
var i,j,n,sv=[false,false]; // sv=new Array(mx+1)
for (n=2; n<=mx; n++) sv[n]=true; // 0,1 non premier
for (j=2; j<=mx; j++)
for (i=2; i<=mx; i++) if (i*j<=mx) sv[i*j]=false;
return sv;
}
Test direct
Ne mettez pas max > 10000, car ce code est très très lent.
Ne vaudrait-il pas mieux d'abandonner ici le développement de cet algorithme ?
Mult1.js
Malgré tout, nous allons continuer et y apporter plusieurs améliorations "évidentes":
a) La boucle extérieure peut finir avec j*j<=mx.
b) Nous pouvons enlever tous les "doublons" en commençant la boucle interne depuis la diagonale (i=j).
c) La boucle interne peut se terminer lorsque i*j>mx.
// i 2 3 4 5 6 7 8 9 ...
// j
// 2 4 6 8 10 12 14 16 18
// 3 9 12 15 18 21 24 27
// 4 16 20 24 28 32 36
// 5 25 30 35 40 45
// 6 36 42 48 54
// 7 49 56 63
// 8 64 72
// 9 81
// ...
function Mult(mx) { // Mult1.js:
var i,j,n,sv=[false,false]; // sv=new Array(mx+1)
for (n=2; n<=mx; n++) sv[n]=true;
for (j=2; j*j<=mx; j++)
for (i=j; i*j<=mx; i++) sv[i*j]=false;
return sv;
}
Test direct
Comme l'amélioration est fulgurante, nous continuons !
Mult2.js
La table de multiplication ne doit pas forcément être calculée avec des multiplications; on peut le faire avec des additions:
function Mult(mx) { // Mult2.js:
var j,n,sv=[false,false]; // sv=new Array(mx+1)
for (n=2; n<=mx; n++) sv[n]=true;
for (j=2; j*j<=mx; j++)
for (n=j*j; n<=mx; n+=j) sv[n]=false;
return sv;
}
Test direct
Mult3.js
Ne traitons la ligne j que lorsque j est premier et remplaçons le test j*j<=mx par j<=sqrt(mx):
function Mult(mx) { // Mult3.js:
var j,n,r=Math.floor(Math.sqrt(mx)),sv=[false,false]; // sv=new Array(mx+1)
for (n=2; n<=mx; n++) sv[n]=true;
for (j=2; j<=r; j++)
if (sv[j]) for (n=j*j; n<=mx; n+=j) sv[n]=false;
return sv;
}
Test direct
Mult4.js
Le prochain but est de ne consirérer que les
nombres impairs.
Comme pour "Sieve2.js", faisons correspondre à un élément du tableau sv non pas
l'index p lui-même, mais, pour p>0, la valeur impaire pp=1+2*p (donc p=(pp-1)/2):
// p: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 ...
// pp: 2 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 ...
//
// Adaptons la table de multiplication pp, et les indexes p correspondants:
// ii 3 5 7 9 11 13 15 17 19 ... i 1 2 3 4 5 6 7 8 9 ...
// jj j
// 3 9 15 21 27 33 39 45 51 57 ... 1 4 7 10 13 16 19 22 25 28 ...
// 5 25 35 45 55 65 75 85 95 2 12 17 22 27 32 37 42 47
// 7 49 63 77 91 105 119 133 3 24 31 38 45 52 59 66
// 9 81 99 117 135 153 171 4 40 49 58 67 76 85
// 11 121 143 165 187 209 5 60 71 82 93 104
// 13 169 195 221 247 6 84 97 110 123
// 15 225 255 285 7 112 127 142
// 17 289 223 8 144 161
// 19 361 9 180
// ... ...
On constate qu'on peut
directement "travailler" sur les indexes, ils augmentent de 1+2*j dans les lignes, et de 4*(1+j)=n dans la diagonale d:
function Mult(mx) { // Mult4.js:
var h=1+Math.floor((mx-1)/2),d,j,jj,p,sv=[]; // sv=new Array(h)
for (k=0; k<h; k++) sv[k]=true;
for (j=1,d=4; d<h; ++j,d+=4*j) // traitement par lignes
if (sv[j]) {
jj=1+2*j;
for (p=d; p<h; p+=jj) sv[p]=false;
}
return sv;
}
Test direct
Cette fonction semble être plus rapide que celles traitées sous
Nombres premiers: Crible d'Eratosthène.
Mult5.js
Comme annoncé au début, étudions l'algorithme du
Crible de Sundaram.
Soit In = 2*n + 1 un entier impair quelconque >= 3, avec n >= 1.
Soit Ip = In + 2*k un entier impair quelconque >= In, avec k >= 0.
On peut obtenir tous les produits de nombres impairs
In*Ip = In*(In+2*k) = (2*n+1)*(2*n+1+2*k)
en incrémentant n et k.
Utilisons comme précédemment:
// p: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 ...
// pp: 2 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 ...
//
// Voici le début du tableau obtenu pour pp=In*Ip, et les indexes p correspondants:
// nn 3 5 7 9 11 13 15 17 19 n 1 2 3 4 5 6 7 8 9
// kk k
// 3 9 1 4
// 5 15 25 2 7 12
// 7 21 35 49 3 10 17 24
// 9 27 45 63 81 4 13 22 31 40
// 11 33 55 77 99 121 5 16 27 38 49 60
// 13 39 65 91 117 143 169 6 19 32 45 28 71 84
// 15 45 75 105 135 165 195 225 7 22 37 52 67 82 97 112
// 17 51 85 119 153 187 221 255 289 8 25 42 59 76 93 110 127 144
// 19 57 95 133 171 209 247 285 323 361 9 28 47 66 85 104 123 142 161 180
On obtient les tableaux transposés décrits sous Mult4.
Sur le web, les "meilleurs" codes que j'ai trouvés correspondent à peu près à:
function Mult(mx) { // Mult5.js: Sundaram du web
var h=1+Math.floor((mx-1)/2),k,kk,n,nn,p,pp,sv=[]; // sv=new Array(h)
for (p=0; p<h; p++) sv[p]=true;
for (n=1; n<h; n++) { // traitement par colonnes
nn=2*n+1;
for (k=n; k<h; k++) {
kk=2*k+1;
pp=nn*kk;
p=(pp-1)/2;
if (p>=h) break;
sv[p]=false;
}
}
return sv;
}
Test direct Sundaram
Mult6.js
Limitons la boucle extérieure à nn*nn<=mx et traitons la ligne n que lorsque sv[n]==true (nn est premier).
function Mult(mx) {
var h=1+Math.floor((mx-1)/2),k,kk,n,nn,p,pp,sv=[]; // sv=new Array(h)
for (p=0; p<h; p++) sv[p]=true;
for (n=1,nn=3; nn*nn<=mx; n++,nn+=2) { // nn=2*n+1
if (sv[n])
for (k=n,kk=nn; k<h; k++,kk+=2) { // kk=2*k+1
p=(nn*kk-1)/2;
if (p>=h) break;
sv[p]=false;
}
}
return sv;
}
Puis simplifions les calculs de la boucle intérieure, ce qui revient à faire les calculs directement sur les indexes.
function Mult(mx) { // Mult6: Sundaram optimisé
var h=1+Math.floor((mx-1)/2),n,nn,p,sv=[]; // sv=new Array(h)
for (p=0; p<h; p++) sv[p]=true;
for (n=1,nn=3; nn*nn<=mx; n++,nn+=2) // nn=2*n+1
if (sv[n]) for (p=(nn*nn-1)/2; p<h; p+=nn) sv[p]=false;
return sv;
}
Test direct Sundaram optimisé
Mult7.js
Revenons à Mult4 et éliminons maintenant tous les nombres divisibles par 2 et 3, (voir aussi Cycle0.js sous Nombres premiers: Crible avec cycle additif)
A partir de p=5, on a p+=c[i%2], avec le "cycle" c=[2,4].
Adaptons également les tables de multiplication:
// n: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 ...
// p: 2 3 5 7 11 13 17 19 23 25 29 31 35 37 41 43 47 49 53 55 59 61 65 67 71 ...
//
// ii 5 7 11 13 17 19 23 25 29 31 35 37 41 43 47 ...
// jj
// 5 25 35 55 65 85 95 115 125 145 155 175 185 205 215 235
// 7 49 77 91 119 133 161 175 203 217 245 259 287 301 329
// 11 121 143 187 209 253 275 319 341 385 407 451 473 517
// 13 169 221 247 299 325 377 403 455 481 533 559 611
// 17 289 323 391 425 493 527 595 629 697 731 799
// 19 361 437 475 551 589 665 703 779 817 893
// 23 529 575 667 713 805 851 943 989 1081
// 25 625 725 775 875 925 1025 1075 1175
// 29 841 899 1015 1073 1189 1247 1363
// 31 961 1085 1147 1271 1333 1457
// 35 1225 1295 1435 1505 1645
// 37 1369 1517 1591 1739
// 41 1681 1763 1927
// 43 1849 2021
// 47 2209
// ...
//
// i 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 u
// j v [x,y]
// 2 9 12 19 22 29 32 39 42 49 52 59 62 69 72 79 8 [ 3, 7]
// 3 17 26 31 40 45 54 59 68 73 82 87 96 101 110 24 [ 5, 9]
// 4 41 48 63 70 85 92 107 114 129 136 151 158 173 16 [ 7,15]
// 5 57 74 83 100 109 126 135 152 161 178 187 204 40 [ 9,17]
// 6 97 108 131 142 165 176 199 210 233 244 267 24 [11,23]
// 7 121 146 159 184 197 222 235 260 273 298 56 [13,25]
// 8 177 192 223 238 269 284 315 330 361 32 [15,31]
// 9 209 242 259 292 309 342 359 392 72 [17,33]
// 10 281 300 339 358 397 416 455 40 [19,39]
// 11 321 362 383 424 445 486 88 [21,41]
// 12 409 432 479 502 549 48 [23,47]
// 13 457 506 531 580 104 [25,49]
// 14 561 588 643 56 [27,55]
// 15 617 674 120
// 16 737
// ...
Remarque: j'ai bien sûr fait un petit programme pour afficher ces deux tableaux.
Le nombre de candidats<=mx qui restent est est un peu difficile à calculer:
h = ((mx-5-a)/3+((a>1)? 4: 3)) , où a = (mx-5)%6.
Dans le tableau des indexes, après avoir noté sous u,v l'augmentation des éléments de la diagonale, observons que u,v augmentent alternativement de 8,16.
Puis, notons l'augmentation cyclique [x,y] des éléments dans les lignes: nous constatons que x=2*j-1 et que y vaut alternativement 4*j-1, 4*j-3.
function Mult(mx) { // Mult7.js: mx>9
var a=(mx-5)%6,h=((mx-5-a)/3+((a>1)?4:3)),d,j,k,sv=[],u,v,xy=[]; // sv=new Array(h)
for (k=0; k<h; k++) sv[k]=true;
for (j=2,d=9,u=0,v=8; d<h; d+=(++j&1)? u+=8: v+=16) {
if (sv[j]) {
if (j&1) {xy[0]=2*j-1; xy[1]=4*j-3;} else {xy[0]=4*j-1; xy[1]=2*j-1;}
for (k=d; k<h; k+=xy[k&1]) sv[k]=false; // k&1 équivalent à k%2
}
}
return sv;
}
Test direct
Une fois le crible calculé, on peut en extraire la liste des nombres premiers p (à partir de 5) avec la boucle
for (i=2,p=5; i<h; p+=cyc[i&1],i++) if (sv[i]) { /* p is prime */ }
Conclusions
Malgré un début très décevant, nous arrivons finalement à un code tout à fait "honnête", qui supporte la comparaison avec les "meilleurs".
Mesures
Le tableau ci-dessous n'a qu'une valeur comparative.
Les mesures dépendent de l'ordinateur, du navigateur et de leur charge.
Sur le même ordinateur, j'ai chaque fois effectué une série de tests et indiqué le meilleur résultat.
Voir aussi
Nombres premiers (fichier Summary.rtf)
Max: 1000 10000 100000 1000000 10000000 100000000
n: 168 1229 9592 78498 664579 5761455 memory:
Div0 2 63 3365 - - - 0
Div1 1 42 1770 145128 - - 0
Div2 0 2 35 309 7579 - 0
Div3 0 2 31 305 6799 - 0
Div4 0 1 30 296 6564 - 0
List0 0 1 26 323 7301 - int[n]
List1 0 1 19 172 2557 - int[n]
Sieve0 0 1 11 79 659 7104 bool[Max]
Sieve1 0 1 9 63 466 5003 bool[Max]
Sieve2 0 1 2 20 142 1495 bool[Max/2]
Sieve3 0 0 1 23 138 1455 bool[Max/2]
Atkin0 0 0 5 49 407 4847 bool[Max]
Atkin1 0 0 5 47 384 4777 bool[Max]
Atkin2 0 0 5 37 360 4154 bool[Max]
Mult0 12 456 140532 - - - bool[Max]
Mult1 0 0 4 48 787 10139 bool[Max]
Mult2 0 0 5 36 776 10091 bool[Max]
Mult3 0 0 3 25 226 2616 bool[Max]
Mult4 0 0 2 17 111 1307 bool[Max/2]
Mult5 Su 0 3 51 491 5356 bool[Max/2]
Mult6 Su-Opt 0 1 19 119 1363 bool[Max/2]
Mult7 0 0 1 8 77 850 bool[Max/3]
Testé avec:
Firefox 23.0 OK
Google Chrome 28.0 OK
Internet Explorer 10.0.7 OK
Opera 12.16 OK
Safari 5.1.7 OK
Vous n'êtes pas encore membre ?
inscrivez-vous, c'est gratuit et ça prend moins d'une minute !
Les membres obtiennent plus de réponses que les utilisateurs anonymes.
Le fait d'être membre vous permet d'avoir un suivi détaillé de vos demandes et codes sources.
Le fait d'être membre vous permet d'avoir des options supplémentaires.