Crible d'eratosthène optimisé

Description

Voici un calcul des nombres premiers compris entre 1 une limite donnée. Ce calcul utilise une version optimisée de l'algorithme tradionnel d'Eratosthène. Sur mon ordinateur, déjà un peu ancien, je calcule les 50847534 nombres premiers compris entre 1 et 1000000000 en 19.5 secondes. Ce logiciel est pévu pour aller à la limite des "unsigned int" du PC : (2^32)-1 Cependant pour des raisons de programmation expliquées la limite utilisable de l'intervalle de calcul est : 4294709602. Une documentation supplémentaire est incluse dans l'envoi. Et le source est largement commenté. L'optimisation porte sur la taille mémoire utilisée et sur le temps de calcul : les nombres pairs sont traités à part et les multiples de 3 ne sont pas repris ensuite. Il me semble que les nombreux autres programmes disponibles sur Internet qui utilisent l'algorithme d'Eratosthène n'ont pas une efficacité semblable. Sur cette base, il serait possible de développer une variante qui dépasserait la capacité limite des "unsigned int".

Source / Exemple :


/*  Crible d'Eratosthène par pgl10  */

#include <stdio.h>   // pour printf
#include <stdlib.h>  // pour system
#include <string.h>  // pour memset
#include <time.h>    // pour clock
#include <math.h>    // pour sqrt

typedef unsigned long ulong;

#define ulmax 4294967295  // (2^32)-1 : valeur maximum des unsigned long 

ulong nmax;    // borne supérieure de l'intervalle [1,nmax] de calcul utilisé
char* crible;  // crible[1+(nmax-1)/16] : une suite de (1+nmax)/2 bits utiles
//   quand le calcul sera fini les nombres pairs ne seront pas marqués
//   et le (1+i)/2-ième bit de crible mis à 1 quand i est premier et 0 si non.
//   crible[0] pour les impairs de 1 à 15, crible[1] pour les impairs de 17 à 31, etc.

void pause() {
    printf("\n");
    system("pause");    
}
void comp(ulong i) { // mettre à 0 le (1+i)/2-ième bit du crible ( i impair et composé )
    ulong q, r;
    char c[9]={'\xFF', '\xFE', '\xFD', '\xFB', '\xF7', '\xEF', '\xDF', '\xBF', '\x7F'};
//  if(i%2==0) return;         // inutile : tous les nombres pairs ne sont pas marqués
    q=i>>4;                    // q=i/16=(x-1)/8 avec x=(1+i)/2 ou x=1+i/2 pour i impair
    r=((1+i)>>1)-(q<<3);       // r=(1+i)/2-8*q=x-8*q avec x=(1+i)/2  ( r de 1 à 8 )
    crible[q]=crible[q]&c[r];  // on met à 0 le bit (1+i)/2 du crible quand i est impair
}
bool prem(ulong i) { // retourner true si le nombre i est premier et false si non
    ulong q, r;
    if(i<2 || i>nmax) return false;
    if(i==2) return true;
    if(i%2==0) return false;
    q=i>>4;                    // q=0 si i de 1 à 15, q=1 si i de 17 à 31, etc.
    r=((1+i)>>1)-(q<<3);
    if(crible[q]>>(r-1)&1) return true;  // si le bit (1+i)/2 de crible est à 1
    return false;
}
void Eratosthene(ulong n) {       // crible d'Eratosthène
    ulong i, j, m, q, r, s, t;
    char c[9]={'\xFF', '\xFE', '\xFD', '\xFB', '\xF7', '\xEF', '\xDF', '\xBF', '\x7F'};
    memset(crible, '\xFF', 1+(n-1)/16);      // initialisation de crible[]
    comp(1);                                 // 1 n'est pas premier
    for(j=9; j<=n; j=j+6) {                  // les multiples de 2 ne sont pas marqués
        q=j>>4;                              // pour marquer les multiples impairs de 3
        r=((1+j)>>1)-(q<<3);                 // on évite d'appeler comp(j) pour gagner
        crible[q]=crible[q]&c[r];            // un peu de temps             
    }
    m=1+(ulong)sqrt((double)n);
    s=4;
    for(i=5; i<m; i=i+s) {                   // on marque les multiples des autres 
        if(prem(i)) {                        // nombres premiers connus en évitant
            if(((i+2)*i)%3) t=4; else t=2;   // les multiples de 2 et 3
            for(j=i*i; j<=n ; j=j+t*i ) { 
                q=j>>4;                      // ici j est toujours impair
                r=((1+j)>>1)-(q<<3);
                if(crible[q]>>(r-1)&1) crible[q]=crible[q]&c[r];
                t=6-t;
//              if(j>ulmax-t*i) {  // pour éviter un dépassement de capacité des ulong
//                  printf("\nLa limite de l'intervalle de calcul est trop grande\n"); 
//                  pause();       // pour gagner du temps on peut éviter ce contrôle, 
//                  exit(1);       // mais il faut veiller à utiliser nmax < 4294709603
//              }
            }
        }
        s=6-s;
    }
}
int main (int argc, char *argv[]) { 
    ulong i, np=0, p1=0, pn;
    char limite[21];
    if(argc!=2) {
        printf("\nPour calculer les nombres premiers de 1 \205 nmax\n");
        printf("\nlancez ce programme en faisant : crible nmax\n");
        pause();
        return 1;
    }
    nmax=atol(argv[1]);
    sprintf(limite,"%lu",nmax);      // on doit y retrouver argv[1]
    printf("\nCalcul des nombres premiers entre 1 et %lu\n", nmax);
    crible=(char *)malloc((1+(nmax-1)/16)*sizeof(char));
    if(strcmp(argv[1],limite)!=0 || crible==NULL) {
        printf("\nLa limite de l'intervalle de calcul est trop grande\n");
        pause();
        return 1;
    }
    clock_t start = clock();
    Eratosthene(nmax);
    printf("\nTemps elapse du calcul : %f s\n", (double)(clock()-start)/CLOCKS_PER_SEC);
    for(i=1; i<=nmax; i=i+1) if(prem(i)){np++; if(p1==0)p1=i; pn=i;} // calcul de np, p1 et pn
    printf("\nNombre de nombres premiers calcul\202s : %lu\n", np);
    printf("\nLe premier nombre premier : %lu\n",p1);
    printf("\nLe %lu-i\212me nombre premier : %lu\n", np, pn);
    free(crible);
    pause();
    return 0;
}

Conclusion :


Je suis preneur de tous les commentaires et améliorations pouvant permettre d'effectuer ce calcul plus rapidement, tout en conservant l'algorithme d'Eratosthène.

Codes Sources

A voir également

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.