Prime Numbers

, par MiKaël Navarro

L’objectif de cet article est la mise au point d’un programme pour déterminer les n premiers nombres premiers. (merci Patrice pour le challenge ;)

Définition

Commençons par rappeler la définition d’un nombre premier :

Un nombre premier est un entier naturel qui admet exactement deux diviseurs distincts entiers et positifs : 1 et lui-même.

Cf. Nombre_premier

Ainsi :
 0 n’est pas premier car divisible par tout nombre non nul
 1 n’est pas premier car il n’a qu’un seul diviseur (lui-même)
 Au delà de 2 (qui est premier), un nombre premier est forçément impaire
 Autre curiosité, tous les nombres premiers (hormis 2 et 5) se terminent tous par 1, 3, 7 ou 9

Approche naïve

La première approche « naïve » pour déterminer la liste des nombres premiers <= n, est alors :
 une fois le nombre 2 identifié comme premier,
 de parcourir tout les nombres impaires,
 et de tester leurs primalité.

Test de primalité

L’approche « naïve » est alors de tester la primalité par divisions successives :

Pour un nombre num donné :
 si c’est 0 ou 1 ce n’est pas un nombre premier
 si c’est 2 c’est un nombre premier
 s’il est paire ce n’est pas un nombre premier
 sinon tester tous les diviseurs de 2 à num-1 :
 si l’on trouve un nombre qui divise num sans reste alors il n’est pas premier
 si on a parcouru toute la plage sans trouver un diviseur sans reste alors il est premier

On peut optimiser un peu le parcours en s’arrêtant à n/2, en effet si on n’a pas trouvé de diviseur avant n/2 on en trouvera pas au-delà. Mais la complexité restera en o(n) !

L’autre optimisation est de ne tester que les nombres de 2 à sqrt(n), en effet :

si num = pq alors soit p <= sqrt(n) soit q <= sqrt(n)

On passe ici à une complexité de o(sqrt(n)).

P.-S. on peut encore améliorer le parcours :
 On observe que tout nombre premier supérieur à 3 peut être écrit sous la forme 6k +/- 1 (avec k > 0).
 En effet tout entier peut être écrit sous la forme 6k + i (avec i = -1, 0, 1, 2, 3 ou 4).
 Comme 2 divise (6k + 0), (6k + 2) et (6k + 4) et 3 divise (6k + 3), il est plus efficient de tester les divisions par 2 et 3 puis tous les nombres de la forme 6k +/- 1. C’est 3x plus rapide !

À noter qu’en Python3 range est un itérateur [1].

prime_numbers-simple.png

Parallélisme

L’application de la méthode précédente permet d’obtenir les 5 761 455 premiers nombres premiers inférieurs à 100 000 000 en 21 minutes sur mon PC i7-8750H CPU @ 2.20GHz
Mais cet algo n’exploite qu’1 seule unité de calcul, mon PC ayant 12 cœurs l’idée est ici de distribuer les différents tests de primalité sur plusieurs cœurs.

La première approche est d’utiliser le multithreading. Pour rappel le multithreading consiste à exécuter plusieurs threads simultanément au sein d’un même processus.

La création de threads est un processus remarquablement léger, rapide et a une faible empreinte mémoire (ils partagent l’espace mémoire avec leur processus parent).

Malheureusement, en raison du GIL (ou Global Interpreter Lock) en CPython (l’implémentation la plus répendue), l’exécution du bytecode Python dans plusieurs threads est exécuté séquentiellement et non parallèlement.

Le GIL est une exclusion mutuelle (mutex) qui rend les choses plus saines au niveau des tâches permettant d’intégrer facilement des librairies externes qui ne sont pas à fil sécurisé, et exécute le code non-parallèle plus rapidement.

Il y a quelques exceptions :
 Les autres implémentations de Python tels que Jython ou IronPython ne sont pas concernés par le GIL
 Tout ce qui se passe en dehors du GIL (l’appel à des lib C via Cython par exemple)
 Les tâches qui font des I/O (ouvrir/écrire des fichiers, des requêtes http,…)
 Ou l’usage de la librairies comme NumPy.

Pour les autres cas, si on veut exécuter des tâches en parallèle en Python, il faudra se tourner vers le multitraitement.

Multiprocessing

Nous allons ici paralléliser la boucle for de test de primalité.

Paralléliser la boucle signifie répartir tous les processus en parallèle en utilisant plusieurs cœurs, chaque calcul n’attendant pas la fin du précédent.

Pour rappel un processus est une abstraction simple du système d’exploitation, c’est un ensenble d’instructions qui s’exécute indépendamment.
Ainsi, par rapport aux threads (appelé aussi processus léger), le lancement d’un nouveau processus demande l’allocations de ressources (I/O, ports, espace d’adressage) ce qui peut représenter un certain coût et notamment un temps au démarrage.

En Python le package multiprocessing prend en charge la création d’un processus enfant à la demande d’un autre processus en cours.
Assez similaire au module threading, il offre une interface quasi identique à celui-ci, mais privilégie les processus au lieu des threads.

On tombe à 6’ pour déterminer les 5 761 455 premiers nombres premiers.

prime_numbers-multi.png

Joblib

Un autre module, joblib utilise le multitraitement pour exécuter sur de multiples cœurs en parallèle un traitement dans une boucle for.
Il fournit un pipeline léger, simple et direct.

Il convient de mentionner que joblib n’est pas magique, le backend ’threading’ souffre du même goulot d’étranglement de GIL et le backend ’multiprocessing’ génère une surcharge importante en raison de la sérialisation de tous les paramètres et valeurs de retour.

Avec la joblib, les 5 761 455 premiers nombres premiers sont déterminé en 6’ sur 10 cœurs.

prime_numbers-joblib.png

Numba

Numba utilise un compilation JIT (Just-In-Time) à l’exécution : lors de la première exécution d’une fonction, Numba va générer un code optimisé pour la machine et c’est ce code qui va être utilisé lors des appels successifs de cette fonction.

Particulièrement efficace sur les fonctions, boucles et tableaux NumPy, Numba promet une exécution x100 avec seulement quelques lignes de code !

Appliqué à notre premier algo « naïf » :

Il réduit le temps de calcul à 1’36" au lieu des 21’ !!!

Crible d’Eratosthène

Le crible d’Ératosthène procède par élimination : il s’agit de supprimer d’une table des entiers de 2 à N tous les multiples d’un entier (autres que lui-même).

À la fin il ne restera que les entiers qui ne sont multiples d’aucun entier à part 1 et eux-mêmes, et qui sont donc les nombres premiers.

On commence par rayer les multiples de 2, puis les multiples de 3 restants, puis les multiples de 5 restants, et ainsi de suite en rayant à chaque fois tous les multiples du plus petit entier restant.

À la fin du processus, tous les entiers qui n’ont pas été rayés sont les nombres premiers inférieurs à N.

Très efficace, on trouve les 5 761 455 premiers nombres premiers en seulement 41".
Cependant, au-delà cet algo requirt beaucoup de mémoire (ici 2 Gio) !

prime_numbers-eratosthene.png

P.-S. quelques optimisations peuvent encore être faites :
 Une fois le nombre 2 identifié comme premier, on ne testera plus que les nombres impaires.
 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.
 Pour contourner le problème de mémoire, on peut créer un itérateur pour éviter de contenir en mémoire toute la liste des nombres premiers.
 On peut même utiliser un bitarray qui ne prend qu’un seul bit pour stocker chaque valeur booléenne.

Perfs

Un résumé ci-dessous des temps de calcul sur un Dell G5 5587 :
 i7-8750H CPU @ 2.20GHz x 12 ; 16Gio Mem
 Ubuntu 22.04 LTS / Python 3.10 (CPython)

méthode \ limite 1 000 000 jit 10 000 000 jit 100 000 000 jit
simple 1.67" 49" 5" 21’ 1’36"
multi(1000) 0.31" 6’ 2’
multi(10000) 0.28" 6’ 44"
multi(100000) 0.14" 11" 2" 6’ 19"
multi(1000000) 0.13" 6’ 14"
joblib(10) 3" 20" 18" 6’ 4’
Ératosthène 0.33" 4" 41" 5"

P.-S.

P.-S. le programme a été mis-à-jour pour tenir compte des dernières optimisations :

  • L’affichage des nombres premiers a aussi été tronqué pour ne pas impacter les perfs.
  • Le crible d’Ératosthène optimisé :
    • Avec l’itérateur et Numba trouve les 50 847 534 nombres premiers inférieurs à 1 000 000 000 en seulement 18" avec 2.65G de mémoire ;
    • Et les 455 052 511 nombres premiers inférieurs à 10 000 000 000 en 16’ au prix de 20G de mémoire (et en mettant en œuvre le swap) ! Au-delà ça core dump !
    • Dernière optimisation avec le bitarray (mais sans Numba) ça trouve les 455 052 511 nombres premiers inférieurs à 10 000 000 000 en 5’25" avec seulement 2.86G de mémoire !
    • Et les 4 118 054 813 nombres premiers inférieurs à 100 000 000 000 en 55’ avec 13G de mémoire ! Au-delà on a MemoryError !
  • Avec le multiprocessing et en réglant le chunksize à 10 000 000 on arrive à limiter l’utilisation mémoire (12x2G/500M) tout en exploitant les 12 cœurs :
    • Et permet de déterminer les 50 847 534 nombres premiers inférieurs à 1 000 000 000 en 5’25".
    • Plus lent que le crible d’Ératosthène, mais permet d’aller au-delà des fatidiques 10 000 000 000 en prenant son temps ;)