Hack en C
Cet article est consacré à ce que certains appellent le hack. Il ne s'agit pas pénétrer une n-ième fois dans le poste de travail d'une secrétaire du Pentagone mais de tirer le maximum de performance des processeurs actuels.
Ces processeurs sont devenus très complexes. Notamment le core P6 (le c
ur du processeur lui-même) à la base des Pentium Pro, des PII avec le MMX et des PIII avec le MMX et le SSE (MMX pour nombres flottants), et des Celeron. Le Pentium 4 utilise un nouveau core encore plus spécial. Ces conseils s'appliquent également pour les autres processeurs comme l'Alpha ou le PowerPC.
Tout le monde utilise aujourd'hui gcc. Malheureusement, il n'est plus vraiment adapté pour garantir un maximum de performance. Par exemple, il ne sait pas encore gérer les instructions comme le MMX et le SSE, qui utilisent un format de données spécial (des paquets de nombres alignés en mémoire). Mais les structures de ces processeurs font que certaines manières de programmer permettent une accélération énorme des performances. Pour qu'un compilateur accélère le code, il faudrait qu'il comprenne ce que l'on cherche à faire (la sémantique). C'est ce que l'on appelle de la réécriture et c'est encore au stade de la recherche. Pour l'instant, il faut y aller à la main.
De ce fait, le code devient incroyablement complexe et tordu et donc parfaitement illisible. Il y a un certain nombre de règles à respecter pour que gcc s'en sorte le mieux possible. Il ne faut utiliser les autres 'trucs' qu'en dernier recours pour garder le code maintenable au maximum et compréhensible, voire portable, si on utilise du code assembleur.
Il ne faut pas oublier non plus que l'algorithme retenu détermine les performances de l'application en premier lieu. Un meilleur codage peut faire gagner un facteur 10, un meilleur algorithme un facteur 1000.
Options de gcc
Il faut déjà que gcc soit poussé dans ses retranchements avant d'aller plus avant. Voici toutes les options trouvées dans le manuel de GCC :
gcc -mcpu=i686 -O3 -g -o t test.c -fstrength-reduce -frerun-loop-opt -fexpensive-optimizations -fschedule-insns2 -funroll-loops -fomit-frame-pointer -malign-double -fno-strict-aliasing -pipe -malign-loops=2 -jumps=2 -malign-functions=2 -DCPU=686 -ffast-math
Il s'agit des options qui ont l'air le plus efficace. Toutes les explications peuvent être trouvées sur : http://gcc.gnu.org/onlinedocs/gcc-2.95.3/gcc_2.html.
-funroll-loops est utilisé pour une histoire de rendement d'instructions, assez simple à comprendre. Imaginons une bête boucle for :
label: add
i++
comparaison avec n
jump label
Donc, une instruction sur quatre fait un calcul, le reste est seulement de la simple gestion de boucle. On peut donc voir un rendement de 25%... Si on a :
label: add
add
add
add
add
add
add
add
i+=8
comparaison
jump label
Donc, 8 sur 11, cela donne ainsi un rendement de 72 %. Mais il faut pouvoir dérouler la boucle. Il semble que gcc aime bien la dérouler 16 fois.
Tous les "align" de la ligne de compilation servent à aligner les données sur la largeur du bus du processeur. Intel précise qu'un mauvais alignement peut dans le pire des cas créer une pénalité de plus de 100 cycles ! Le core P6 dispose de 3 décodeurs. Ils reçoivent ainsi le flot d'instructions par paquet. Lors d'une boucle, si l'adresse de saut n'est pas alignée, au pire, un seul décodeur peut servir pour décoder l'instruction suivante. D'où l'intérêt de l'alignement pour garantir que les 3 décodeurs fonctionnent en même temps le plus souvent.
Il existe le même problème pour les données. En effet, le bus faisant 64 bits de large, lors d'accès sur des données plus petites, le processeur doit faire un décalage après un load pour récupérer les bonnes données. C'est pour cela que l'on force l'alignement pour éviter ce décalage de bits qui fait perdre du temps.
-pipe ne sert qu'à accélérer la compilation en n'écrivant pas de fichiers temporaires sur le disque, ceci à condition que vous ayez suffisamment de mémoire.
Pour le reste, allez faire un tour du côté du manuel (en bon anglais, RTFM ;-)!
Pour être sûr des augmentations de performance, il faut créer un cadre strict et bien défini.
Après avoir effectué une suite de tests, on s'aperçoit vite que les temps mesurés pour un même programme, sur une même machine, changent beaucoup, surtout s'ils sont courts. En fait, ces variations sont dues, le plus souvent, au temps passé dans le noyau. Toutes les 10 ms, l'ordonnanceur est lancé. Donc, si une performance tourne autour de 10 ms, le calcul peut se finir avant ou après le passage dans l'ordonnanceur.
Dans le temps mesuré pour une fonction donnée, il n'y a pas que le temps passé dans le code du noyau, mais aussi ce que l'on appelle la pollution du cache. C'est-à-dire que les données contenues dans les caches sont remplacées par celles du noyau au détriment de celles de l'application, d'où une perte de temps qui, de plus, est assez aléatoire.
En fait, il existe deux schémas d'ordonnancement dans le noyau Linux qui sont dits "temps réel", SCHED_FIFO et SCHED_RR. Le modèle habituel est le SCHED_OTHER. Il y a une notion de priorité statique mais qui ne nous intéresse pas ici (man sched_setscheduler pour en savoir plus). Pour faire simple, un processus SCHED_FIFO garde toujours la main sauf pour la rendre expressément (par une entrée/sortie, par l'arrivée d'un processus de priorité statique supérieur ou par sched_yield() ). SCHED_RR est un peu plus gentil puisqu'il partage le temps entre les processus de même priorité.
Donc, SCHED_FIFO nous garantit que l'ordonnanceur n'est jamais appelé. Évidemment, avec une boucle infinie et ce genre de priorité, vous perdez définitivement la main. On peut s'en sortir en donnant la même priorité au shell que vous utilisez. Une autre solution consiste à faire de temps en temps des entrées-sorties pour rendre un peu la main au système. Le bon vieux printf() fera l'affaire. Mais je vous aurai prévenu ! Il n'est plus question de ctrl-alt-F2 pour retrouver une petite console.
#include<sched.h>
const struct sched_param *schedParam;
sched_setscheduler(0,SCHED_FIFO,schedParam);
0 désigne le processus courant. On peut utiliser le PID pour en désigner un autre. SCHED_FIFO désigne le style d'ordonnancement voulu. schedParam renvoie les paramètres effectivement utilisés par le processus.
Pour des raisons de sécurité, tout cela ne peut s'exécuter qu'en tant que root.
Pour garantir que la mémoire de l'application ne sera pas "swappée" sur le disque pour faire de la place, on peut utiliser :
#include <sys/mman.h>
mlockall(MCL_CURRENT); /* to avoid swaping page to disk */
munlockall();
Ces techniques permettent de vérifier que le code s'exécute un poil plus vite qu'un autre. Mais il ne faut pas oublier les applications réelles où les processus ne devraient jamais tourner avec l'identité root.
Mesure du temps
Bon d'accord, pour mesurer un temps, on pourrait commencer par faire plus simple. Par exemple, on peut boucler plusieurs fois autour de la même fonction pour en mesurer le temps effectif. Mais en fait, on ne mesure plus tout à fait la même chose, encore à cause de la mémoire cache !
En fait, à la première exécution, on remplit ce cache que l'on utilise par la suite. Il semble que le temps d'exécution se stabilise après 3 passes. Cela dépend également de l'application. Dans tous les cas, c'est une méthode à oublier car c'est trop loin de la réalité !
Il faut pouvoir mesurer le temps de la façon la plus précise. L'outil de base est :
#include <sys/times.h>
struct tms timesApres;
times(×Apres);
Cette fonction renvoie donc une information sur le temps écoulé depuis le début d'exécution du processus. Pour plus d'information, faites un man 2 times.
struct tms
clock_t tms_utime; /* durée utilisateur */
clock_t tms_stime; /* durée système */
clock_t tms_cutime; /* durée utilisateur des fils */
clock_t tms_cstime; /* durée système des fils */
;
On peut aussi utiliser la fonction time du bash. Par exemple, "time cat toto.txt" rend le temps écoulé à lister l'hypothétique fichier toto.txt.
Ces fonctions ont l'avantage de séparer le temps passé dans le noyau et celui passé dans l'application. Il n'y a pas que l'ordonnanceur qui prend du temps, un malloc ou un printf en prenne également. Il faut d'ailleurs savoir que le temps d'un malloc est quasiment le même quelle que soit la taille de mémoire demandée (avec un minimum de 16 octets). Donc, à utiliser le moins souvent possible !
Mais en fait, ces fonctions sont inutiles la plupart du temps. La résolution temporelle de ces fonctions est bloquée à 10 ms. Impossible de descendre en dessous !
Avec une résolution de 10 ms et l'impossibilité de lancer plusieurs fois la même fonction, est-ce donc la fin du monde ? Non ! On va directement attaquer le compteur interne du core P6. Ce n'est évidemment pas portable mais cela concerne 95% du public. C'est un compteur 64 bits qui s'incrémente à chaque impulsion d'horloge et qui démarre au lancement de la machine. Il doit faire un tour complet en un millénaire.
inline unsigned long long int GetTick()
unsigned long long int x;
/* __asm__ volatile (".byte 0x0f, 0x31" : "=A" (x));*/
__asm__ volatile ("RDTSC" : "=A" (x));
return x;
On utilise ici une autre propriété non portable offerte par GCC, les long long. Ce sont en fait des nombres entiers sur 64 bits. On utilise ici un code assembleur à l'intérieur du code C. RDTSC (Read Time-Stamp Counter) est le vrai mnémonique assembleur pour lire ce compteur. Si l'assembleur utilisé est trop vieux, on peut utiliser la directive .byte qui écrit directement la valeur binaire dans le code. Donc, 0f31 est le code de l'instruction RDTSC.
L'instruction __asm__ de gcc est bien plus puissante que l'habituelle inclusion de code dans le fichier assembleur. Le "=A" (x) renseigne le compilateur sur les variables à manipuler. Ici, cela signifie que l'on utilise les registres du processeur de type A, le registre 64 bits (EDX:EAX). Le = signifie que l'on écrit dedans. (x) indique que l'on utilise la variable x. (Plus d'informations peuvent être trouvées dans la documentation de gcc dans C Extensions >Extended Asm).
Cette fonction renvoie un nombre qui correspond au nombre d'impulsions d'horloge. Bien sûr, le fait de récupérer cette information coûte un peu de temps puisqu'il faut bien sauvegarder EDX:EAX quelque part. Cela correspond à 32 cycles, ce qui est souvent complètement négligeable. Donc, avec une machine à 300 Mhz, en divisant par 300 000 000, on obtient le temps en secondes (n'oubliez pas de passer en float ;p).
Mais on peut faire plus précis en allant chercher les infos là où elles sont :
$cat /proc/cpuinfo | grep cpu MHz
cpu MHz : 300.691
J'ai donc ici une fréquence de 300 691 000 Hz.
En fait, on peut aller encore plus loin. Intel a disposé dans ces processeurs des compteurs d'événements, 2 dans le P6 et 18 dans le Pentium 4. Il existe des dizaines d'événements que l'on peut surveiller. Par exemple, on peut regarder le nombre de cycles où le processeur ne fait rien (à cause de dépendances ou d'attente de données), ou encore le taux de cache miss (voir un peu plus bas).
Malheureusement, si on peut lire les compteurs en mode utilisateur, on doit être en mode kernel pour dire quel événement doit être surveillé. Les instructions en question sont : RDMSR (Read Model Specific Counter), WRMSR (write Model Specific Counter) et RDPMC (read performance-monitoring counter).
En fait, avec le kernel 2.4, on peut accéder aux msr au travers de /dev/cpu/%d/msr. A vous de voir comment ! Il n'existe pas encore de documentation. Il n'y a pas de 'man msr', il faut lire les sources du noyau : 'locate msr' pour les trouver.
Les arcanes des processeurs
Le core P6 dispose d'un pipeline de longueur 10, l'Athlon, 14, le Pentium 4, 20. On peut parler aussi des 14 étages de l'UltraSparc III. Cela fait donc au moins 7 instructions en cours d'exécution à chaque instant pour le P6 (c'est le principe du pipeline : découper chaque tâche en morceau pour pouvoir augmenter la fréquence d'horloge, car il y a moins de choses à faire entre 2 coups d'horloge ; c'est une sorte de parallélisme car 7 (un peu moins que 10) instructions pour le P6 sont en cours d'exécution par pipeline, mais il faut 10 coups d'horloge pour finir une opération).
Dans ces nouveaux processeurs, plusieurs instructions peuvent être exécutées en même temps (5 instructions pour le P6 et jusqu'à 20 pour l'Itanium).
Le P6 dispose de 5 pipelines : 2 pour les unités entières dont un utilisé aussi pour les opérations flottantes, 2 pour faire des stores (store address calculation unit avec un buffer de 12 entrées et store data unit avec aussi un buffer de 12 entrées ), un pour le load avec 16 load buffer.
Donc, 7*5=35 instructions peuvent tourner au même instant. A l'extrême, l'alpha 21264 peut exécuter jusqu'à 80 instructions simultanément ! Le problème consiste à remplir complètement et en permanence ce pipeline et éviter les insertions de "bulles" (pseudo-instructions dans le pipeline qui ne servent à rien).
Dans les calculs, il faut considérer les dépendances entre les instructions. Cela signifie qu'une instruction ne peut pas s'exécuter tant que l'instruction précédente ne soit terminée. Par exemple,
1) x=a+b
2) y=x*a
Ici, les deux instructions ne peuvent être exécutées en parallèle, donc l'addition doit se finir avant de passer à la multiplication (car x doit être calculé avant). On parle de dépendance vraie ou dépendance de flot. Les autres dépendances sont les antidépendances (write after read) ou les dépendances de sortie (write after write) qui peuvent être évitées. Les vraies dépendances sont les dépendances dites RAW (read after write). Les dépendances RAR (read after read) ne gênent pas vraiment.
1) x = a + b
2) y = x * a
3) x = c + d
Ici, il y a toujours une vraie dépendance entre 1) et 2) sur x. Mais il y a une fausse dépendance entre 2) et 3) sur x. En effet, dans 3), il y a l'utilisation de x qui sert de variable temporaire. Or, on pourrait utiliser une autre variable et casser cette dépendance.
1) x = a + b
2) y = x * a
3) x' = c + d
Ainsi, 2) et 3) peuvent s'effectuer en parallèle. Heureusement, la plupart des processeurs sont là pour accélérer les vieux binaires. Ils sont donc capables de renommer x en x' et réorganiser les calculs s'il y a une fausse dépendance.
1) x = a + b
2) y = x * a
3) x' = c + d
4) w = e + f
Ici, le processeur peut bouger la 3ème et la 4ème instruction entre 1) et 2). Alors, pourquoi s'ennuyer à y faire attention ? Simplement parce que la fenêtre de travail est forcément limitée. De plus, on peut souvent ré-écrire le programme pour éviter les vrais dépendances.
Au niveau du C, le fait de prendre un maximum de variables temporaires plutôt que de les réutiliser aide le compilateur à mieux optimiser. Car il est pour l'instant incapable de comprendre la sémantique du programme. Par exemple :
s = a + b + c + d
correspond en fait à s = ((a + b) + c) + d , or s = (a + b) + (c + d) et est en réalité plus efficace car deux additions peuvent être effectuées en parallèle.
Et cela devient :
tmp1 = a+b
tmp2 = tmp1+c
s = tmp2+d
versus
tmp1 = a+b || tmp2 = c + d
s = tmp1 + tmp2
|| signifie "s'exécute en parallèle avec".
Les progrès des compilateurs font que ces optimisations pourraient être appliquées automatiquement. Pour faire cela, les compilateurs devraient savoir que l'addition est commutative (a+b=b+a), transitive ( (a+b)+c = a+(b+c) ) pour mettre en oeuvre la réécriture. Mais, selon la norme IEEE, le compilateur n'a pas le droit de toucher à l'arithmétique flottante à moins de le préciser par une option (-ffast-math pour gcc) et... qu'il puisse faire quelque chose. Donc, tout ce qui touche aux flottants ne sera, a priori, jamais optimisé à cause des problèmes d'arrondi.
Voici un programme exemple, compilé avec les options précédentes, mais avec un '-S' au lieu du '-o t'. Ainsi, la compilation s'arrête au fichier assembleur.
gcc -mcpu=i686 -O3 -g -S test.c -fstrength-reduce -frerun-loop-opt -fexpensive-optimizations -fschedule-insns2
-funroll-loops -fomit-frame-pointer -malign-double -fno-strict-aliasing -pipe -malign-loops=2 -jumps=2
-malign-functions=2 -DCPU=686 -ffast-math
#include<stdio.h>
#include<sched.h>
main()
int s;
int a,b,c,d;
long long t1,t2;
const struct sched_param *schedParam;
sched_setscheduler(0,SCHED_FIFO,schedParam);
scanf("%i %i %i %i",&a,&b,&c,&d);
__asm__ volatile ("RDTSC" : "=A" (t1));
s = ((a+b)+(c+d));
__asm__ volatile ("RDTSC" : "=A" (t2));
printf("opt : %Ldn(%i)n",t2-t1,s);
scanf("%i %i %i %i",&a,&b,&c,&d);
__asm__ volatile ("RDTSC" : "=A" (t1));
s = a+b+c+d;
__asm__ volatile ("RDTSC" : "=A" (t2));
printf("Usuel : %Ldn(%i)n",t2-t1,s);
Les printf et les scanf assurent que l'optimiseur n'optimisera "pas trop", en enlevant les calculs qu'il juge inutiles car reliés à aucune sortie. On remarque l'utilisation de %Ld dans le printf pour imprimer un long long.
movl %eax,120(%esp)
movl %edx,124(%esp)
movl 140(%esp),%eax
movl 132(%esp),%edx
addl 128(%esp),%edx
addl 136(%esp),%eax
addl %eax,%edx
movl %edx,100(%esp)
versus
movl 132(%esp),%esi
movl %eax,120(%esp)
addl 128(%esp),%esi
movl %edx,124(%esp)
addl 136(%esp),%esi
addl 140(%esp),%esi
On voit donc plus d'instructions dans le premier cas, mais l'exécution est en fait plus rapide ! 140(%esp) est une indirection vers la pile, donc le temps passé dans les movl est passé également dans les addl correspondants. On remarque aussi qu'il y a beaucoup moins de dépendances vraies. Dans le deuxième cas, il y a une dépendance entre presque chaque instruction (%esi), donc même une exécution dans le désordre ne pourra pas accélérer le code contrairement au premier cas.
Pour exécuter le programme, vous devez entrer des valeurs pour les variables. Pour 1 2 3 4 écrit deux fois, l'exécution du programme rend :
Opt : 38
(10)
Usuel : 59
(10)
Cela donne un avantage pour la version avec parenthèses mais vous remarquerez que les chiffres changent à chaque lancement. De plus, si l'on augmente le nombre de variables, le nombre d'accès à la mémoire augmente (pour stocker les variables temporaires), car le x86 a un tout petit nombre de registres. On peut, par exemple, pour 8 variables, faire quelque chose comme ceci : s=(a+b)+(c+d)+(e+f)+(g+h).
Toujours pas convaincu ? On va voir ce que l'on appelle de la réécriture de boucle. Par exemple, il s'agit simplement d'additionner un long vecteur.
int scatter1(int n, int A[n])
int i,ret=0;
for(i=0;i<n;i++)
ret +=A[i];
return ret;
Ce n'est pas compliqué, hein ? Comme vous vous en doutez, on peut faire mieux :
int scatter2(int n, int A[n])
int i,ret=0;
int tmp1=0,tmp2=0,tmp3=0,tmp4=0;
for(i=0;i<n;i+=4)
tmp1 +=A[i];
tmp2 +=A[i+1];
tmp3 +=A[i+2];
tmp4 +=A[i+3];
ret = (tmp1+tmp2)+(tmp3+tmp4);
return ret;
Donc, ma toute petite boucle vient d'être éclatée en quatre. Cela ressemble à comment faire simple quand on peut faire compliqué. Et pourquoi quatre ? Bah, parce que ! En fait, par essais/erreurs, 4 a été élu le meilleur chiffre. C'est un rapport entre le nombre d'instructions exécutables en parallèle, la profondeur du pipeline et les dépendances associées, le nombre d'accès à la mémoire en attente consécutive, la latence des instructions et le rapport entre la vitesse du core (le c
ur du processeur) et celle de la mémoire.
Si les vecteurs sont trop petits, on peut perdre des performances. Mais on va voir plus loin comment faire mieux.
Les sauts
Un autre problème pour garder le pipeline rempli provient des sauts dans le programme.
Lorsqu'il y a une rupture du flot d'instructions, le contenu du pipeline doit être vidé avant de commencer l'exécution des instructions en provenance d'une autre adresse. Certaines vieilles études assuraient qu'il y avait un saut toutes les trois instructions en moyenne. Avec une vingtaine d'instructions minimum possibles en cours d'exécution, le processeur passerait son temps à ne rien faire.
Heureusement, les processeurs essayent de préparer un éventuel saut vers des adresses qu'ils croisent. Ils font de la prédiction de branchement, ce qui assurent un cycle de pénalité, au pire, si la prédiction est juste. Sinon, plus d'une dizaine de cycles peuvent être perdus.
Il faut donc minimiser l'utilisation de saut (if, la plupart du temps) car la taille des tables de prédictions est forcement limitée. Des performances catastrophiques peuvent être obtenues si la probabilité de faire un saut ou non (on dit prendre un saut, on parle de taken or not taken en anglais) est de 50%. Dans ce cas, le système va se tromper plus de 50% du temps.
Il existe plusieurs méthodes pour tenter de prédire un branchement. Une étude rapide considère que dans le cas de boucle, 90 % du temps, le branchement est pris. Dans le cas d'une clause if, on ne peut pas savoir. De fait, une prédiction qui donne toujours taken a raison 60-70% du temps.
Mais si la pénalité est de 10 cycles, en cas d'erreur dans 30% des cas, 30% du code prend, grossièrement, 10 fois plus de temps à s'exécuter que le reste, qui devient du même coût presque négligeable en temps d'exécution par comparaison.
La méthode la plus courante utilise une fsm (finite state machine ou machine à état fini) à quatre états. Ceux-ci s'appellent : strong taken, taken, not taken, strong not taken. Imaginons que l'état initial soit taken ; selon la première étude, la fsm prédit que le branchement est pris. S'il l'est effectivement, elle passe en strong taken, sinon en not taken, et ainsi de suite. Chaque "bonne" prédiction (deux de suite en fait) renforce la croyance du système dans ses convictions.
On comprend bien ce qui se passe avec une boucle for ou while. La plupart du temps, la boucle est prise un certain nombre de fois avant de continuer. La prédiction reste ainsi à taken. Par contre, si une clause 'if' est "rarement" effectuée, elle est sautée par défaut.
Dans les processeurs les plus avancés, il existe une hiérarchie dans les tables qui conservent ces valeurs. La question est de savoir, lorsque l'on n'a plus de place, quelle ancienne valeur liée à une adresse de saut doit être remplacée.
Le P6 utilise encore un système plus complexe avec quatre niveaux de profondeur au lieu de deux. La table (dite BTB) stocke 512 adresses. Une erreur coûte en moyenne entre 10 et 15 cycles et au pire 26 cycles ! Intel utilise aussi un prédicat statique qui est utilisé lorsque le processeur voit le saut pour la première fois. Le processeur prévoit de prendre le saut s'il repart en arrière car il s'agit certainement d'une boucle. Si le saut va vers l'avant, il prévoit de ne pas le prendre considérant qu'il s'agit d'une clause if.
Il existe une nouvelle instruction dans les P6 qui permet d'éviter les vrais sauts. Il s'agit de l'instruction CMOV. Il s'agit de copier le contenu d'un registre dans un autre sous la condition de la valeur d'un 3ième. L'avantage est de ne pas toucher au flux des instructions. Cette instruction peut être utilisée dans tous les cas où des résultats de calcul sont nécessaires pour faire un choix. Souvent, la prédiction de branchement est très mauvaise sur ce genre de cas. Connaissant la pénalité due aux mauvaises prédictions, il vaut mieux faire quelques calculs en plus, puis faire un choix avec CMOV.
L'influence des instructions du processeur
Par ailleurs, les instructions disponibles influent sur les performances que l'on peut obtenir.
Par exemple, il est conseillé d'utiliser des entiers non signés pour les compteurs de boucles et les index des tableaux. Notamment, les divisions et le modulo non signé sont plus rapides. La fameuse multiplication, ou division par deux, peut être remplacée par un shift uniquement en non signé.
A l'inverse, si vous voulez convertir un entier en flottant, il vaut que cela soit un entier signé car il existe une instruction dédiée, contrairement aux entiers non signés pour les codes x86.
Il vaut mieux parfois faire certaines optimisations à la main plutôt que de faire confiance au compilateur et être sûr des performances. Il s'agit entre autres des déroulements des boucles. Toutes les boucles de moins d'une vingtaine d'itérations méritent d'être déroulées (cela correspond aux différentes latences et pénalités en cas de mauvaises prédictions).
Il faut aussi réduire au maximum le contenu interne des boucles. Cela paraît évident puisque le contenu d'une boucle est fait pour être répété. Mais on ne pense pas toujours qu'il vaut mieux sortir un if et dupliquer la boucle dans le code que de laisser le if à l'intérieur. for() if()... else... doit être remplacé par : if()for()... elsefor()...
Il vaut mieux éviter les divisions qui sont les instructions les plus lentes (18 cycles horloges en flottant 32 bits contre un pour une multiplication). On peut parfois les remplacer par des multiplications.
m =i /j /k;
m =i /(j * k); est plus rapide.
Il vaut mieux extraire à la main les sous-expressions des calculs et utiliser des variables temporaires autant que possible. Cela évite de refaire des calculs inutiles, surtout avec les variables flottantes où il est interdit de changer l'ordre des calculs à cause des problèmes d'arrondis. Mais il ne faut pas non plus introduire plus de "choses à faire" qu'auparavant. En effet, l'accès à la mémoire peut être plus coûteux que de refaire une addition.
Avec des nombres flottants :
a1=b/c
a2=d/c
En mieux :
t=1/c
a1=b*t
a2=d*t
il vaut mieux déclarer les variables dans l'ordre croissant de taille pour favoriser l'alignement 64 bits et éviter de futurs mauvais alignements et les pénalité associées (comme gcc aligne les données, cela fera gagner de la place mémoire). Il faut mettre les doubles avant les int, qui doivent être avant les short, etc.
Lorsque l'on écrit un switch, les conditions sont testées les unes après les autres. Donc, il vaut mieux mettre la clause la plus probable en premier et ainsi de suite.
On peut faire aussi ce genre de choses avec les clauses if. Cela se complique avec les optimisations du C. En effet, lors d'une opération booléenne avec un || (OR), si la première clause est vraie, quel que soit le résultat de la seconde, le résultat sera vrai. Il est donc inutile d'exécuter la deuxième clause. Il se passe la même chose avec un && (ET) si la première expression est fausse. La première clause doit donc être la clause qui rend le résultat qui a le plus d'importance dans la décision (vrai avec un OU, faux avec un ET), ou qui prendra le moins de temps de calcul.
Mais si on prend en compte la prédiction de branchement, cela se complique encore un peu plus. A chaque exécution de clause correspond un saut et la prédiction associée. Il faut donc utiliser au maximum des clauses faciles à prédire (revoir plus haut). S'il reste encore des clauses impossibles à prédire, il vaut mieux les mettre le plus loin possible dans le test du if.
Pour l'instant, on a essayé de se passer des optimisations du compilateur. Mais il y a deux moyens de l'aider. D'abord, l'utilisation de 'const' pour qualifier les paramètres d'entrée d'une fonction. Cela autorise plus d'optimisation. On peut aussi utiliser le préfixe 'static' devant les fonctions d'un même fichier. Cela signifie que la fonction n'est visible que dans le fichier correspondant et cela permet de faire de "l'inlining agressif" (le compilateur copie/colle le code de la fonction là où elle est appelée). L'intérêt est de gagner le code d'appel des fonctions et de supprimer les pénalités associées. Cela concerne toutes les petites fonctions. Ainsi, elles n'existent plus dans le binaire et ne peuvent pas être appelées par un autre objet (provenant d'un autre fichier, par exemple).
Avant de continuer la description des processeurs, je vais terminer par un dernier conseil, qui me fait particulièrement plaisir, car les conseils m'ont traumatisé quand j'étais petit. Il concerne les pointeurs. Pour faire bref, supprimez-les ! En fait, les compilateurs ne savent pas trop comment optimiser les accès aux données les utilisant.
Par exemple,
int a ;
int * p ;
....
p = & a ;
...
f(a)
g(*p)
Cela s'appelle le problème d'aliasing. La valeur contenue dans a est dotée de deux noms. A moins de faire une étude complète sur la vie du pointeur, le compilateur n'a pas d'autre choix que d'y accéder en passant systématiquement par la mémoire, alors qu'il pourrait mettre la valeur dans un registre.
D'une manière générale, un compilateur ne saura pas quoi faire d'un pointeur sans de gros calculs. Il vaut mieux essayer d'utiliser un tableau, il est beaucoup plus facile d'y extraire du parallélisme. Il semble que gcc 3.0 gère mieux le problème.
Les accès mémoires
Le principal goulot d'étranglement des processeurs actuels est la bande passante avec la mémoire.
Un core P6 dispose de 3 décodeurs, c'est-à-dire qu'il peut lire trois instructions en même temps. Or, une instruction x86 fait entre 2*8 et 15*8 bits. Si l'on considère une valeur moyenne, 16 bits, un processeur à 1 Ghz a besoin de 3*2*1=6 Go/s de bande passante... Or, le bus système à 133 Mhz propose une bande passante de 133*64/8=1,064 Go/s en crête (cela n'est vrai que lors de lecture d'adresse consécutive (mode burst), un cycle réellement aléatoire (qui correspond à la latence plus le cycle d'un accès lui-même) prend aujourd'hui autour de 50 ns, soit une bande passante de 20*64/8=160 Mo/s. D'ailleurs, si les RDRAM ont si mauvaise presse, c'est parce que leur latence tourne plutôt autour de 60 ns.).
En fait, je n'ai compté que les accès pour le fetch ("pour aller chercher") des instructions. Or, si on regarde le code généré plus haut, il y a quasiment un accès mémoire (sur 32 bits) par instruction. Il faut donc 1*32/8= 4 Go/s en plus. Si on utilise des instructions MMX ou SSE, les données qui transitent sont respectivement sur 64 et 128 bits. Il faut donc encore doubler, voire quadrupler le chiffre ! On arrive ainsi à 4*4+6 = 22Go/s dans le cas extrême...
Le pauvre sous-système mémoire ne peut fournir en burst que 1Go/s, voire maintenant 2Go/s avec les mémoires DDR-SDRAM. Vous comprenez mieux ainsi l'utilité de la mémoire cache.
Le premier niveau de mémoire cache doit donc absolument tourner à la même vitesse que le processeur, sinon celui-ci ne pourra jamais aller à pleine vitesse. La quantité de mémoire de ce cache L1 tourne autour de 64 Ko-128 Ko. Le cache L2, 256 Ko pour l'Athlon, est un peu plus lent. Il a besoin de quelques cycles pour être accédé. De plus, il est "direct mapped". En gros, cela signifie que les données ont un peu moins de chances de s'y trouver.
Comment fonctionne une mémoire cache ? Un cache de 256 Ko dispose de 18 bits d'adresse (2^18/1024 = 256). On sait que les lignes de cache font 32 octets sur les P6. Les 32 octets s'adressent avec 5 bits (l'adresse est toujours relative à l'octet, il n'y a guère que les DSP qui font autrement). On peut voir les 32 octets comme le paquet de base. L'adresse des processeurs est en gros sur 32 bits (les X86 seraient prévus jusqu'à plus de 40 bits de mémoire physique, mais avec toujours uniquement 32 bits pour chaque processus).
Si, dans le cache, on veut accéder aux données, on accède à la ligne ; il ne reste donc que 13 bits d'adressage. Lors d'un accès mémoire, l'adresse est découpée en trois morceaux. Les lsb (least significant bit ou bit de poids faible; 7 bits) sont mis de côté pour pouvoir choisir le bon mot dans une ligne. Le deuxième paquet (13 bits) accède directement au cache et à une mémoire de tag (de drapeaux ou d'informations). Les msb (most significant bit ou bit de poids fort), les 14 bits restants, sont comparés avec la sortie de la mémoire tag. Si les valeurs sont identiques, cela signifie que le cache contient bien la donnée recherchée. Sinon, il faut accéder à la mémoire principale.
Vous n'avez rien compris ? Essayez de comprendre le schéma avant de relire le paragraphe précédent, vous verrez, ce n'est pas si complexe.
On se rend bien compte de la limite de la méthode (si, si), puisque chaque adresse de la mémoire principale multiple de la taille du cache se reflète dans la même adresse de celui-ci. Dit autrement, une case mémoire ne se reflète que dans une seule et unique case de la cache. On peut contourner le problème en mettant l'équivalent de plusieurs caches en parallèle. S'il y en a quatre, on appelle cela des 4-ways set associative cache. Donc, il a quatre emplacements possibles pour chaque adresse mémoire. S'il y a un mécanisme pour chaque ligne, on dit que le cache est full associative. C'est rare d'en trouver ; plus le nombre de voix augmente et moins l'efficacité augmente. 16-way est peut-être ce que l'on trouve de mieux.
Malheureusement, les caches multivoies ont beaucoup plus de contrôle, notamment la mémoire de tag qui devient assez grosse. Donc, le cache L1 est souvent multivoie pour augmenter les chances de hits ; ce sont les mémoires les plus rapides. Pour le cache L2, on utilise la mémoire la plus dense possible, qui est 1-way pour limiter la taille de la mémoire de tag (tout en gardant un temps d'accès acceptable). Ici, la mémoire tag est de 14 Ko pour 256 Ko de mémoire.
Un cache miss signifie que la donnée (ou l'instruction) ne se situe pas dans le cache et qu'il faut aller la chercher dans le niveau de cache supérieur ou dans la mémoire principale. La pénalité tourne alors autour de 150 cycles, où le processeur ne fait absolument rien s'il y a une dépendance de flot.
Il faut donc absolument rester en mémoire cache. Ces petites mémoires ultra-rapides n'ont de sens que si les programmes exécutés respectent le concept de localité. Cela signifie que la prochaine information que l'on va chercher est proche de la position de l'information en cours de traitement. Ainsi, lors d'un accès à la mémoire centrale, une lecture burst est effectuée et au moins 32 octets (pour le P6, 64 pour l'Athlon) sont lus consécutivement pour remplir "une ligne de cache". Ainsi, à la prochaine lecture, on pense avoir plus de chance de trouver la donnée suivante. L'utilisation de burst permet aussi de mieux utiliser les accès à la mémoire.
Il existe aussi un autre petit cache constitué de 12 "store buffer" (core P6) où les données sont stockées avant d'être écrites vers la mémoire pour rendre la main beaucoup plus rapidement. C'est une sorte d'écriture asynchrone mais où l'ordre est préservé. Ces buffers servent à masquer la latence de la mémoire. Pour gagner du temps, il faut toujours éviter d'avoir plus de 12 writes dans le même bloc d'instructions. Dans le cas contraire, le cpu doit attendre l'écriture effective dans la mémoire.
Mais attention, ce temps peut être long. Cela dépend du mode dans lequel fonctionne le cache : Write-thought and write-back. L'un écrit d'abord dans le cache qui retransmet effectivement dans la mémoire quand il doit évincer une ligne de cache. Donc, le contenu de la SDRAM n'est pas cohérent avec le contenu du cache. L'autre écrit directement dans la mémoire centrale. C'est lent mais la mémoire est cohérente.
Le strip-mining est une technique qui consiste à casser des boucles pour augmenter la localité des données. Imaginons que l'on manipule un grand tableau de N*M valeurs, int tab[n][m]. On sait que le déplacement d'une adresse se fait sur le dernier indice pour gcc (*(tab+i*n+j)=tab[i][j]). Mais il semble que cela ne soit pas le cas de tous les compilateurs.
for(k=0;k<K;k++)
for(j=0;j<m;j++)
for(i=0;i<n;i++)
..
tab[i][j]
..
En imaginant que n et m soient assez grands, tab ne peut pas rentrer dans le cache en entier. Sachant qu'à chaque cache miss, une ligne de 32 octets est chargée, il faut en profiter et travailler sur les données alignées !
Il faut déjà inverser les deux compteurs, pour que se soit J qui soit au plus profond de la boucle pour respecter au mieux l'alignement des données.
De plus, on parcourt K fois ce tableau. Or, à chaque fois, on pollue le cache (c'est-à-dire que l'on remplace des données par d'autres que l'on ne va finalement pas réutiliser ; bref, ici, le cache ne sert presque à rien). Il serait beaucoup plus intéressant de travailler K fois sur un bout du tableau et passer au suivant. Pour cela, on doit réécrire l'algorithme en rajoutant un niveau.
for(ii=0;ii<m;ii+=II)
for(k=0;k<K;k++)
for(i=ii;i<ii+II;i++)
for(j=0;j<n;j++)
..
tab[i][j]
..
On considère ici que II fait en sorte que les II*n éléments rentrent dans le cache, ainsi les trois boucles les plus profondes tourneront 90% du temps à l'intérieur de celui-ci.
Un exemple pratique
Je vous ai beaucoup parlé. Maintenant, passons à un cas beaucoup moins théorique, où l'on voit tout l'intérêt de cet article.
Un problème ultra classique est la multiplication de matrice. Habituellement, cela s'écrit :
void mulMatrix1(int n,int A[n][n] ,int B[n][n],int C[n][n])
int i,j,k;
for(i=0; i<n ;i=i+1)
for(j=0; j<n ;j=j+1)
for(k=0; k<n ;k=k+1)
C[i][j] = C[i][j] + A[i][k]*B[k][j];
L'ordre des trois boucles n'a pas l'air d'avoir d'importance. Or, on voit que i et j interviennent sur C et alternativement sur A, puis B. Ici, k est à l'intérieur des boucles. On casse une dépendance entre deux itérations de boucles (puisque les C[i][j] sont des cases différentes). Par contre, il faut aller chercher le contenu de deux cases mémoires.
void mulMatrix2(const unsigned int n, const int A[n][n] , const int B[n][n],int C[n][n])
register unsigned int i,j,k;
for(i=0; i<n ;i=i+1)
for(k=0; k<n ;k=k+1)
for(j=0; j<n ;j=j+1)
C[i][j] = C[i][j] + A[i][k]*B[k][j];
Ici, j est à l'intérieur des boucles. Il y a une dépendance entre C[i][j], mais ici on accède au même A[i][k]. De plus, B[k][j] court ainsi sur des adresses alignées. On peut compter aussi sur les write buffers pour récupérer plus facilement les données.
void mulMatrix3(const unsigned int n, const int A[n][n], const int B[n][n],int C[n][n])
register unsigned int i,j,k;
register int tmp1,tmp6;
for(i=0; i<n ;i=i+1)
for(j=0; j<n ;j=j+1)
tmp1=0;
tmp6 = C[i][j];
for( k=0; k<n ;k=k+1 )
tmp1 += A[i][k] * B[k][j];
C[i][j] = tmp6+tmp1;
Ici, tout a été réécrit pour casser les dépendances le plus possible. On retarde au maximum tout accès à la mémoire. Mais, c'est surtout pour la forme : vous comprendrez en regardant les performances.
void mulMatrix4(const unsigned int n, const int A[n][n], const int B[n][n],int C[n][n])
register int tmp6;
register unsigned int i,j,k;
for(i=0;i<n;i=i+1)
for(k=0; k<n ;k=k+1)
tmp6 = A[i][k];
for(j=0;j<n;j=j+1)
C[i][j] = C[i][j] + tmp6*B[k][j];
La même chose, mais en utilisant l'autre optimisation d'ordre de boucle possible, qui est nettement plus intéressante. Ici, les trois accès sont optimisés : C[i][j] avec les write buffers, A[i][k] avec ce que l'on appelle un preload, et pour B[k][j], les données sont utilisées alignées.
Comment faire mieux ? Avec le strip-minning ! Let's go !
void mulMatrix5(const unsigned int n, const int A[n][n], const int B[n][n],int C[n][n])
register int tmp6;
register unsigned int i,j,k,ii,kk;
for(ii=0; ii<n ;ii=ii+K)
for(kk=0; kk<n ;kk=kk+K)
for(i=ii;i<ii+K;i=i+1)
for(k=kk; k<kk+K ;k=k+1)
tmp6 = A[i][k];
for(j=0;j<n;j=j+1)
C[i][j] = C[i][j] + tmp6*B[k][j];
Ici, on découpe i et k, on laisse j car on veut profiter au maximum de l'alignement des données. Cela fonctionne bien uniquement si une ligne complète de la matrice rentre dans le cache. Sinon, on perdrait l'intérêt de ce découpage.
Qu'est-ce qu'il reste à faire ? On a optimisé en cassant au maximum les dépendances et en optimisant la gestion du cache. On peut encore essayer de faire mieux en "prefetchant" les données et en déroulant la boucle interne à la main.
Le prefetch existe dans la plupart des processeurs récents. Tout se passe comme si une lecture de la mémoire était faite. Mais rien ne change dans le flot du programme. Alors, à quoi cela sert-il ? Comme d'habitude, à remplir les caches !
void mulMatrix11(const unsigned int n, const int A[n][n], const int B[n][n], int C[n][n])
register int tmp6,tmp7,tmp8,tmp9;
register int tmp10,tmp11,tmp12,tmp13;
register unsigned int i,j,k,ii,jj,kk,jjj,kkk;
for(ii=0; ii<n ;ii=ii+K)
for(kk=0; kk<n ;kk=kk+K)
for(i=ii;i<ii+K;i=i+1)
for(kkk=kk;kkk<kk+K;kkk+=4)
__asm__ __volatile__ ("prefetchnta 16(%0) " : : "S" (&(A[i][kkk+4])) );
tmp6 = A[i][kkk];
tmp7 = A[i][kkk+1];
tmp8 = A[i][kkk+2];
tmp9 = A[i][kkk+3];
for(j=0;j<n;j=j+128)
__asm__ __volatile__ ("prefetchnta 512(%0)" : : "D"(&(B[kkk][j+128])));
__asm__ __volatile__ ("prefetchnta 512(%0)" : : "S"(&(B[kkk+1][j+128])));
__asm__ __volatile__ ("prefetchnta 512(%0)" : : "D"(&(B[kkk+2][j+128])));
__asm__ __volatile__ ("prefetchnta 512(%0)" : : "S"(&(B[kkk+3][j+128])));
for(jjj=j;jjj<j+128;jjj+=2)
C[i][jjj] += tmp6*B[kkk][jjj]+tmp7*B[kkk+1][jjj]
+ tmp8*B[kkk+2][jjj]+tmp9*B[kkk+3][jjj];
C[i][jjj+1] += tmp6*B[kkk][jjj+1]+tmp7*B[kkk+1][jjj+1]
+ tmp8*B[kkk+2][jjj+1]+tmp9*B[kkk+3][jjj+1];
Sur une distribution pas trop récente, il est possible que as (le programme d'assemblage avec gcc) ne reconnaisse pas ce mnémonique (sur une Mandrake 8.0, cela marche : gcc version 2.96).
Le 128 qui traîne ne correspond à rien de précis. En fait, c'est un savant mélange entre le nombre d'octets nécessaires au calcul qui suit (le memory stride, ms), le nombre de cycles passés dans le c
ur de la boucle(nc). AMD donne la formule suivante pour calculer la distance à donner au prefetch : 200*ms/nc. Le 200 est dépendant de l'implémentation du chip (ici, l'Athlon). Eh oui, la vie est mal faite parfois. Cela correspond, entre autre, à la vitesse de la mémoire en rapport avec la vitesse du core.
Ici, pour faire des tests, on n'utilise que des matrices multiples de puissance de 2. C'est facile d'imaginer des petits calculs préliminaires pour pouvoir utiliser une taille de matrice quelconque. Ou encore, on peut couper le calcul en deux, avec une série de boucles simples pour finir le calcul.
Pour parler performances, sur ma machine d'une fréquence de 4.5*66=300, j'obtiens une augmentation de 45 % de performance sur des "petites" matrices en comparant le pire et le meilleur algorithme. Sur les matrices qui ne rentrent pas dans le cache, j'obtiens 400 % d'augmentation de performance (voire *10 dans un cas bien précis avec des matrices 512*512)... Qui dit mieux ? Je suppose qu'avec les nouveaux processeurs qui ont un rapport entre vitesse de core et vitesse de la mémoire encore plus grand, cela doit être encore plus important (j'ai vu un x25 avec un Athlon 1.2Ghz...). Convaincu ? NON ! Et bien compilez maintenant...
Sur notre exemple précédent, on peut écrire :
int scatter3(int n, int A[n])
int i,ii,ret=0;
int tmp1=0,tmp2=0,tmp3=0,tmp4=0;
for(ii=0;ii<n;ii=ii+32)
__asm__ __volatile__ ("prefetcht0 128(%0)" : : "S"(&(A[ii])));
for(i=ii;i<ii+32;i+=4)
tmp1 +=A[i];
tmp2 +=A[i+1];
tmp3 +=A[i+2];
tmp4 +=A[i+3];
ret = (tmp1+tmp2)+(tmp3+tmp4);
return ret;
inline int scatter5(int n, int A[n])
unsigned int i,j;int ret=0;
int tmp1=0,tmp2=0,tmp3=0,tmp4=0;
for(i=0;i<n;i+=4)
__asm__ __volatile__ ("prefetchnta 256(%0)" : : "r"(&(A[i])));
tmp1 +=A[i];
tmp2 +=A[i+1];
tmp3 +=A[i+2];
tmp4 +=A[i+3];
ret = (tmp1+tmp2)+(tmp3+tmp4);
return ret;
inline int scatter1b(int n, int A[n])
int i,ii,ret=0;
int tmp1=0,tmp2=0,tmp3=0,tmp4=0;
for(ii=0;ii<n;ii=ii+16)
__asm__ __volatile__ ("prefetchnta 64(%0)" : : "r"(&(A[ii+16])));
for(i=ii;i<ii+16;i+=1)
tmp1 +=A[i];
ret = (tmp1+tmp2)+(tmp3+tmp4);
return ret;
On peut toujours écrire des versions différentes du code. L'intérieur de la boucle étant très petit, les "prefetch" n'ont pas vraiment le temps de s'effectuer avant la demande effective des données. C'est pour cela que l'on rajoute un offset (128) pour aller "prefetcher" les 128/4=32 entiers plus loin. Le "prefetch" est asynchrone et on ne sait pas quand il sera effectivement terminé. Dans le cas où il ne serait pas fini, il occupe l'unité de load pour rien. Prévoir la bonne distance est très difficile. Intel conseille d'éviter cette instruction si la donnée rentre dans le cache.
Selon la taille des matrices et le processeur (et le rapport entre vitesse du core et celle de la mémoire), l'algorithme optimal change. Celui qui semble être assez constant est le 1b. Ici, on est en présence d'un traitement 'memory bound' (opposé à CPU bound). Cela signifie qu'il est complètement limité par la bande passante mémoire, donc le temps passé dans les calculs est négligeable devant l'attente des données venant de la mémoire.
Il existe 4 types de "prefetch" : prefetcht0, prefetcht1, prefetcht2 et prefetchnta. Les trois premiers dépendent des caches vers lesquels doivent être transférées les données (L1 et L2, uniquement L2). Avec le core P6, prefetcht1 et prefetcht2 sont identiques. Le dernier ne remplit que L1 et ne touche pas à L2, "pour minimiser la pollution du cache", sachant que le cache L2 contient aussi du code. Pour le P4, le prefetch ne concerne que le cache L2 (les instructions sont donc identiques).
Pour vous convaincre de bannir à jamais les pointeurs de votre code, testez :
void mulMatrix14(const unsigned int n, const int *A, const int *B, int *C)
register int tmp6,tmp7,tmp8,tmp9,*p,*p2,*p3;
register int tmp10,tmp11,tmp12,tmp13;
register unsigned int i,j,k,ii,jj,kk,jjj,kkk;
long long tickv,tickp;
for(ii=0; ii<n ;ii=ii+K)
for(kk=0; kk<n ;kk=kk+K)
for(i=ii;i<ii+K;i=i+1)
for(kkk=kk;kkk<kk+K;kkk+=4)
p = (A+i*n+kkk);
__asm__ __volatile__ ("prefetchnta 16(%0) " : : "S" ( (p+4) ));
tmp6 = *(p);
tmp7 = *(p+1);
tmp8 = *(p+2);
tmp9 = *(p+3);
for(j=0;j<n;j=j+128)
p2=B+kkk*n+j+128;
__asm__ __volatile__ ("prefetchnta 512(%0)" : : "D"(p2) );
__asm__ __volatile__ ("prefetchnta 512(%0)" : : "S"(p2+n) );
__asm__ __volatile__ ("prefetchnta 512(%0)" : : "D"(p2+2*n) );
__asm__ __volatile__ ("prefetchnta 512(%0)" : : "S"(p2+3*n) );
for(jjj=j;jjj<j+128;jjj+=2)
p3=B+kkk*n+jjj;
*(C+i*n+jjj) += tmp6**(p3)+tmp7**(p3+n)
+ tmp8**(p3+2*n)+tmp9**(p3+3*n);
*(C+i*n+jjj+1) += tmp6**(p3+1)+tmp7**(p3+n+1)
+ tmp8**(p3+2*n+1)+tmp9**(p3+3*n+1);
Malgré les apparentes simplifications des calculs nécessaires à l'accès en mémoire, j'obtiens une perte de 25 % des performances par rapport à la version optimale.
Nous venons de survoler la mécanique interne des processeurs en vue de mieux les programmer en C. J'espère ainsi que certaines applications pourront gagner en performance.
Nicolas Boulay
Ingénieur en Microélectronique
Membre du projet f-cpu
URL
AMD :
http://www.amd.com/products/cpg/athlon/techdocs/index.html (le document sur l'optimisation est à lire absolument ; par contre, impossible d'avoir le document en entier avec les lecteurs de pdf disponible pour Linux (Ghostview, et même Acroread !) à l'exception de xpdf)
Intel :
http://developer.intel.com/design/PentiumIII/manuals/
Sun :
http://www.sun.com/microelectronics/manuals/index.html
http://www.sun.com/sparc/UltraSPARC-III/index.html
IBM :
http://www-3.ibm.com/chips/techlib/techlib.nsf/productfamilies/PowerPC
Digital :
www.alphapowered.com
fcpu :
www.f-cpu.org
Code du testeur
Tout le code est évidemment sous GPL.
#include<stdio.h≶
#include<sched.h≶
#include<stdlib.h≶
#include <sys/mman.h≶
/*
CPUHZ doit être récupéré par :
cat /proc/cpuinfo | grep cpu MHz
*/
#define CPUHZ 300689000
#define K 32
/*Pour vérifier que l'on n'a pas fait de bêtises*/
void affiche(int n, int S[n][n])
int tmp=0,i,j;
for(i=0;i<n;i=i+1)
for(j=0;j<n;j=j+1)
tmp += S[i][j];
printf("(%d)",tmp);
/*Pour être à peu près sûr de vider les caches si size > 1 Mo*/
void TrashCache(int * trash, int size)
int i;
for(i=0;i<size;i++)
(*trash)++;
void result(int size , long long tick,int *C)
/*
Je compte les accès mémoire 'théoriques' 2*n read + un write par case,
une version +réelle avec 2*n read et n write, le nombre de mops correspondant (n multiplication et n-1 addition par case ),
et le rendement face au maximum théorique du processeur.
Si cela ne vous plaît pas, faites votre propre tambouille.
*/
float flot = (size*size*(size*2+1)*4/(((float)(tick))/CPUHZ))/1000000;
float flot2 = (size*size*size*3*4/(((float)(tick))/CPUHZ))/1000000;
float mops = (size*size*(size+(size-1))/(((float)(tick))/CPUHZ))/1000000;
float fracThMops=100*1000000*(mops/(CPUHZ*1.5));/*rendement moyen 2 add/clk, 1 mul/clk mais il peut faire
UNE addition en même temps que 1 mul*/
affiche(size,C);
printf(" : %.2f Mo/s (%.2f Mo/s) - %.2f mops(%.0f%)n",flot,flot2,mops,fracThMops);
/*pour les scatter*/
void resultAdd(int size , long long tick )
float flot = (size*4/(((float)(tick))/CPUHZ))/1000000;
printf(" : %.2f Mo/sn",flot);
int main()
unsigned int tmp1,vect;
unsigned int i, j,ret;
int *A,*B,*C,*trash;
long long int t1,t2;
float ref;
const struct sched_param *schedParam;
sched_setscheduler(0,SCHED_FIFO,schedParam);
mlockall(MCL_CURRENT); /* to avoid swaping page to disk */
printf("Fréquence de base utilisée pour les calculs : %.1fMhzn",(float)CPUHZ/1000000);
printf("Fréquence définie par #define CPUHZn");
printf("Pour l'avoir : more /proc/cpuinfo | grep cpu MHzn");
printf("[fonction][%% en comparaison du 1er][vérification]:nt[débit théorique](débit 'réelle') - [nb d'ops](rendement P3)n");
for (tmp1=128;tmp1 < 1025;tmp1+=128 )
printf("matrice[%d][%d] (%.2f Ko)n",tmp1,tmp1,((float)(tmp1*tmp1*sizeof(int)))/1024);
A=(int*)malloc(sizeof(int[tmp1][tmp1]));
B=(int*)malloc(sizeof(int[tmp1][tmp1]));
C=(int*)malloc(sizeof(int[tmp1][tmp1]));
if( A == NULL || B == NULL || C == NULL )
printf("ERR: Plus de mémoiren");exit(0);
memset((A),0,sizeof(int[tmp1][tmp1]));
memset((B),0,sizeof(int[tmp1][tmp1]));
memset((C),0,sizeof(int[tmp1][tmp1]));
for(i=0;i<tmp1;i=i+1)
for(j=0;j<tmp1;j=j+1)
*(A+i*tmp1+j)=i+2*j;
for(i=0;i<tmp1;i=i+1)
for(j=0;j<tmp1;j=j+1)
*(B+i*tmp1+j)=i+3*j;
trash=(int*)malloc(sizeof(int)*1000000);
TrashCache(trash,1000000);
t1 = GetTick();
mulMatrix1(tmp1,A,B,C);
t2 = GetTick();
printf("mulMatrix1 (+%.1f%)t",((ref-(t2-t1))*100)/(t2-t1));
result(tmp1,t2-t1,C);
ref = t2-t1;
TrashCache(trash,1000000);
memset((C),0,sizeof(int[tmp1][tmp1]));
t1 = GetTick();
mulMatrix2(tmp1,A,B,C);
t2 = GetTick();
printf("mulMatrix2(futur ref) (+%.1f%)t",((ref-(t2-t1))*100)/(t2-t1));
result(tmp1,t2-t1,C); /*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
memset((C),0,sizeof(int[tmp1][tmp1])); /*reset de C*/
ref = t2-t1; /*on est gentil : algo simple le - con (déjà x2 en perf)*/
t1 = GetTick();
mulMatrix3(tmp1,A,B,C);
t2 = GetTick();
printf("mulMatrix3 (+%.1f%)t",((ref-(t2-t1))*100)/(t2-t1));
result(tmp1,t2-t1,C); /*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
memset((C),0,sizeof(int[tmp1][tmp1])); /*reset de C*/
t1 = GetTick();
mulMatrix4(tmp1,A,B,C);
t2 = GetTick();
printf("mulMatrix4 (+%.1f%)t",((ref-(t2-t1))*100)/(t2-t1));
result(tmp1,t2-t1,C); /*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
memset((C),0,sizeof(int[tmp1][tmp1])); /*reset de C*/
t1 = GetTick();
mulMatrix5(tmp1,A,B,C);
t2 = GetTick();
printf("mulMatrix5 (+%.1f%)t",((ref-(t2-t1))*100)/(t2-t1));
result(tmp1,t2-t1,C); /*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
memset((C),0,sizeof(int[tmp1][tmp1])); /*reset de C*/
t1 = GetTick();
mulMatrix11(tmp1,A,B,C);
t2 = GetTick();
printf("mulMatrix11 (+%.1f%)t",((ref-(t2-t1))*100)/(t2-t1));
result(tmp1,t2-t1,C); /*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
memset((C),0,sizeof(int[tmp1][tmp1])); /*reset de C*/
t1 = GetTick();
mulMatrix11p(tmp1,A,B,C);
t2 = GetTick();
printf("mulMatrix11p (+%.1f%)t",((ref-(t2-t1))*100)/(t2-t1));
result(tmp1,t2-t1,C); /*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
memset((C),0,sizeof(int[tmp1][tmp1])); /*reset de C*/
free(A);
vect =1024*tmp1;
A=(int*)malloc(sizeof(int[vect]));
printf("vecteur[%d] (%.2f ko)n",vect,(((float)vect*sizeof(int)))/1024);
if(A==0)
perror("Pas assez de mémoire !n");
exit(0);
for(i=0;i<vect;i++)
(*(A+i))=i;
t1 = GetTick();
ret=scatter1(vect,A);
t2 = GetTick();
printf("scatter1 (00%)t(%d)t",ret); /*De quoi on cause*/
resultAdd(vect,t2-t1);/*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
ref =t2-t1;
t1 = GetTick();
ret=scatter1b(vect,A);
t2 = GetTick();
printf("scatter1b(+%.1f%)t(%d)t",((ref-(t2-t1))*100)/(t2-t1),ret);
resultAdd(vect,t2-t1);/*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
t1 = GetTick();
ret=scatter2(vect,A);
t2 = GetTick();
printf("scatter2 (+%.1f%)t(%d)t",((ref-(t2-t1))*100)/(t2-t1),ret);
resultAdd(vect,t2-t1);/*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
t1 = GetTick();
ret=scatter3(vect,A);
t2 = GetTick();
printf("scatter3 (+%.1f%)t(%d)t",((ref-(t2-t1))*100)/(t2-t1),ret);
resultAdd(vect,t2-t1);/*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
t1 = GetTick();
ret=scatter4(vect,A);
t2 = GetTick();
printf("scatter4 (+%.1f%)t(%d)t",((ref-(t2-t1))*100)/(t2-t1),ret);
resultAdd(vect,t2-t1);/*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
t1 = GetTick();
ret=scatter5(vect,A);
t2 = GetTick();
printf("scatter5 (+%.1f%)t(%d)t",((ref-(t2-t1))*100)/(t2-t1),ret);
resultAdd(vect,t2-t1);/*l'affichage*/
TrashCache(trash,1000000); /*le zigouillage de cache*/
free(A);free(B);free(C);free(trash);
munlockall();
return 0;