TD 3-4. Introduction à CUDA
par Eric Goubault et Sylvie Putot
Utiliser CUDA dans les salles machines de l'X
Quelles cartes NVIDIA et ou sont-elles?
- En salles 30, 35 et 36 (liste des machines)
- Cartes GeForce GTX 260:
- Total amount of global memory: 938803200 bytes
- Number of multiprocessors: 24
- Number of cores: 192
- Total amount of constant memory: 65536 bytes
- Total amount of shared memory per block: 16384 bytes
- Total number of registers available per block: 16384
- Warp size: 32
- Maximum number of threads per block: 512
- Maximum sizes of each dimension of a block: 512 x 512 x 64
- Maximum sizes of each dimension of a grid: 65535 x 65535 x 1
- Maximum memory pitch: 262144 bytes
- Texture alignment: 256 bytes
- Clock rate: 1.24 GHz
Ecrire, compiler et exécuter un programme CUDA
- La première fois, il faut positionner les variables d'environnement. Créer ou ouvrir le fichier ~/.bashrc ou ~/.bash_profile et y ajouter:
PATH=$PATH:/usr/local/cuda/bin:/usr/local/cuda/open64/bin
LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/usr/local/cuda/lib
export PATH
export LD_LIBRARY_PATH
Puis recharger le fichier (source ~/.bashrc par exemple).
- NVIDIA propose un ensemble d'exemples, accompagné d'une librairie d'utilitaires (cutil), qui contient des timers et toutes les fonctions
en "CUT", dont en particulier des fonctions de récupération/affichage d'erreur.
- Pour compiler vos fichiers indépendamment (en dehors donc) de la SDK, vous pouvez, par exemple pour le répertoire template, partir du Makefile suivant
- Vous pouvez aussi tout mettre dans un seul fichier mon_essai.cu et compiler puis exécuter
par
nvcc mon_essai.cu -I/usr/local/cuda/include/ -L/usr/local/cuda/lib -o mon_essai
./mon_essai
Si vous voulez utiliser la librairie cutil de la sdk (divers utilitaires, dont les timers et tout ce qui est en cut), il faut ajouter des include et des librairies, par exemple:
nvcc mon_essai.cu -I/usr/local/cuda/include/ -I/users/profs/info/putot/NVIDIA_GPU_Computing_SDK/C/common/inc/ -L/usr/local/cuda/lib/
-L/users/profs/info/putot/NVIDIA_GPU_Computing_SDK/C/lib/ -lcutil_i386 -o mon_essai
La documentation
Une première mise en oeuvre : calcul de PI en parallèle
Implémenter le calcul de π en parallèle avec la formule suivante (cf aussi le TD1) :
\[ \pi = \int_0^1 \frac{4}{1+x^2} dx \approx \sum_{i=1}^n \frac{1}{n} \frac{4}{1+((i-\frac{1}{2})\frac{1}{n})^2}. \]
Principe général de calcul
Une solution est le paradigme Maître/Esclave.
Un maître va lancer N esclaves chargés de calculer les sommes
partielles
\[ P_k = \sum_{i=k*n/N+1}^{(k+1)*n/N} \frac{1}{n} \frac{4}{1+((i-\frac{1}{2})\frac{1}{n})^2} \, , \; k=0,\ldots,N-1, \]
puis faire la somme des résultats partiels.
Quelques pistes et suggestions pour l'implémentation CUDA
Vérification des résultats et affichage des performances
- Comparer le résultat obtenu sur GPU avec le calcul de Pi pour la même valeur de $n$
- Ajouter des mesures des temps de calcul et afficher l'accélération par rapport au programme séquentiel;
plusieurs solutions pour mesurer le temps de calcul:
- les fonctions habituelles C clock() ou time(),
- les timers CUDA (cudaEvent, cf le manuel Cuda C best practice),
- les timers de la SDK (voir template par exemple)
- Ne pas hésiter à récupérer les messages d'erreur éventuels en utilisant les fonctions CUT_SAFE_CALL/cutilSafeCall(attention sinon, l'exécutable peut très bien s'être planté sans rien dire...) ; pour les problèmes d'accès mémoire (voire pour avoir plus d'info sur ce qui est réellement calculé),
vous pouver utiliser cuda_memcheck (compiler avec -g -G)
Précision des calculs
- Si les calculs sont effectués en simple précision et si on prend n trop grand, un phénomène d'absorption fait que la précision des résultats devient mauvaise,
- Or avec les options de compilation par défaut, les double sont interprétés sur le GPU comme des float, et il y a
des problèmes de transmission entre CPU et GPU, donc appeler nvcc avec l'option -arch=sm_13
- Si vous utilisez des double, attention aux performances, un seul coeur par carte calcule en double, et les double sont
plus lents! Attention aussi, avec l'option -arch=sm_13, seuls les unités double seront utilisées, donc même si votre programme
utilise des float et non des double, les calculs seront lents...
Faire varier taille des grilles et des blocs
Faire varier n, le nombre de blocs, de threads par blocs et observer les performances relatives.
Quelques outils
Visual profiler
Occupancy calculator
Pour essayer de comprendre les performances: mode d'emploi inclus dans le fichier excel.
Le debugger
cuda_gdb (seulement en mode console), cuda_memcheck (pas de pb de mode console; compiler avec -g -G)
Simulation différences finies
Equation de Laplace
Nous voulons calculer la distribution de température dans une pièce carrée, sachant que la température sur les bords est fixée: la température sur les murs de la pièce est de 20 degrés, sauf sur l'un des murs (disons la cheminée...), ou elle est à 100 degrés.
En régime stationnaire, c'est-à-dire lorsque la température est stabilisée, la température vérifie l'équation de Laplace (bi-dimensionnelle ici) \[ \Delta T = \frac{\partial^2 T}{\partial x^2} + \frac{\partial^2 T}{\partial y^2}= 0.\]
Nous allons calculer la température sur une grille régulière par une méthode de différences finies classique,
en approchant les dérivées partielles au point (x,y) par
\[ \frac{\partial^2 T}{\partial x^2}(x,y) \approx \frac{1}{h^2} (f(x+h,y) - 2 f(x,y) + f(x-h,y) ), \]
\[ \frac{\partial^2 T}{\partial y^2}(x,y) \approx \frac{1}{h^2} (f(x,y+h) - 2 f(x,y) + f(x,y-h) ), \]
ou h est un pas de discrétisation suffisamment petit. En supposant qu'on a discrétisé le rectangle par une grille de N par N
points, alors la température en chaque point (i,j) de la grille peut
s'exprimer directement en fonction de ses voisins immédiats, par
\[ T_{i,j} = \frac{T_{i-1,j}+T_{i+1,j}+T_{i,j-1}+T_{i,j+1}}{4}, \; \forall i,j \in \{1,\ldots,N-2\}, \]
les températures pour i=0, i=N-1, j=0 et j=N-1 étant données par les conditions limites.
Calcul séquentiel de la température en chaque point
Pour calculer la température en chaque point, on va donc itérer ce calcul jusqu'à ce que la différence entre deux itérées successives devienne inférieure à un certain seuil. On va dans un premier temps implémenter la version de référence séquentielle et afficher les résultats.
Vous pouvez faire plusieurs versions:
- une première, la plus simple, qui pour un nombre d'itérations n fixé, calcule par l'itération de Jacobi:
\[ T^{k+1}_{i,j} = \frac{T^k_{i-1,j}+T^k_{i+1,j}+T^k_{i,j-1}+T^k_{i,j+1}}{4}, \; \forall i,j \in \{1,\ldots,N-2\},\; \forall k \in \{0,\ldots,n\}\]
et qui sera celle à implémenter sur le GPU dans un premier temps, car elle est naturellement parallélisable;
- une deuxième qui utilise une itération de Gauss-Seidel: on utilise les valeurs nouvellement calculées (donc pour le même k) pour accélérer la convergence de la méthode; elle sera la version séquentielle de réference pour
calculer l'accélération due à l'utilisation du GPU;
- enfin, on modifiera chacune de ces deux versions pour s'arrêter non plus après un nombre fixé d'itérations, mais lorsque le calcul a convergé (la différence max entre la température aux pas k et k+1 est inférieure à un seuil donné)
Une possibilité pour l'affichage
Pour l'affichage des résultats, vous pouvez par exemple:
- insérer votre code dans le fichier suivant Chaleur.cu, en mettant également le fichier pixmap_io.h dans le répertoire; la fonction main vous montre comment afficher le tableau T;
- l'exécution produit un fichier de résultats out.pnm, que vous pouvez afficher par exemple par gimp out.pnm (ou bien, en salle machine, simplement en double cliquant sur le fichier)
- dans un second temps, pour éventuellement visualiser l'évolution de la température au cours du temps, dans le cas de l'équation de la chaleur,
vous pourrez sauver plusieurs fichiers à différents instants et utiliser ImageMagick pour l'affichage
Versions parallèles
- Implémenter et comparer les performances en temps de calcul des versions parallèles suivantes de l'algo précédent:
- Version pour un nombre fixé d'itérations, et dans laquelle tous les threads sont synchronisés (par le CPU)
à chaque itération k, version couteuse, qui nécessite de relancer un "kernel" à chaque itération
- Si on relache la contrainte de synchronisation, on se retrouve avec une itération de style Gauss-Seidl beaucoup plus efficace mais suivant les cas elle peut ne pas converger, expérimenter des versions de ce type
- Version avec critère d'arrêt sur la convergence du calcul, et non plus nombre d'itérations fixé: comment faire efficacement et correct?
Quelques suggestions:
- Pour chaque version, bien vérifier que les résultats sont bien comparables sur GPU et CPU, et comparer les temps d'exécution lorsque l'on fait varier les tailles des blocs et grilles, ainsi que le nombre de points N de discrétisation,
- Vous pouvez essayer de tirer parti de la mémoire shared, des fonctions atomic
Equation de la chaleur
On s'intéresse maintenant à l'évolution de la température dans la pièce, modélisée par l'équation de la chaleur
\[ \frac{\partial T}{\partial t} = \alpha \Delta T, \]
en partant par exemple d'une température égale à 20 degrés partout sauf sur un mur ou elle est égale à 100,
et avec comme conditions aux limites 20 degrés sur tous les murs sauf celui ou elle vaut 100.
Calculer la température en fonction du temps par une méthode de différences finies. Cette fois, la
synchronisation stricte est nécessaire. Comme précédemment, implémenter et comparer plusieurs versions
en essayant d'améliorer progressivement les performances:
- synchronisation par le CPU comme précédemment,
- synchronisation par blocs et inter blocs en utilisant les fonctions atomic,
- mémoire shared