Pour la définition de ce qu'est un nombre premier, voir ici: http://fr.wikipedia.org/wiki/Nombre_premier

Les fonctions présentées ici renverront la liste des nombres premiers inférieurs à un nombre donné. Par exemple, pour ceux inférieurs à 100:

premiers(100) ⇒ [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97]

Combien y a-t-il de nombres premiers? Une infinité. Voilà quelques exemples:

nombres de nombres premiers inférieurs ou égaux à 10**n:
10**1  4
10**2  25
10**3  168
10**4  1229
10**5  9592
10**6  78498
10**7  664579
10**8  5761455
10**9  50847534
10**10 455052511  <= ok
10**11 4118054813
10**12 37607912018
10**13 346065536839
10**14 3204941750802
10**15 29844570422669
10**16 279238341033925
10**17 2623557157654233
10**18 24739954287740860
10**19 234057667276344607
10**20 2220819602560918840

Avec le code le plus rapide parmi ceux de cette page, on peut trouver les 455052511 nombres premiers inférieurs ou égaux à 10 milliards en une heure environ.

Au delà, il faudra une grosse mémoire et beaucoup de patience…

Le principe du calcul est indiqué ici: http://fr.wikipedia.org/wiki/Crible_d'%C3%89ratosth%C3%A8ne

Commençons par un code simple et très proche du principe de base:

def eratosthene(n):
    """retourne la tableau des nombres premiers <= n (crible d'eratosthene)"""
    n += 1
    tableau = [0,0] + [i for i in xrange(2, n)]
    for i in xrange(2, n):
        if tableau[i] != 0:
            # c'est un nombre 1er: on garde, mais on neutralise ses multiples
            for j in xrange(i*2, n, i):
                tableau[j] = 0
    return [p for p in tableau if p!=0]        

Le code se lit bien:

• on crée une liste de nombres comme [0,0,2,3,4,5,6,7,8,9,10,…]
• le nombre d'indice 1 est 0, parce que 1 n'est pas un nombre premier
• on progresse dans la liste, et à chaque fois qu'on trouve un nombre différent de 0: c'est un nombre premier, et on met à 0 tous ses multiples
• à la fin, on retourne le liste de tous les nombres qui ne sont pas à 0

C'est déjà un code rapide, mais on peut l'améliorer:

• une fois le nombre 2 identifié comme premier, on ne testera plus que les nombres impairs
• dans la boucle de test, une fois qu'on a atteint un nombre premier supérieur à racine de n, il n'y a plus de multiple à supprimer
• dans le tableau, on va plutôt contenir des booléens (True/False) que des nombres

Voilà un tel code amélioré:

def eratosthene(n):
    """retourne la liste des nombres premiers <= n (crible d'eratosthene)"""
    if n<2:
        return []
    n += 1
    tableau = [False,False] + [True]*(n-2)
    tableau[2::2] = [False]*((n-2)//2 + n%2) # sup. des nb pairs
    premiers = [2] # initialisation de la tableau des nb 1ers (2 est 1er)
    racine = int(n**0.5)
    racine = racine + [1,0][racine%2] # pour que racine soit impair
    for i in xrange(3, racine+1, 2):
        if tableau[i]:
            # on enregistre le nouveau nb 1er
            premiers.append(i)
            # on élimine i et ses multiples
            tableau[i::i] = [False]*((n-i)//i + int((n-i)%i>0)) 
    for i in xrange(racine, n, 2):
        if tableau[i]:
            # on enregistre le nouveau nb 1er (pas de multiples à supprimer)
            premiers.append(i)
    return premiers

Cette version est de plus très rapide. Par exemple, on trouve les 664579 nombres premiers inférieurs à 10 millions en une seule seconde!!! Et les 5761455 nombres premiers inférieurs à 100 millions en une dizaine de secondes.

Mais au delà, malheureusement, cette méthode consomme beaucoup de mémoire et génère une exception “MemoryError”.

Pour contourner ce problème de mémoire, commençons par créer un itérateur pour éviter d'avoir à contenir en mémoire toute la liste des nombres premiers:

def eratosthene_iter(n):
    """Itérateur retourne tous les nb premiers <= n (crible d'Eratosthene)"""
    if n<2:
        pass # il n'y a aucun nb 1er en dessous de 2!
    else:
        n += 1 # pour avoir les nb 1ers <=n et pas seulement <n
        tableau = [False, False] + [True]*(n-2)
        tableau[2::2] = [False]*((n-2)//2 + n%2) # supr. des nb pairs
        yield 2 # 2 est un nombre premier
        racine = int(n**0.5)
        racine = racine + [1,0][racine%2] # pour que racine soit impair
        i, fin, pas = 3, racine+1, 2
        while i<fin: # on ne traite que les nb impairs
            if tableau[i]:
                yield i # i est un nombre premier
                # on élimine i et ses multiples
                tableau[i::i] = [False]*((n-i)//i + int((n-i)%i>0)) 
            i += pas
        i, fin, pas = racine, n, 2
        while i<fin:  # on ne traite que les nb impairs
            if tableau[i]:
                yield i # i est un nombre premier
            i += pas

Comme c'est un itérateur, le plus naturel est de l'utiliser dans une boucle for:

for p in eratosthene_iter(100):
    print p

Mais on peut, à la rigueur, l'utiliser comme suit:

gen = eratosthene_iter(100)
while True:
    try:
        print gen.next()
    except StopIteration:
        break

Comme c'est un itérateur, on peut toujours retrouver la liste complète des nombres renvoyés:

premiers = [p for p in eratosthene_iter(1000)]

On peut même créer une fonction lambda:

premiers = lambda n: [p for p in eratosthene_iter(n)]

Ce code est très rapide, mais le tableau consomme encore beaucoup de mémoire.

Alors, on peut utiliser le module “bitarray' qui ne prend qu'un seul bit pour stocker chaque valeur booléenne.

On trouve ce module ici: http://pypi.python.org/pypi/bitarray. Sous Windows, il faut l'installer avec setuptools.

Voilà le nouveau et dernier code, toujours sous forme d'itérateur:

def eratosthene_ba_iter(n):
    """Itérateur retourne tous les nb premiers <= n (crible d'Eratosthene)
       on utilise ici le module 'bitarray' pour stocker les booléens
    """
    if n<2:
        pass # il n'y a aucun nb 1er en dessous de 2!
    else:
        n += 1 # pour avoir les nb 1ers <=n et pas seulement <n
        tableau = bitarray(n)
        tableau.setall(True)
        tableau[0], tableau[1] = False, False # 1 n'est pas un nb 1er
        yield 2  # 2 est un nombre premier
        tableau[2::2] = False # on élimine tous les nombres pairs
        racine = int(n**0.5)
        racine = racine + [1,0][racine%2] # pour que racine de n soit impair
        i, fin, pas = 3, racine+1, 2
        while i<fin: # on ne traite que les nb impairs
            if tableau[i]:
                yield i # on a trouvé un nouveau nb 1er: on le renvoie
                tableau[i::i] = False # on élimine i et ses multiples
            i += pas
        i, fin, pas = racine, n, 2
        while i<fin:  # on ne traite que les nb impairs
            if tableau[i]:
                yield i # on a trouvé un nouveau nb 1er: on le renvoie
            i += pas

Comme c'est un itérateur, on peut toujours retrouver la liste complète des nombres renvoyés:

premiers = [p for p in eratosthene_ba_iter(1000)]

On peut même créer une fonction lambda:

premiers = lambda n: [p for p in eratosthene_ba_iter(n)]

Cette dernière fonction est vraiment rapide et permet d'obtenir, par exemple, les 455052511 nombres premiers inférieurs à 10 milliard en une heure environ.

Il faut bien se rendre compte de ce que représentent une telle quantité de nombres: un fichier de 5,03 Go (5 403 267 048 octets)!

Le principe du calcul est assez simple. Pour le nombre n donné, on teste chacun des nombres de 2 à n pour savoir s'il est premier. Pour chaque test, le nombre est premier si on ne lui a pas trouvé de diviseur. Quand le nombre est premier, on le conserve dans une liste. A la fin, on renvoie cette liste.

Pour gagner du temps, on utilise plusieurs astuces:

- A part “2” qui est premier, on ne teste pas les nombres pairs qui ne sont forcément pas premiers.

- Pour chaque nombre k à tester, On arrête les divisions à racine de k parce que si on n'a pas encore trouvé de diviseur, on n'en trouvera plus après, et le nombre est donc premier.

- Comme on doit conserver tous les nombres premiers trouvés, on n'utilise qu'eux comme diviseur dans les tests de division!

Pour vérifier, vous pouvez voir ici pas mal de nombres premiers déjà calculés: http://nombrespremiersliste.free.fr/

Voilà le code proposé:

def premiers(n, p=[2,3,5]):
    """Retourne la liste des nombres premiers <= n (méthode=division)"""
    k = p[-1]+2
    if n < k:
        return [x for x in p if x<=n]
    while k <= n:
        i = 0
        while i < len(p):
            if p[i]*p[i] > k:
                p.append(k)
                break
            if (k % p[i]) == 0:
                break
            i += 1
        k += 2
    return p
 
# Exemple d'utilisation
print premiers(100)  # affiche [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97]

Performance: Cette fonction est très rapide. Par exemple, toujours avec l'accélérateur psyco, la fonction trouve les 5.761.455 nombres premiers inférieurs à 100 millions en une centaine de secondes. A noter que l'affichage de cette très longue liste sera plus longue que son calcul…

Avec cette fonction, on peut avoir n=1 milliard, et on trouve alors les 50.847.534 nombres premiers en à peu près 40 minutes. Mais il ne faut pas rêver: le temps de calcul grimpe très très vite, et ce n'est pas avec ça que vous allez casser le cryptage RSA: le dernier cryptage factorisé, le RSA-640, a demandé plus de 5 mois de calcul avec 30 CPU … Et on en est maintenant au RSA-2048 qui correspond à un nombre d'environ 600 chiffres!

Cas particulier: si n=1, la fonction renvoie une liste vide, parce que “1” n'est pas un nombre premier.

Vous pouvez bien entendu ajouter un test pour empêcher que la fonction ne soit appelée avec n<=0.

Vous pouvez tester cette fonction avec la Calculext ici: http://calculext.jpvweb.com, mais soyez raisonnable: au delà de premiers(100000), vous allez dépasser le temps maxi (8s) de calcul autorisé sur le serveur (qui, de plus, n'a pas psyco).