CUDA approfondi


précédentsommairesuivant

X. CUDPP

X-A. Qu'est-ce que CUDPP ?

CUDA Data Parallel Primitives.

Il s'agit d'une autre API C. Elle est donc orientée vers les algorithmes de traitement parallèle des données, comme la réduction parallèle et le tri parallèle des listes. Ces algorithmes, primitifs, sont nécessaires à bien d'autres algorithmes, comme la compression de données, sans être forcément très faciles à implémenter efficacement.

CUDPP se base, comme CuFFT et FFTW pour les transformées de Fourier, sur le principe des plans.

La majorité des développeurs utilise un framework (Qt, GTK+), ou une librairie pluripotente (Boost) pour s'aider dans leurs développements. Ces librairies permettent d'effectuer rapidement, efficacement et simplement des pans entiers de nos programmes, comme le tri. CUDPP effectue les mêmes opérations, mais sur un GPU, à la manière de CuBLAS et de CuFFT : l'utilisateur ne voit presque pas la différence entre GPU et CPU.

Cette librairie est aussi disponible pour Java.

X-B. Compilation

CUDPP n'est pas développé par NVIDIA, ni supporté, par conséquent, non distribué avec le SDK de CUDA. Il vous faudra donc le télécharger puis le compiler. Cette opération prend un peu de temps, environ 15 minutes sur une configuration assez musclée.

Cependant, CUDPP dépend de CUTIL, une librairie développée par NVIDIA, mais non supportée, qui est censée aider au développement avec CUDA. Cette librairie est incluse et doit être compilée avant utilisation.

Commencez par récupérer la librairie sur le site officiel : GPGPU CUDPP. La version actuelle est la 1.0 alpha.

Une fois l'archive récupérée, décompressez-la dans un répertoire, peu importe où, tant qu'il n'y a pas d'espace dans le chemin.

La suite dépend de votre environnement. Vous devrez avoir installé CUDA correctement et configuré l'environnement pour son utilisation.

CUDPP est censé compiler avec CUDA 2.0 et être incompatible avec les autres versions. Or, j'ai pu le compiler sans problème avec CUDA 2.2. Si vous avez des problèmes, il vous suffit de revenir à une version plus ancienne de CUDA.

X-B-1. Visual Studio

Dans le dossier common, ouvrez les solutions cutil.sln et paramgl.sln. Compilez les configurations Debug et Release.

Ensuite, dans le dossier cudpp, ouvrez la solution cudpp.sln. Sélectionnez la configuration qui vous intéresse (un conseil : compilez au moins une version Debug et une version Release, pour faciliter le débogage de votre application), puis lancez la compilation, comme pour toute solution.

Ensuite, il faut installer ces librairies. Copiez les fichiers .lib du dossier lib dans votre dossier CUDA\lib, les .h du dossier inc dans CUDA\include et les fichiers .dll du dossier lib dans votre dossier CUDA\bin.

Si vous n'arrivez pas à compiler correctement ces librairies, téléchargez cette archive, qui les contient toutes (vous n'avez qu'à l'extraire dans le dossier Visual Studio 9.0\VC). Les fichiers compilés sont prévus pour Visual Studio 9 et CUDA 2.2, toute utilisation avec d'autres versions n'est pas garantie de fonctionnement.

L'utilisation des compilateurs Intel est déconseillée dans ce cas, car nvcc n'a pas encore été testé avec lui ! Le résultat de la compilation peut être très aléatoire !

Pendant la compilation, vous verrez beaucoup d'avertissements : ils ne sont que les signes de l'optimiseur du compilateur, qui optimise tout ce qu'il peut, surtout du côté GPU.

X-B-2. GCC

Dans un shell Linux (sh, par exemple), lancez ces quelques commandes.

 
Sélectionnez
( cd common ; make cuda-install=/usr/local/cuda )
( cd cudpp  ; make cuda-install=/usr/local/cuda )

Elles compilent et installent le tout automatiquement, dans le dossier de CUDA.

X-C. La découverte par l'exemple

Cet exemple est celui livré avec la librairie CUDPP, sous le nom de simpleCUDPP.

 
TéléchargerCacherSélectionnez

L'entièreté du fichier est disponible ici et les fonctions utilisées non définies dedans le sont dans ce fichier.

X-D. Décryptage de l'exemple

X-D-1. Initialisation

Tout d'abord, qu'est-il censé effectuer ? Un scan, aussi dénommé somme partielle, somme des préfixes ou réduction des préfixes. Mais qu'est-ce donc ?

Prenez une liste. Appliquez-lui un scan. Il en résultera une autre liste, dont chaque élément est la somme des précédents.

Liste d'origine
Sélectionnez
1
2
4
8
16
14
2004
0
Liste scannée
Sélectionnez
1
3
7
15
31
45
2049
2049
 
TéléchargerSélectionnez
2.
3.
4.
5.
6.
7.
8.
9.
10.
11.
12.
13.
CUT_DEVICE_INIT();

unsigned int numElements = 32768;
unsigned int memSize = sizeof(float) * numElements;


float* h_idata = (float *) malloc(memSize);

for (unsigned int i = 0; i < numElements; ++i) 
{
    h_idata[i] = (float) (rand() & 0xf);
}

CUTIL initialise un périphérique, on spécifie le nombre d'éléments à traiter puis la mémoire utilisée pour stockerces éléments. On initialise le tableau qui va les contenir et on le remplit de données aléatoires (entre 0 et 15).

 
TéléchargerSélectionnez
17.
18.
19.
20.
21.
22.
23.
24.
25.
float* d_idata;
CUDA_SAFE_CALL( cudaMalloc( (void **) & d_idata, memSize));

CUDA_SAFE_CALL( cudaMemcpy( d_idata, h_idata, memSize,
                            cudaMemcpyHostToDevice) );


float* d_odata;
CUDA_SAFE_CALL( cudaMalloc( (void **) & d_odata, memSize));

On crée un tableau semblable sur le GPU et on y transfère les données. On prépare aussi un emplacement pour le résultat.

Maintenant, nous pouvons enfin commencer à parler de CUDPP.

X-D-2. Les plans

Nous devons commencer par configurer CUDPP correctement pour scanner notre tableau. La configuration des algorithmes de CUDPP repose sur le concept des plans.

Un plan est, tout comme dans CUFFT, une structure qui contient les informations sur les algorithmes à employer pour le calcul.

Ici, nous allons effectuer un scan sur des flottants, exclusivement vers l'avant, d'un nombre d'éléments connu d'avance.

 
TéléchargerSélectionnez
27.
28.
29.
30.
31.
32.
33.
34.
35.
36.
37.
38.
39.
40.
CUDPPConfiguration config;
config.op = CUDPP_ADD;
config.datatype = CUDPP_FLOAT;
config.algorithm = CUDPP_SCAN;
config.options = CUDPP_OPTION_FORWARD | CUDPP_OPTION_EXCLUSIVE;

CUDPPHandle scanplan = 0;
CUDPPResult result = cudppPlan(&scanplan, config, numElements, 1, 0);

if (CUDPP_SUCCESS != result)
{
    printf("Error creating CUDPPPlan\n");
    exit(-1);
}

X-D-3. Le calcul

 
TéléchargerSélectionnez
43.
cudppScan(scanplan, d_odata, d_idata, numElements);

CUDPP effectue le scan, enfin.

 
TéléchargerSélectionnez
46.
47.
48.
49.
50.
51.
52.
53.
54.
55.
56.
float* h_odata = (float*) malloc( memSize);  

CUDA_SAFE_CALL( cudaMemcpy( h_odata, d_odata, memSize,  
                            cudaMemcpyDeviceToHost) );  

float* reference = (float*) malloc( memSize);  
computeSumScanGold( reference, h_idata, numElements, config);  


CUTBoolean res = cutComparef( reference, h_odata, numElements);  
printf( "Test %s\n", (1 == res) ? "PASSED" : "FAILED");  

On récupère le résultat du calcul, on effectue le même calcul en local (sur le CPU) puis on compare les résultats : s'ils sont identiques, le test est réussi.

 
TéléchargerSélectionnez
58.
59.
60.
61.
62.
63.
64.
65.
66.
67.
68.
69.
70.
result = cudppDestroyPlan(scanplan);  
if (CUDPP_SUCCESS != result)  
{  
    printf("Error destroying CUDPPPlan\n");  
    exit(-1);  
}  

free(h_idata);  
free(h_odata);  
free(reference);  
CUDA_SAFE_CALL(cudaFree(d_idata));  
CUDA_SAFE_CALL(cudaFree(d_odata));  
} 

On détruit le plan, puis on libère la mémoire allouée.

Maintenant, vous pouvez voir la simplicité de CUDPP : on alloue de la mémoire, on la remplit, on configure CUDPP en lui disant ce qu'il doit faire, puis il le fait et on récupère le résultat. C'est tout. Et c'est optimisé.

X-E. Gestion des erreurs

Chaque fonction de CUDPP renvoie un CUDPPResult, qui peut prendre quelques valeurs, en fonction de l'erreur qui s'est déroulée.

Code Signification
CUDPP_SUCCESS Pas d'erreur.
CUDPP_ERROR_INVALID_HANDLE CUDPPHandle spécifié invalide.
CUDPP_ERROR_ILLEGAL_CONFIGURATION La configuration demandée est impossible (illogique ou invalide).
CUDPP_ERROR_UNKNOWN Erreur inconnue ou impossible à tracer.

X-F. Options de l'algorithme

Tout algorithme de CUDPP doit être préalablement configuré, grâce à une structure de type CUDPPConfiguration. Cette structure est ainsi définie.

Définition de CUDPPConfiguration
Sélectionnez
struct CUDPPConfiguration
{
    CUDPPAlgorithm algorithm;
    CUDPPOperator  op;
    CUDPPDatatype  datatype;
    unsigned int   options;
};

Nous allons étudier les valeurs de chacune de ces variables.

X-F-1. Algorithme

Premièrement, nous utilisons une variable de type CUDPPAlgorithm pour configurer un algorithme. Il s'agit, justement, de l'algorithme à employer.

  • CUDPP_SCAN
  • CUDPP_SEGMENTED_SCAN
  • CUDPP_COMPACT
  • CUDPP_REDUCE
  • CUDPP_SORT_RADIX
  • CUDPP_SORT_RADIX_GLOBAL
  • CUDPP_SPMVMULT
  • CUDPP_SORT_INVALID

Évidemment, cet algorithme devra correspondre à la fonction appelée plus tard avec cette configuration.

X-F-2. Opérateur

Deuxièmement, vous pouvez spécifier un opérateur, l'opération qui sera effectuée, uniquement lors d'un scan. Pour les autres algorithmes, laissez ceci vide.

Tous ces opérateurs sont binairement associables.

  • CUDPP_ADD
  • CUDPP_MULTIPLY
  • CUDPP_MIN
  • CUDPP_MAX

X-F-3. Type

Troisièmement, vous devez définir le type des données en entrée et en sortie. Voici un tableau récapitulatif de ces valeurs, avec leur équivalent en C.

CUDPPDatatype C
CUDPP_CHAR char
CUDPP_UCHAR unsigned char
CUDPP_INT int
CUDPP_UINT unsigned int
CUDPP_FLOAT float

X-F-4. Configuration

La configuration d'un algorithme se fait, finalement, aussi grâce à un entier. Cet entier est constitué d'une combinaison de paramètres, grâce à l'opérateur pipe (|).

Nom Action Valeur
CUDPP_OPTION_FORWARD L'algorithme se déroulera du début vers la fin des données. 0x1
CUDPP_OPTION_BACKWARD L'algorithme se déroulera de la fin vers le début des données. 0x2
CUDPP_OPTION_EXCLUSIVE L'algorithme inclura tous les éléments précédents sauf l'actuel (uniquement pour des scans). 0x4
CUDPP_OPTION_INCLUSIVE L'algorithme inclura tous les éléments précédents et l'actuel (uniquement pour des scans). 0x8
CUDPP_OPTION_CTA_LOCAL L'algorithme sera effectué sans communication entre les blocs (ne fonctionne que pour les tris, pour le moment). 0x10

X-G. Fonctions de plans

 
Sélectionnez
CUDPPResult cudppPlan(CUDPPHandle        * planHandle, 
                      CUDPPConfiguration config, 
                      size_t             n, 
                      size_t             rows, 
                      size_t             rowPitch);
                      
CUDPPResult cudppDestroyPlan(CUDPPHandle plan);

cudppPlan() crée un plan planHandler, qui respectera la configuration config, qui traitera n éléments répartis en rows colonnes de rowPitch éléments chacune.

Le plan sera détruit par la fonction cudppDesroyPlan() dès qu'il ne sera plus utilisé.

X-H. Fonctions de calcul

X-H-1. Scan

 
Sélectionnez
CUDPPResult cudppScan(CUDPPHandle   planHandle,
                      void        * d_out, 
                      const void  * d_in, 
                      size_t        numElements);

Les données à traiter se situent sur le GPU, à l'emplacement pointé par d_in ; la sortie se situe aussi sur le GPU, à l'emplacement pointé par d_out. La fonction utilisera les options prévues dans le plan planHandle et scannera numElements éléments.

 
Sélectionnez
CUDPPResult cudppSegmentedScan(CUDPPHandle          planHandle,
                               void               * d_out, 
                               const void         * d_idata,
                               const unsigned int * d_iflags,
                               size_t               numElements);

La fonction traitera les segments précisés dans d_iflags de d_idata.

 
Sélectionnez
CUDPPResult cudppMultiScan(CUDPPHandle planHandle,
                           void        *d_out, 
                           const void  *d_in, 
                           size_t      numElements,
                           size_t      numRows);

La fonction effectuera numRows scans sur les données avant d'écrire le résultat.

 
Sélectionnez
CUDPPResult cudppCompact(CUDPPHandle          planHandle,
                         void               * d_out, 
                         size_t             * d_numValidElements,
                         const void         * d_in, 
                         const unsigned int * d_isValid,
                         size_t               numElements);

Cette fonction compacte un tableau, avec les éléments d_isValid définis comme valides, les autres étant à éjecter. d_numValidElements sera le nombre d'éléments effectivement gardés.

Par exemple :

Exemple de compactage
Sélectionnez
d_in        = [ a b c d e f ]
deviceValid = [ 1 0 1 1 0 1 ]
d_out       = [ a c d f ]

X-H-2. Tri

 
Sélectionnez

CUDPPResult cudppSort(CUDPPHandle planHandle,
                      void        *d_out, 
                      const void  *d_in,
                      size_t      numElements);

Trie le tableau.

XI. Thrust

XI-A. Qu'est-ce que Thrust ?

Il s'agit d'une bibliothèque de templates à utiliser avec CUDA. En tant que telle, elle ne nécessite pas de compilation : pour l'installer, téléchargez la dernière version sur son site web et extrayez l'archive dans votre dossier d'includes de CUDA.

Thrust dispose quand même de fonctions : des fonctions inline. Cela permet de gagner un certain temps à l'exécution : un appel de fonction produit un remplissage des registres du CPU puis un saut vers la fonction, qui resautera vers l'emplacement précédent, ce qui prend un certain temps. En contre partie, votre exécutable sera plus lourd (comme si vous aviez inclus une nouvelle librairie en statique).

Thrust se base sur la STL pour fournir ses fonctions : il arrive assez fréquemment que des fonctions aient le même nom. C'est pourquoi toutes les fonctions de Thrust sont dans le namespace thrust.

XI-B. Un premier exemple : la version

 
TéléchargerCacherSélectionnez

Ce code montre comment récupérer la version de Thrust, tout simplement.

XI-C. Les vecteurs

 
TéléchargerCacherSélectionnez

Le terme vecteur est ici à comprendre au sens de la STL et non des mathématiques : il ne s'agit pas d'un ensemble de nombres, mais bien d'un tableau dynamique.

Il existe deux types de vecteurs : l'un pour l'hôte (host_vector), l'autre pour le périphérique (device_vector).

Les copies de l'hôte vers le périphérique sont simplifiées à l'extrême : on initialise un vecteur de chaque type, et on affecte la valeur de l'un à l'autre. C'est tout.

Chaque modification d'un vecteur sur le périphérique fait appel à cudaMemcpy(), ils sont donc fort coûteux. Il faut les limiter au maximum.

 
TéléchargerCacherSélectionnez

Ici, on initialise un vecteur sur le périphérique, on remplit les 7 premiers éléments avec 9. On crée un autre vecteur, sur l'hôte, qui contient les 5 premières valeurs du premier vecteur. L'hôte est rempli avec des valeurs incrémentant à partir de 1. On copie l'hôte de H.begin() à H.end() sur D.begin()

XI-D. Itérateurs et pointeurs

Vous n'êtes probablement pas sans savoir que les fonctions begin() et end() que nous venons d'utiliser sont des itérateurs : ils reprennent le nom utilisé pour la STL. Ces fonctions sont implémentées comme des templates : une seule et unique implémentation pour gérer l'hôte et le périphérique. Leurs types de renvoi sont quand même différents selon le milieu où elles doivent être exécutées.

Vous pouvez aussi utiliser des pointeurs créés directement par CUDA avec Thrust.

 
TéléchargerSélectionnez
size_t N = 10;

int * raw_ptr;
cudaMalloc( (void **) & raw_ptr, N * sizeof(int));

thrust::device_ptr<int> dev_ptr(raw_ptr);

thrust::fill(dev_ptr, dev_ptr + N, (int) 0);

Vous pouvez aussi récupérer un pointeur CUDA depuis un pointeur Thrust.

 
Sélectionnez
size_t N = 10;

thrust::device_ptr<int> dev_ptr = thrust::device_malloc<int>(N);

int * raw_ptr = thrust::raw_pointer_cast(dev_ptr);

Thrust peut aussi être compatible avec les conteneurs de la STL, quand il n'en dispose pas.

 
TéléchargerCacherSélectionnez

XI-E. Algorithmes

Thrust dispose d'un certain nombre d'algorithmes parallèles, dont une certaine majorité possède un équivalent dans la STL. Les deux versions (STL et Thrust) ont alors le même nom.

Chacun des algorithmes de Thrust est implémenté pour l'hôte et pour le périphérique. Suivant le conteneur ou itérateur qui est passé à l'algorithme, la version hôte ou périphérique sera utilisée.

La majorité des algorithmes de Thrust requièrent que les itérateurs soient dans le même milieu (hôte ou périphérique), à quelques exceptions logiques près, comme copy().

XI-E-1. Transformations

Ces algorithmes appliquent une même opération à tous les éléments d'une entrée et les stockent à un autre endroit.

 
TéléchargerCacherSélectionnez

thrust::modulus et thrust::negate sont implémentés comme des foncteurs, dans <thrust/functional.h>.

Nous allons poursuivre avec une implémentation personnelle d'une fonction de BLAS : SAXPY, qui effectue les opérations de la forme y = a * x + y.

 
TéléchargerCacherSélectionnez

Vous pouvez voir deux implémentations : une rapide et une lente. La première effectue 2x lectures et x écritures. La seconde effectue 4x lectures et 3x écritures.

Vu que SAXPY est limité par la mémoire et non par la puissance de calcul, la première est plus rapide que la seconde.

La technique utilisée s'appelle la réduction des kernels, car, pour obtenir de meilleures performances, dans ce cas, vous devez effectuer le moins possible de lectures et d'écritures en mémoire.

Cette manière de faire (avec thrust::transform) limite les arguments : un ou deux en entrée. Si vous avez besoin de plus d'arguments, voyez thrust::counting_iterator et thrust::for_each tels qu'utilisés dans ce nouvel exemple.

 
TéléchargerCacherSélectionnez

XI-E-2. Réductions

Un algorithme de réduction utilise une opération binaire pour réduire une entrée à une valeur. Voici la manière d'obtenir la somme d'un tableau avec Thrust.

 
Sélectionnez
int sum = thrust::reduce(D.begin(), D.end(), (int) 0, thrust::plus<int>());

Cette manière de procéder convient généralement, mais pas toujours, par exemple pour compter. Thrust, comme la STL, fournit quelques autres fonctions, pour votre facilité.

 
Sélectionnez
#include <thrust/count.h>
#include <thrust/device_vector.h>

thrust::device_vector<int> vec(5,0);
vec[1] = 1;
vec[3] = 1;
vec[4] = 1;

int result = thrust::count(vec.begin(), vec.end(), 1);
// result == 3

On peut utiliser ces réductions pour calculer la norme d'un vecteur.

 
TéléchargerCacherSélectionnez

En fusionnant la multiplication avec l'addition, nous avons un kernel très optimisé, autant qu'un kernel écrit à la main et patiemment optimisé.

XI-E-3. Scans

 
Sélectionnez
#include <thrust/scan.h>
int data[6] = {1, 0, 2, 2, 1, 3};
thrust::inclusive_scan(data, data + 6, data);
// data == {1, 1, 3, 5, 6, 9}
 
Sélectionnez
#include <thrust/scan.h>
int data[6] = {1, 0, 2, 2, 1, 3};
thrust::exclusive_scan(data, data + 6, data);
// data == {0, 1, 1, 3, 5, 6}

XI-E-4. Tri

Thrust, une fois de plus, se calque sur la STL pour fournir ses algorithmes de tri.

 
TéléchargerCacherSélectionnez

De plus, Thrust propose le tri par clé.

 
Sélectionnez
const int N = 6;
int    keys[N] = {  1,   4,   2,   8,   5,   7};
char values[N] = {'a', 'b', 'c', 'd', 'e', 'f'};
thrust::sort_by_key(keys, keys + N, values);
// keys   == {  1,   2,   4,   5,   7,   8}
// values == {'a', 'c', 'b', 'e', 'f', 'd'}

Et selon un critère défini par l'utilisateur.

 
Sélectionnez
const int N = 6;
int A[N] = {1, 4, 2, 8, 5, 7};
thrust::stable_sort(A, A + N, thrust::greater<int>());
// A == {8, 7, 5, 4, 2, 1}

précédentsommairesuivant

Vous avez aimé ce tutoriel ? Alors partagez-le en cliquant sur les boutons suivants : Viadeo Twitter Facebook Share on Google+   

CUDA et le GPGPU
Introduction à CUDA
CUDA approfondi
La FAQ GPGPU & CUDA
  

Copyright © 2009 Thibaut Cuvelier. Aucune reproduction, même partielle, ne peut être faite de ce site et de l'ensemble de son contenu : textes, documents, images, etc. sans l'autorisation expresse de l'auteur. Sinon vous encourez selon la loi jusqu'à trois ans de prison et jusqu'à 300 000 € de dommages et intérêts.