Intel Acceler8: CPrimes

Envoyer un nouvel article

retweet
July 24, 2011 3:00 PM PDT


$2^{nd}$ Problème Acceler8 : cprime

GOOG_ICEANDMAGIC

22/05/2011

 

Introduction

Nous expliquons, par le présent, notre solution au problème cprime du concours acceler8. Vous allez remarqué que nous avons fait un effort de recherche pour chercher les implémentations les plus fiables et les plus rapides. Nous aurons peut dû faire nos propres implémentations, mais nous considérons qu'il ne faut pas réinventer la roue, et qu'il y a des gens qui accordent une importance capitale à ce problème. Autant exploiter leurs efforts.

Le travail de recherche fait équivaut à un vrai effort d'implémentation. Du moins, à notre, sens!

 

Définitions

Soient:

 

  • $P$ : l'ensemble des nombres premiers.
  • $C_k$ : l'ensemble des séquences de $k$ nombres premiers consécutifs. Nous appelons une telle séquence clique.
    $c = (p_1, p_2, ... ,p_k), c \in C_k \Leftrightarrow \forall i, 1 \leq i \leq k, p_i \in P, \forall np, p_i < np < p_i+1, not (np \in P)$

     

  • $S(c)$ : une fonction de $C_k$ dans $N$
    $C_k \quad \quad \quad \quad \quad \longrightarrow N$
    $ c = (p_1, ... ,p_k) \longrightarrow \sum_{i=1}^k p_i$

     

  • $P(c)$ : une fonction de $C_k$ dans $N$
    $C_k \quad \quad \quad \quad \quad \longrightarrow N$
    $ c = (p_1, ... ,p_k) \longrightarrow \prod_{i=1}^k p_i$

     

  • Couple de cliques : un pair $(c_1,c_2)$ de deux éléments de $C_k$. $(c_1,c_2) \in C_k \times C_k$.
  • $l$ est un entier désignant une borne inférieure.
  • $u$ est un entier désignant une borne supérieure.

     

 

Problématique

Étant quatre entiers, $l,u,\alpha,k$. Il faut trouver tous les couples $(c_1,c_2)$ dites $\alpha$-proches, $ c_1=(p_1, ..., p_k), c_2=(q_1, ..., q_k)$, de cliques vérifiant:

$ l \leq p_k, q_k \leq u, \vert S(c1) - P(c_2) \vert \leq \alpha $

Il est évident que le gros travail consiste à calculer les 'séquences' de nombres premiers entre $l$ et $u$. Cela avait été confirmé par mes statistiques. Le reste du rapport sera divisé en deux parties, une conscarée à la recherche des nombres premiers, l'autre au calcul des cliques $\alpha$-proches.

Nombres Premiers en Pratique

Mâme si notre effort de documentation n'est certainement pas à la hauteur des travaux faits, nous résumons ce dont nous étions capable de trouver et de comprendre.

 

Tests Probabilistes & Deterministes

Il est clair que, pour le moment, les tests probabilistes dominent les implémentations 'sérieuses'. À l'image de la très connue bibliothèque GMP. Ce genre de tests de primalité renvoie un booléen flou 1. En fonction d'une précision paramétrable, pour un entier donné, ce booléen peut avoir les valeurs suivantes:
  • 0 : l'entier est composée.
  • 2 : l'entier est premier.
  • 1 : l'entier est probablement premier.

La fiabilité de l'algorithme, et sa lourdeur, dépend essentiellement du degré de confiance qu'on cherche à atteindre. Cela est fait en calibrant le paramètre de précision.

Ce genre de tests peut âtre utilisé pour repérer, rapidement, les nombres premiers probables avant de faire un test, naif peut âtre, beaucoup plus rigoureux.

La littérature semble se congratuler pour l'invention d'un test de primalité, AKS en occurrence 2, le plus rapide (complexité polynomiale). Nous n'avons trouvé aucune implémentation sérieuse, et l'algorithme nous semble assez flou pour faire un effort d'implémentation.

 

La fenâtrisation : Wheel Factorization

l'Idée est d'éviter de balayer bâtement tous les nombres impair d'un espace de recherche, mais de considérer que les éléments qui sont suspectibles d'âtre premiers. Il s'agit de la fenâtrisation (wheel vectorization). John MOYER 3 nous explique que seuls les nombres suivants, étant un $N > 0$, sont suspectibles d'âtres premiers :

$N \times 30 + 1$
$N \times 30 + 7$
$N \times 30 + 11$
$N \times 30 + 13$
$N \times 30 + 17$
$N \times 30 + 19$
$N \times 30 + 23$
$N \times 30 + 29$

Appelons le $30$, dans $N \times 30$, fenâtre de recherche.

Évaluer les décalages $d$ (dans l'exemple: 1, 7, 11, ... etc) par rapport à $N \times 30$ consiste à retenir ceux qui ne font pas sortir un facteur commun avec la valeur de la fenâtre $30 = 2 \times 3 \times 5$, et donc d'éviter le: 2, 3, 4, 5, 6, 8, 9 ... Du coup, faire $8$ tests de primalités au lieu de $30$.

La largeur de la fenâtre peut âtre personnalisée, d'ailleurs, John MOYER 4 dans son implémentation choisit une largeur $2310 = 2 \times 3 \times 4 \times 5 \times 7 \times 11$. Ce qui lui permet de ne faire que $480$ tests au lieu de $2310$. Le choix de la largeur de fânetre est simple, le produit des premiers nombres premiers.

 

Les Cribles

Étant le coût très elevé des tests de primalités unitaires, il est plus adapté (dans notre cas), d'utiliser une méthode dite par .

La Crible d'd'Ératosthène

Consiste à éliminer tous les nombres divisibles par l'un des nombres premiers plus petits que lui. Les inconvignients sont multiples car:
  • Considérer tous les nombres d'un intervalle, ce qui peut âtre très gourment en mémoire.
  • Commencer par 2 quelque soit la borne inférieur de l'intervalle.

L'implémentation la plus aboutie de cette technique est PrimeSieve, écrite en C++, 5 de Kim WALISCH. Comparée à celle de MOYER, PrimeSieve est beaucoup plus rapide, et beaucoup plus simple à utiliser.

 

La Crible d'Atkin par Bernstein

Nous nous sommes pas aller plus loin dans les détails, mais la crible d'Atkin semble âtre une amélioration de la crible d'Ératosthène pour réduire sa consommation de mémoire et améliorer la convergence de l'algorithme.

La seule implémentation que nous ayons trouvé est celle de l'inventeur de la méthode, BERNSTEIN, primegen 6 écrite en C.

 

PrimeSieve vs primegen

Tous nos tests ont montré que PrimeSieve est sensiblement plus rapide que primegen. Nous l'avons retenu dans le release officiel de la compétition.

Présentation de la Solution

Discussion

Notre première intuition nous a orienté vers une approche où l'on ne calcule qu'une seule fois les nombres premiers. Cela revient à estimer une borne supérieure pour les nombres premiers, dits utiles (qui seront utilisé par nos calculs), puis les stocker dans une structure.

Cette dernière doit favoriser un accès direct et rapide aux contenu, doit être également scalable pour réagir d'une manière flexible au nombre grandissant de nombres premiers d'un espace trop grand. Nous avons pensé au arbres balancés et au STL map, mais nous nous sommes vite rendu compte que la structure n'est pas scalable et entraine un surcoût considérable.

Un bon compromis aurait être l'utilisation d'une table (STL vector) indexé par un arbre balancé (STL map). Verdict, même les tables dynamiques génèrent des exception liées à l'incapacité de stocker tous les nombres. Et puis, même les tables engendrent un surcoût assez considérable.

Nous nous sommes, finalement, restreints à ne précalculer que les nombres premiers qui appartiennent à l'espace de recherche donné par l'utilisateur, et puis calculer au besoin tous ceux qui ne le sont pas.

 

Borne liée à la capacité : BC

Il nous est demandé de représenter un nombre premier par un entier de 64bits. Étant donné la fonction $P$, et en considérant le dépassement de capacité, il parait clair que tous les éléments d'une clique $c$ doivent être inférieurs à $BC(k) = (2^{63}-1)^{1/k}$. Voici l'évolution de cette borne en fonction de $k$ :



$k$ $BC(k)$
2 $3.0370005 \times 10^9$
3 2097152
4 55108.9875
5 6208.37506
6 1448.15469
7 512
8 234.753035
9 128
10 78.7932425

L'ensemble des valeurs de $k$ qui permettent un calcul qui ne provoquera pas un dépassement de capacité est :
$ K = \{ k / BC(k) \leq \prod_{i=1}^k p_i , p_i \in P , p_1 \geq l \} $.
$l$ étant la borne inférieure de l'espace de recherche.

Considérant le meilleur cas : i.e. $p_1 = 2$, cet ensemble est très restreint $K = \{2, 3, 4, 5\}$.

Pour un $k$ donné, calculer des nombres premiers pour une clique en vue d'un calcul par la fonction $P$ provoquera un dépassement de capacité et donc des valeurs incorrectes.
$u_{reel} = min (BC(k), u_{utilisateur})$

ici $u_{utilisateur}$ correspond à ce que l'utilisateur donne comme $2^e$ argument au binaire.

 

Chercher les Nombres Premiers

Le pseudo-code de la recherche des nombres premiers comprisentre $l$ et $u$ est le suivant :

 

void	chercher_premiers(l,u, premiers){
Initialiser_recherche(l);
p = suivant();

for(;;){
   if( p > u)
         return;
   inserer(premiers, p);
   p=suivant();
  }

}

Il est, à noter, qu'il y a une dépendance de données inter-itérations forte sur la variable p. La seule issue de parallelisation, à notre sens, est de faire un blocking de la sorte :

 

void	chercher_premiers_bloc(l,u, premiers){
   // blocking
for(i=l; i + BLOCK_SIZE< u ; i+=BLOCK_SIZE)
    chercher_premiers(i,i+BLOCK_SIZE-1,premiers);
   
   // epilogue du blocking
chercher_premiers(i,i+BLOCK_SIZE-1,premiers);
}

Cette version maintient une dépendance sur premiers, mais elle n'est qu'une dépendance de sortie, donc, facile à enlever. La parallélisation de la boucle i est envisageable selon la rentabilité de cette transformation, qui dépend essentiellement de la valeur BLOCK_SIZE. Car nous nous payons les frais de :

  • u/BLOCK_SIZE initialisations.
  • le double calcul de u/BLOCK_SIZE nombres premiers.

 

Chercher les Cliques Proches

Une fois la liste des nombres premiers calculée, la recherche des cliques proches devient un calcul à coût dérisoire. Nos statistiques ont montré que cette partie ne consomme, dans tous les cas, que $20\%$ du temps global. Voici notre algorithme :

 

void    chercher_cliques_proches(premiers,k,alpha){
vecteur_circulaire 	c1, c2;


for(i=0; i < taille(premiers)-k+1; i++){

   // remplir c1 par les k premiers elements en commençant \`a partir de i
     remplir(c1,i,k,premiers);

      // estimer une borne inf\'erieure pour les \'el\'ements de la clique c2
    bi = estimer_debut_pour_c2(k, alpha, P(c1));
    bs = estimer_fin_pour_c2(k, alpha, P(c1));

      // calculer les nombres premiers compirs entre bi et bs
    chercher_premiers(bi,bs, autres_premiers);


      // remplir la table circulaire
    remplir(c2,0,k,autres_premiers);
    j=k;
    while( S(c2) + alpha < P(c1) ){
       if( abs( S(c2) - P(c1) ) <= alpha ){
           afficher(c1);
           afficher(c2);
           afficher(abs( S(c2) - P(c1) ));
           }
       j++;
       inserer(c2,autres_premiers[j]);
       }

    inserer(c1,premiers[i]);
    }
}

Les boucles itc1, itc2 sont complètement parallèles.

Il est à noter qu'il y a un risque de recalculer plusiers fois les mêmes nombres premiers (dans autres_premiers). Mais cela est très dérisoire car:

  • Pour le nombres petits, le calcul est assez rapide. Plus rapide même que le coût.
  • Il y a une très faible probabilité de recalcul pour les nombres assez grands.

Cependant, cet algorithme peut être amélioré dans ce sens là. Quand nous aurons le temps !

Mise en Oeuvre du Parallélisme & Résultats

Parallélisation

La collecte des nombres premier est éxecutée en parallèle pour tout intervalle de longeure supérieure à $2 \times 10^8$

La recherche des cliques proches est lancée en parallèle pour les très grands nombres ou un intervalle de recherche relativement large.

 

Quelques Chiffres

Tous les tests ont été réalisés 100 fois. Les performances medianes ont été retenues. La machine utilisée est une 8-multicoeurs Intel Xéon 3Ghz avec 8Mo de chache. La MTL étant trop chargée pour l'utiliser.

Pour un $\alpha = 1$, $ k = 2$ à partir de 1 :

bourne supérieure $u$ tout (us) Threads
$10$ 84 1
$10^2$ 1378 1
$10^3$ 1457 1
$10^4$ 1098 1
$10^5$ 85192 1
$10^6$ 9318984 5
$10^7$ 8s90 7

D'autres tests ont été faits, malheureusement, nous ne connaissons pas ceux qui sont les plus pertinents. Nous nous contentons de précédents pour se comparer aux chiffres annoncés sur le forum.

Conclusion

Nous avons la conviction que notre travail, élaboré dans la précipitation due au manque du temps, peut être grandement amélioré en utilisant des tests de primalités unitaire couplés à un mécanisme de fenêtrisation (wheel factorization) au lieu des cribles.

Mais cela ne garantit, malheureusement, pas la scalabilité de la solution. étant donné que même le test de primalité, prouvé le plus rapide (AKS), est polynomial. Ce qui n'arrange que peu notre problème.

 

Téléchargement


Footnotes

... flou1
Voir : Logique Floue
... occurrence2
M. Agrawal, N. Kayal, N. Saxena. Prime is in P
...Moyer 3
J. Moyer, http://www.rsok.com/ jrm/printprimes.html.
...Moyer 4
sieve2310 de John Moyer, ftp://ftp.rsok.com/pub/source/sieve2310.c
... C++,5
PrimeSieve de Kim Walisch, http://code.google.com/p/primesieve/
... primegen6
primegen de Bernstein, http://cr.yp.to/primegen.html