Mise en œuvre de votre GPU avec CuPy et PyCUDA

Dans cette seconde partie, nous décou­vrons comment utili­ser votre GPU avec la librai­rie CuPy, qui est l’équi­valent de NumPy et vous permet­tra de calcu­ler des matrices NumPy sur votre GPU sans chan­ger plus d’une ligne à votre code !

Puis; nous rentre­rons plus en profon­deur dans l’écri­ture de code natif avec PyCUDA. Je vous rassure, dans la vraie vie on peut s’en passer. Mais pour comprendre les fonde­ments, c’est une étape appro­priée.

Un exemple avec CuPy

On ne peut pas se lancer dans le calcul sur GPU sans comprendre un mini­mum l’or­ga­ni­sa­tion logi­cielle/maté­rielle de ces systèmes. Ils ne se programment pas comme un CPU et, par exemple, avec numba, votre code Python devra être adapté en consé­quence.

Mais, comme vous êtes certai­ne­ment très impa­tients de toucher du doigt ce nouveau monde, voici un petit exemple avec CuPy qui va nous permettre de décou­vrir à la fois la simpli­cité de la librai­rie et la puis­sance de votre carte graphique.

Nous vous propo­sons de dessi­ner la frac­tale de l’en­semble de Mandel­brot.

Pour cela nous allons repar­tir de l’exemple de matplot­lib dont le code présente l’avan­tage de bien se prêter à l’exer­cice.

L’ins­tal­la­tion de la librai­rie CuPy est très simple via pip :

$ pip install cupy-<cuda version>

Dans notre exemple avec CUDA 12.6, ce sera :

$ pip install cupy-cuda12x

Si vous utili­sez Mini­forge/Conda/Mamba :

$ mamba install cupy

CuPy implé­mente globa­le­ment la même API que NumPy. On peut donc rempla­cer la librai­rie NumPy par CuPy et le tour est joué:

En remplaçant le seul import :

import numpy as np

Par ce dernier :

import cupy as np

Votre programme exécu­tera ses calculs sur GPU sans autre effort.

Mais cela c’est la théo­rie dans un monde parfait. La pratique n’est pas aussi simple si l’on veut des résul­tats à la hauteur des capa­ci­tés de votre GPU.

Dans l’exemple proposé ci-dessous, 2 fonc­tions mandel­brot_set ont été défi­nies.

• mandelbrot_set_np : corres­pond à l’im­plé­men­ta­tion avec NumPy du calcul de la frac­tale
• mandelbrot_set_cp : corres­pond à l’im­plé­men­ta­tion avec CuPy du calcul de la frac­tale

Voici, les chan­ge­ments opérés entre les 2 versions de la fonc­tion par rapport au code origi­nal de l’exemple Matplot­lib :

• Dans mandelbrot_set_cp le module np a été remplacé par cp
• Dans les 2 versions le troi­sième para­mètre de linspace a été converti en entier, int(xn):
  • Cela évite les warning NumPy
  • Cela évite une erreur CuPy (TypeError: 'float' object cannot be interpreted as an integer)
• La fonc­tion abs(Z) a été rempla­cée par cp/np.abs(Z)
• La fonc­tion CuPy doit en plus conver­tir son résul­tat en tableau NumPy pour l’af­fi­chage final.
  C’est le trans­fert des données stockées dans la RAM du GPU vers celle du CPU, sinon matplot­lib ne saura pas affi­cher l’image géné­rée.
• Enfin, dans la version CuPy, nous avons ajouté un bloc with permet­tant de choi­sir le GPU sur lequel le calcul sera exécuté lorsque vous en possé­dez plusieurs

La fonc­tion compute prend quant à elle en para­mètre la fonc­tion mandelbrot_set qui sera utili­sée pour le calcul de la frac­tale (avec en option le numéro du GPU pour la version CuPy) et exécute le calcul de la frac­tale puis l’af­fiche en indiquant le temps de calcul mesuré.

import numpy as np
import time 
import matplotlib
from matplotlib import colors
import matplotlib.pyplot as plt
from timeit import default_timer as timer
import cupy as cp

def mandelbrot_set_np(xmin, xmax, ymin, ymax, xn, yn, maxiter, horizon=2.0, dummy=None):
    X = np.linspace(xmin, xmax, int(xn), dtype=np.float32)
    Y = np.linspace(ymin, ymax, int(yn), dtype=np.float32)
    C = X + Y[:, None]*1j
    N = np.zeros(C.shape, dtype=np.int64)
    Z = np.zeros(C.shape, np.complex64)
    for n in range(maxiter):
        I = np.less(np.abs(Z), horizon)
        N[I] = n
        Z[I] = Z[I]**2 + C[I]
    N[N == maxiter-1] = 0
    return Z, N

def mandelbrot_set_cp(xmin, xmax, ymin, ymax, xn, yn, maxiter, horizon=2.0, gpu_number=0):
    with cp.cuda.Device(gpu_number):  # Choosing the GPU to use
        # Just renaming NP with CP
        X = cp.linspace(xmin, xmax, int(xn), dtype=cp.float32)
        Y = cp.linspace(ymin, ymax, int(yn), dtype=cp.float32)
        C = X + Y[:, None]*1j
        N = cp.zeros(C.shape, dtype=cp.int64)
        Z = cp.zeros(C.shape, cp.complex64)
        for n in range(maxiter):
            I = cp.less(cp.abs(Z), horizon)
            N[I] = n
            Z[I] = Z[I]**2 + C[I]
        N[N == maxiter-1] = 0

        return cp.asnumpy(Z), cp.asnumpy(N)  # Convert CP tables (stored in GPU to NP table in CPU RAM)


def compute(compute_set, gpu=0):

    xmin, xmax, xn = -2.25, +0.75, 3000/2
    ymin, ymax, yn = -1.25, +1.25, 2500/2
    maxiter = 400
    horizon = 2.0 ** 40
    log_horizon = np.log(np.log(horizon))/np.log(2)

    start = timer()
    Z, N = compute_set(xmin, xmax, ymin, ymax, xn, yn, maxiter, horizon, gpu)
    stop = timer()

    msg = "Rendered in %08.6f seconds\n" % (stop - start)
    print(msg)

    # Normalized recount as explained in:
    # https://linas.org/art-gallery/escape/smooth.html
    # https://www.ibm.com/developerworks/community/blogs/jfp/entry/My_Christmas_Gift

    # This line will generate warnings for null values but it is faster to
    # process them afterwards using the nan_to_num
    with np.errstate(invalid='ignore'):
        M = np.nan_to_num(N + 1 -
                          np.log(np.log(abs(Z)))/np.log(2) +
                          log_horizon)

    dpi = 72
    width = 4
    height = 4*yn/xn
    fig = plt.figure(figsize=(width, height), dpi=dpi)
    ax = fig.add_axes([0.0, 0.0, 1.0, 1.0], frameon=False, aspect=1)

    # Shaded rendering
    light = colors.LightSource(azdeg=315, altdeg=10)
    M = light.shade(M, cmap=plt.cm.hot, vert_exag=1.5,
                    norm=colors.PowerNorm(0.3), blend_mode='hsv')
    plt.imshow(M, extent=[xmin, xmax, ymin, ymax], interpolation="bicubic")
    ax.set_xticks([])
    ax.set_yticks([])

    # Some advertisement for matplotlib
    year = time.strftime("%Y")
    major, minor, micro = matplotlib.__version__.split('.', 2)
    text = ("The Mandelbrot fractal set\n" +
            msg +
            "Rendered with matplotlib %s.%s, %s - http://matplotlib.org"
            % (major, minor, year))
    ax.text(xmin+.025, ymin+.025, text, color="white", fontsize=12, alpha=0.5)

    plt.show()

compute(mandelbrot_set_np, gpu=None)

Rendered in 2.401172 seconds

Le premier appel à CuPy ne sera pas mira­cu­leux, car il intègre le temps de compi­la­tion du code sur le GPU, cela est assez couteux…

compute(mandelbrot_set_cp, gpu=0)  

Rendered in 0.647558 seconds

Si vous le relan­cez, le code étant déjà compilé, le temps sera meilleur.

compute(mandelbrot_set_cp, gpu=0)

Rendered in 0.374514 seconds

print("Ratio CPU/GPU:", 2.5/0.37)

Ratio CPU/GPU: 6.756756756756757

Avec très peu de modi­fi­ca­tions, nous avons exécuté le calcul sur GPU. Le véri­table chan­ge­ment consiste à avoir renommé numpy en cupy.

C’est trop facile avec Python !

La carte graphique offre un gain d’en­vi­ron x6/7 par rapport à la version NumPy sur CPU. C’est déjà bien.

Mais c’est un succès très mitigé si nous regar­dons le nombre total d’opé­ra­tions que chaque proces­seur peut réali­ser à la seconde.

Compa­rons les carac­té­ris­tiques des 2 proces­seurs pour connaître leurs perfor­mances.

!lscpu | grep -v Drapaux | grep -v Vulnerability

Architecture :                          x86_64
Mode(s) opératoire(s) des processeurs : 32-bit, 64-bit
Address sizes:                          39 bits physical, 48 bits virtual
Boutisme :                              Little Endian
Processeur(s) :                         12
Liste de processeur(s) en ligne :       0-11
Identifiant constructeur :              GenuineIntel
Nom de modèle :                         Intel(R) Core(TM) i7-8700 CPU @ 3.20GHz
Famille de processeur :                 6
Modèle :                                158
Thread(s) par cœur :                    2
Cœur(s) par socket :                    6
Socket(s) :                             1
Révision :                              10
CPU(s) scaling MHz:                     93%
Vitesse maximale du processeur en MHz : 4600,0000
Vitesse minimale du processeur en MHz : 800,0000
BogoMIPS :                              6399,96
Virtualisation :                        VT-x
Cache L1d :                             192 KiB (6 instances)
Cache L1i :                             192 KiB (6 instances)
Cache L2 :                              1,5 MiB (6 instances)
Cache L3 :                              12 MiB (1 instance)
Nœud(s) NUMA :                          1
Nœud NUMA 0 de processeur(s) :          0-11

Le CPU utilisé est un proces­seur Intel i7 de huitième géné­ra­tion :

• Il possède une fréquence de 4,6Ghz.
• Il dispose de 6 cœurs physiques et 12 hyper­threads, soit 2 hyper­threads par cœur.

Les hyper­threads sont des cœurs virtuels/logiques. Tous les proces­seurs n’en possèdent pas. Quand un proces­seur en est doté, un seul hyper­thread fonc­tionne sur le cœur physique à la fois, comme évoqué dans le chapitre de présen­ta­tion des concepts du paral­lé­lisme.
Les hyper­threads sont utili­sés unique­ment quand le proces­seur réalise des opéra­tions d’en­trée/sortie : plutôt que de perdre de précieuses milli­se­condes à patien­ter pour attendre que le disque dur ait fini d’écrire ce qu’il a en cours et lui envoyer la suite, le proces­seur bascule sur le second hyper­thread et peut ainsi exécu­ter d’autres tâches pendant cette attente.

Dans le cas d’un calcul numé­rique comme ici, il n’y a pas de véri­table opéra­tion d’en­trée/sortie, les hyper­threads ne sont pas utili­sés.

La fréquence repré­sente le nombre de cycles qu’un proces­seur peut réali­ser en 1 seconde. Ici 4,6 milliards.

Un cycle est le temps qu’il faut au proces­seur pour chan­ger l’état d’un bit. Dans un ordi­na­teur, en géné­ral un bit est repré­senté par un tran­sis­tor possé­dant une tension de 0 ou 5 volts. Si sa tension est nulle, il vaut 0, sinon, il vaut 1. Quand la valeur du bit passe de 0 à 1, sa tension croit progres­si­ve­ment jusqu’à 5 volts. Tant que la tension lue n’est pas égale soit à 0 soit à 5, son état est instable et sa valeur indé­ter­mi­née. Il ne peut pas être utilisé pour le calcul.

La fréquence repré­sente donc le nombre de chan­ge­ments de bits que le proces­seur peut faire à la seconde.

Certaines opéra­tions ne demandent qu’un cycle machine, comme inver­ser l’état des bits d’un registre. D’autres peuvent deman­der plusieurs cycles. Comme char­ger un registre avec des données issues de la RAM.

Pour simpli­fier, disons qu’un cycle repré­sente 1 instruc­tion du proces­seur.

Nous avons donc un CPU pouvant exécu­ter : 27.6 milliards d’ins­truc­tions à la seconde

nb_coeurs = 6
frequence = 4.6 #  en Giga Hertz

nombre_instructions_cpu = nb_coeurs * frequence
print("Nombre d'instructions à la seconde du CPU en milliards:", nombre_instructions_cpu)

Nombre d'instructions à la seconde du CPU en milliards : 27.599999999999998

Regar­dons main­te­nant du côté de la carte graphique :

!nvidia-smi

Mon Jan 27 15:09:44 2025       
+-----------------------------------------------------------------------------------------+
| NVIDIA-SMI 560.35.05              Driver Version: 560.35.05      CUDA Version: 12.6     |
|-----------------------------------------+------------------------+----------------------+
| GPU  Name                 Persistence-M | Bus-Id          Disp.A | Volatile Uncorr. ECC |
| Fan  Temp   Perf          Pwr:Usage/Cap |           Memory-Usage | GPU-Util  Compute M. |
|                                         |                        |               MIG M. |
|=========================================+========================+======================|
|   0  NVIDIA GeForce RTX 2070        Off |   00000000:01:00.0  On |                  N/A |
| N/A   48C    P0             31W /  115W |    1114MiB /   8192MiB |      8%      Default |
|                                         |                        |                  N/A |
+-----------------------------------------+------------------------+----------------------+

+-----------------------------------------------------------------------------------------+
| Processes:                                                                              |
|  GPU   GI   CI        PID   Type   Process name                              GPU Memory |
|        ID   ID                                                               Usage      |
|=========================================================================================|
|    0   N/A  N/A      1722      G   /usr/lib/xorg/Xorg                            503MiB |
|    0   N/A  N/A      2762      G   cinnamon                                       41MiB |
|    0   N/A  N/A      3392      G   /usr/bin/nextcloud                              1MiB |
|    0   N/A  N/A      3903      G   /usr/lib/firefox/firefox                        0MiB |
|    0   N/A  N/A     21959      C   ...lib/libreoffice/program/soffice.bin         92MiB |
|    0   N/A  N/A     32453      C   ...s/miniforge3/envs/rapids/bin/python        218MiB |
+-----------------------------------------------------------------------------------------+

Il s’agit d’une RTX 2070. Comme le montrent la commande deviceQuery exécu­tée plus tôt ou la docu­men­ta­tion NVIDIA, elle dispose de 2304 cœurs et d’une fréquence maxi­male de 1,71Ghz.

Soit un total de 3 939.94 milliards d’ins­truc­tions à la seconde.

nb_coeurs = 2304
frequence = 1.71 #  en Giga Hertz

nombre_instructions_gpu = nb_coeurs * frequence
print("Nombre d'instructions à la seconde du GPU en milliards:", nombre_instructions_gpu)

Nombre d'instructions à la seconde du GPU en milliards : 3939.84

Nous aurions donc pu espé­rer un rapport de vitesse envi­ron 142 fois plus rapide en faveur du GPU. Nous en sommes loin.

nombre_instructions_gpu / nombre_instructions_cpu

142.74782608695654

Remarque

Nous avons comparé le nombre de cycles des proces­seurs pour compa­rer leur puis­sance respec­tive.
Dans cet exemple, les calculs sont réali­sés avec des réels 32 bits.

Il aurait été plus juste d’uti­li­ser les FLOPS (Floa­ting Point Opera­tion Per Second) comme mesure.
Les flops repré­sentent le nombre d’opé­ra­tions sur des réels que peuvent réali­ser des proces­seurs en 1 seconde.

Dans ce dernier cas, le proces­seur Intel possède des instruc­tions SIMD – Single Instruc­tion Multiple Data, aussi appe­lées vecto­rielles, de type AVX2 qui permettent d’exé­cu­ter des opéra­tions sur 8 réels 32 bits simul­ta­né­ment en 1 cycle.

bash$ lscpu | grep -i drapaux
Drapaux :                               fpu vme de pse tsc msr pae mce cx8 apic sep mtrr pge mca cmov pat pse36 clflush dts acpi mmx fxsr sse sse2 ss 
ht tm pbe syscall nx pdpe1gb rdtscp lm constant_tsc art arch_perfmon pebs bts rep_good nopl xtopology nonstop_tsc cpuid aperfmperf pni pclmulqdq dtes64
monitor ds_cpl vmx smx est tm2 ssse3 sdbg fma cx16 xtpr pdcm pcid sse4_1 sse4_2 x2apic movbe popcnt tsc_deadline_timer aes xsave avx f16c rdrand lahf_lm 
abm 3dnowprefetch cpuid_fault epb pti ssbd ibrs ibpb stibp tpr_shadow flexpriority ept vpid ept_ad fsgsbase tsc_adjust bmi1 avx2 smep bmi2 erms invpcid 
mpx rdseed adx smap clflushopt intel_pt xsaveopt xsavec xgetbv1 xsaves dtherm ida arat pln pts hwp hwp_notify hwp_act_window hwp_epp vnmi md_clear 
flush_l1d arch_capabilities

Il convien­drait donc de divi­ser ce ratio par 8. (soit 17.75 fois plus rapide pour le GPU)

Les cartes NVIDIA possèdent aussi des instruc­tions SIMD pour les modèles ayant les capa­ci­tés de calcul 3.0 jusqu’à 5.0. Cepen­dant, elles ont été reti­rées pour la plupart à partir de la capa­cité 5.0. Mais ceci est un autre sujet.

Nous revien­drons plus tard sur la notion de capa­cité de calcul de la carte qui corres­pond à une nomen­cla­ture interne à NVIDIA. Pour infor­ma­tion, la carte RTX 2070 corres­pond à la nomem­cla­ture 7.5 comme affi­ché par la commande deviceQuery précé­dente ; elle ne possède que peu d’ins­truc­tions SIMD (unique­ment un jeu réduit pour la vidéo et les nombres entiers).

142/8  # Ratio GPU/CPU en tenant compte des instructions SIMD du processeur

17.75

Ques­tion : Pourquoi notre test CuPy n’a-t-il pas montré un tel gain ?

Il est en effet plus de 3 fois infé­rieur au ratio que l’on pour­rait espé­rer ?

Nous allons tenter d’ap­por­ter des éléments de réponses avec un autre test réalisé via la librai­rie PyCuda.

Un autre exemple avec PyCuda

Conti­nuons donc, pour assou­vir votre soif de connais­sances et d’ex­pé­ri­men­ta­tion de votre carte graphique dernier cri, avec la librai­rie PyCuda.

PyCuda est un projet proche de PyOpenCL. Elle permet d’uti­li­ser l’API CUDA (langage C) des cartes NVIDIA depuis Python.
Elle ne fonc­tionne qu’avec les cartes NVIDIA.

PyCuda ne s’ins­talle pas très faci­le­ment si votre tool­kit NVIDIA n’est pas bien confi­guré au niveau du système.

Pensez à véri­fier que le compi­la­teur nvcc est bien dans le PATH système. A défaut, ajou­tez-le comme présenté dans le para­graphe sur la compi­la­tion.
N’ou­bliez pas de redé­mar­rer Jupy­ter pour qu’il ait connais­sance de ce chan­ge­ment.

Ensuite, l’ins­tal­la­tion avec pip sera assez facile :

$ pip install pycuda
Collecting pycuda
...
Successfully built pycuda
Installing collected packages: pycuda
Successfully installed pycuda-<version>

PyCUDA permet d’uti­li­ser la librai­rie CUDA direc­te­ment depuis Python, avec quelques avan­tages :

• Elle gère la mémoire pour vous, vous n’au­rez pas à vous soucier des allo­ca­tions et libé­ra­tions de mémoire
• Elle offre une classe GPUAr­ray qui faci­lite l’échange de données entre le CPU et le GPU
• Elle propose depuis Python un accès à l’in­té­gra­lité de la librai­rie CUDA, mais cela impo­sera souvent d’écrire des bouts de codes en C qui seront tran­mis à la librai­rie

Avant de commen­cer à lancer des calculs sur votre GPU, il convient de char­ger, puis d’ini­tia­li­ser le driver de la carte graphique :

import pycuda.driver as drv
drv.init()

Nous dispo­sons alors de quelques fonc­tions pour iden­ti­fier la confi­gu­ra­tion du maté­riel :

drv.Device.count()  # Nombre de GPU disponibles

1

gpu0 = drv.Device(0)  # Accès au premier GPU

gpu0.name()  # Nom de la carte

'NVIDIA GeForce RTX 2070'

gpu0.total_memory()  # mémoire disponible pour le GPU

8144093184

gpu0.pci_bus_id()  # BUS PCI utilisé

'0000:01:00.0'

La méthode get_attributes four­nira beau­coup d’autres infor­ma­tions :

attributes = gpu0.get_attributes()  # Dictionary containing gpu specifications
for attr, value in attributes.items():
    print(f"{attr} : {value}")

ASYNC_ENGINE_COUNT : 3
CAN_MAP_HOST_MEMORY : 1
CAN_USE_HOST_POINTER_FOR_REGISTERED_MEM : 1
CLOCK_RATE : 1440000
COMPUTE_CAPABILITY_MAJOR : 7
COMPUTE_CAPABILITY_MINOR : 5
COMPUTE_MODE : DEFAULT
COMPUTE_PREEMPTION_SUPPORTED : 1
CONCURRENT_KERNELS : 1
CONCURRENT_MANAGED_ACCESS : 1
DIRECT_MANAGED_MEM_ACCESS_FROM_HOST : 0
ECC_ENABLED : 0
GENERIC_COMPRESSION_SUPPORTED : 0
GLOBAL_L1_CACHE_SUPPORTED : 1
GLOBAL_MEMORY_BUS_WIDTH : 256
GPU_OVERLAP : 1
HANDLE_TYPE_POSIX_FILE_DESCRIPTOR_SUPPORTED : 1
HANDLE_TYPE_WIN32_HANDLE_SUPPORTED : 0
HANDLE_TYPE_WIN32_KMT_HANDLE_SUPPORTED : 0
HOST_NATIVE_ATOMIC_SUPPORTED : 0
INTEGRATED : 0
KERNEL_EXEC_TIMEOUT : 1
L2_CACHE_SIZE : 4194304
LOCAL_L1_CACHE_SUPPORTED : 1
MANAGED_MEMORY : 1
MAXIMUM_SURFACE1D_LAYERED_LAYERS : 2048
MAXIMUM_SURFACE1D_LAYERED_WIDTH : 32768
MAXIMUM_SURFACE1D_WIDTH : 32768
MAXIMUM_SURFACE2D_HEIGHT : 65536
MAXIMUM_SURFACE2D_LAYERED_HEIGHT : 32768
MAXIMUM_SURFACE2D_LAYERED_LAYERS : 2048
MAXIMUM_SURFACE2D_LAYERED_WIDTH : 32768
MAXIMUM_SURFACE2D_WIDTH : 131072
MAXIMUM_SURFACE3D_DEPTH : 16384
MAXIMUM_SURFACE3D_HEIGHT : 16384
MAXIMUM_SURFACE3D_WIDTH : 16384
MAXIMUM_SURFACECUBEMAP_LAYERED_LAYERS : 2046
MAXIMUM_SURFACECUBEMAP_LAYERED_WIDTH : 32768
MAXIMUM_SURFACECUBEMAP_WIDTH : 32768
MAXIMUM_TEXTURE1D_LAYERED_LAYERS : 2048
MAXIMUM_TEXTURE1D_LAYERED_WIDTH : 32768
MAXIMUM_TEXTURE1D_LINEAR_WIDTH : 268435456
MAXIMUM_TEXTURE1D_MIPMAPPED_WIDTH : 32768
MAXIMUM_TEXTURE1D_WIDTH : 131072
MAXIMUM_TEXTURE2D_ARRAY_HEIGHT : 32768
MAXIMUM_TEXTURE2D_ARRAY_NUMSLICES : 2048
MAXIMUM_TEXTURE2D_ARRAY_WIDTH : 32768
MAXIMUM_TEXTURE2D_GATHER_HEIGHT : 32768
MAXIMUM_TEXTURE2D_GATHER_WIDTH : 32768
MAXIMUM_TEXTURE2D_HEIGHT : 65536
MAXIMUM_TEXTURE2D_LINEAR_HEIGHT : 65000
MAXIMUM_TEXTURE2D_LINEAR_PITCH : 2097120
MAXIMUM_TEXTURE2D_LINEAR_WIDTH : 131072
MAXIMUM_TEXTURE2D_MIPMAPPED_HEIGHT : 32768
MAXIMUM_TEXTURE2D_MIPMAPPED_WIDTH : 32768
MAXIMUM_TEXTURE2D_WIDTH : 131072
MAXIMUM_TEXTURE3D_DEPTH : 16384
MAXIMUM_TEXTURE3D_DEPTH_ALTERNATE : 32768
MAXIMUM_TEXTURE3D_HEIGHT : 16384
MAXIMUM_TEXTURE3D_HEIGHT_ALTERNATE : 8192
MAXIMUM_TEXTURE3D_WIDTH : 16384
MAXIMUM_TEXTURE3D_WIDTH_ALTERNATE : 8192
MAXIMUM_TEXTURECUBEMAP_LAYERED_LAYERS : 2046
MAXIMUM_TEXTURECUBEMAP_LAYERED_WIDTH : 32768
MAXIMUM_TEXTURECUBEMAP_WIDTH : 32768
MAX_BLOCKS_PER_MULTIPROCESSOR : 16
MAX_BLOCK_DIM_X : 1024
MAX_BLOCK_DIM_Y : 1024
MAX_BLOCK_DIM_Z : 64
MAX_GRID_DIM_X : 2147483647
MAX_GRID_DIM_Y : 65535
MAX_GRID_DIM_Z : 65535
MAX_PERSISTING_L2_CACHE_SIZE : 0
MAX_PITCH : 2147483647
MAX_REGISTERS_PER_BLOCK : 65536
MAX_REGISTERS_PER_MULTIPROCESSOR : 65536
MAX_SHARED_MEMORY_PER_BLOCK : 49152
MAX_SHARED_MEMORY_PER_BLOCK_OPTIN : 65536
MAX_SHARED_MEMORY_PER_MULTIPROCESSOR : 65536
MAX_THREADS_PER_BLOCK : 1024
MAX_THREADS_PER_MULTIPROCESSOR : 1024
MEMORY_CLOCK_RATE : 7001000
MEMORY_POOLS_SUPPORTED : 1
MULTIPROCESSOR_COUNT : 36
MULTI_GPU_BOARD : 0
MULTI_GPU_BOARD_GROUP_ID : 0
PAGEABLE_MEMORY_ACCESS : 1
PAGEABLE_MEMORY_ACCESS_USES_HOST_PAGE_TABLES : 0
PCI_BUS_ID : 1
PCI_DEVICE_ID : 0
PCI_DOMAIN_ID : 0
READ_ONLY_HOST_REGISTER_SUPPORTED : 0
RESERVED_SHARED_MEMORY_PER_BLOCK : 0
SINGLE_TO_DOUBLE_PRECISION_PERF_RATIO : 32
STREAM_PRIORITIES_SUPPORTED : 1
SURFACE_ALIGNMENT : 512
TCC_DRIVER : 0
TEXTURE_ALIGNMENT : 512
TEXTURE_PITCH_ALIGNMENT : 32
TOTAL_CONSTANT_MEMORY : 65536
UNIFIED_ADDRESSING : 1
WARP_SIZE : 32

Nous avons main­te­nant un aperçu de comment nous pouvons inter­ro­ger le maté­riel dont nous dispo­sons.

Utili­sons cet exemple issu de la page d’ac­cueil du site de PyCuda afin d’éprou­ver une nouvelle fois la puis­sance de ce GPU.

Ici, le programme présenté utilise un petit bout de code C pour défi­nir la fonc­tion qui sera lancée sur le GPU. En effet, CUDA permet à la base d’exé­cu­ter des programmes écrits en langage C. Aussi, il est souvent néces­saire de devoir écrire quelques fonc­tions en C avec PyCuda.

Le bout de code ci-après est une petite astuce dans votre note­book pour vous assu­rer que Python trou­vera bien votre compi­la­teur nvcc requis pour la suite du TP.

import os 
print(os.environ['PATH'])
os.environ['PATH'] += ':/usr/local/cuda/bin'

/usr/local/cuda/bin:/home/aliquis/miniforge3/envs/JBook/bin:/home/aliquis/miniforge3/condabin:/home/aliquis/.local/bin:/home/aliquis/bin:/usr/local/sbin:/usr/local/bin:/usr/sbin:/usr/bin:/sbin:/bin:/usr/games:/usr/local/games:/snap/bin

Dans ce programme la fonc­tion multiply_them multi­plie 2 tableaux de réels, a et b pour stocker le résul­tat dans un tableau dest

import pycuda.autoinit  # Another way to initialize the GPU driver
import pycuda.driver as drv  # Here it is loaded
from pycuda.compiler import SourceModule  # Use to compile C functions for Cuda using nvcc
import numpy as np

# Mod will compile the C function
mod = SourceModule("""
__global__ void multiply_them(float *dest, float *a, float *b)
{
  const int i = threadIdx.x;
  dest[i] = a[i] * b[i];
}
""")

# Then the C function is retrieved this way so it can be called directly from Python 
multiply_them = mod.get_function("multiply_them")

Quelques expli­ca­tions :

• La fonc­tion C, qui est créée, est celle qui sera exécu­tée sur le GPU.
  On l’ap­pelle fonc­tion device.
  Elle prend 3 tableaux de réels 32 bits en entrée (poin­teurs float)
  Le premier contien­dra le résul­tat, les 2 autres sont les vecteurs que l’on souhaite multi­plier.
• Remarquez qu’au­cune boucle n’est réali­sée pour multi­plier tous les éléments des tableaux
  Au lieu de cela, la fonc­tion récu­père l’in­dice de l’élé­ment qu’elle doit calcu­ler via ThreadIdx.x
• La fonc­tion C qui aura été compi­lée dans la foulée par SourceModule est rendue acces­sible en Python via la méthode get_function

Comment la fonc­tion itère-t-elle sur tous les éléments du tableau ?

En effet, elle ne mani­pule qu’un élément à l’in­dice threadIdx.x. Il n’y a pas de boucle.

Le concept est le suivant: quand le GPU va lancer le calcul il va lancer cette fonc­tion sur> diffé­rents cœurs/threads en paral­lèle. Chaque cœur/thread possède un indice – ici threadIdx.x – et s’oc­cu­pera de la case située à cet indice.

Il n’y a donc pas besoin d’une boucle. Mais il faudra indiquer au GPU combien de fois la fonc­tion devra être lancée : autant qu’il y a d’élé­ments dans le tableau !

Réali­sons le même travail avec NumPy :

a = np.random.randn(400).astype(np.float32)  # We create 2 tables of numpy numbers having 400 items
b = np.random.randn(400).astype(np.float32)  # single precision is used as GPU are computing 
                                             # most of the time with float having 32 bits

dest = np.zeros_like(a)  # dest is a table full of zeros and having same size as a !

print(a[:5])
print(b[:5])
print(dest[:5])

[-0.00555495  0.25411418 -0.49110463  0.31549597 -0.17547004]
[0.42375204 0.90399384 0.1391739  0.44877696 1.2170142 ]
[0. 0. 0. 0. 0.]

Mesu­rons le temps d’exé­cu­tion de a * b avec NumPy sur CPU :

%timeit dest = a * b

422 ns ± 5.96 ns per loop (mean ± std. dev. of 7 runs, 1,000,000 loops each)

r = a * b
print(r[:5])

[-0.00235392  0.22971766 -0.06834894  0.14158732 -0.21354952]

Là aussi, il n’y a pas de boucle dans notre code. Mais elle existe bien dans la réalité : c’est l’opé­ra­teur * implé­menté dans la classe NDAr­ray de NumPy qui la réalise pour nous. Numpy ne travaille que sur 1 seul cœur. La librai­rie n’est pas paral­lé­li­sée sur tous les cœurs de la machine.

Dans cet exemple, le résul­tat devrait être envi­ron, en termes de nombre de cycles, (142 * 6) fois plus rapide en faveur de CuPy…

Excepté que nous n’uti­li­sons que 400 cœurs du GPU.

Le ration devrait donc être :

nb_coeurs = 400
frequence = 1.71 #  en Giga Hertz

nombre_instructions_gpu_400 = nb_coeurs * frequence

nb_coeurs = 1
frequence = 4.6 #  en Giga Hertz

nombre_instructions_cpu_1 = nb_coeurs * frequence

nombre_instructions_gpu_400 / nombre_instructions_cpu_1   # 148 fois plus rapide en théorie...

148.69565217391306

Véri­fions cela :

import cupy as cp
ac = cp.random.randn(400).astype(cp.float32)  
bc = cp.random.randn(400).astype(cp.float32)  

destc = cp.zeros_like(ac)  # dest is a table full of zeros and having same size as a !

print(ac[:5])
print(bc[:5])
print(destc[:5])

[-1.6642687   0.7745935   1.227624    0.5488435  -0.55261844]
[-1.9712749  -0.49186102 -0.06115467  0.69422895  0.68600804]
[0. 0. 0. 0. 0.]

%timeit destc = ac * bc

10.6 μs ± 46.9 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)

10.6 / 0.422  # 10us pour cupy contre 422ns soit 0.422us pour numpy...

25.118483412322274

Mais nous sommes 25 fois plus lent que le CPU ! Cela peut sembler n’avoir aucun sens !
Sur 1 seul cœur, le CPU est 25 fois plus rapide que le GPU qui utilise norma­le­ment 400 cœurs dans cet exemple !

Misère ! Rembour­sez-nous cette carte graphique, elle ne semble pas fonc­tion­ner correc­te­ment !

Véri­fions la mesure avec la fonc­tion multiply_them.

Nous allons pouvoir lancer le même calcul sur le GPU avec la fonc­tion précé­dem­ment créée.

• Il convient pour cela de trans­mettre les tableaux Numpy à la carte graphique, car, rappe­lez-vous, le GPU travaille dans sa propre mémoire. Ceci est réalisé via les fonc­tions drv.In ci-après.
• Pour récu­pé­rer le résul­tat, la fonc­tion drv.Out est utili­sée afin de trans­mettre les données de la RAM du GPU vers le CPU
• Les para­mètres block et grid défi­nissent la « matrice » des cœurs de calculs qui seront utili­sés.
  Ces notions seront détaillées dans la section suivante.
  Ici nous réser­vons une grille de 400 cœurs/threads.
• Nous avons un tableau de 400 éléments, les 400 multi­pli­ca­tions seront réali­sées en paral­lèle.
  Sous réserve que la carte dispose de suffi­sam­ment de cœurs (elle en a plus de 2000)

import numpy as np
a = np.random.randn(400).astype(np.float32)  # We create 2 tables of numpy numbers having 400 items
b = np.random.randn(400).astype(np.float32)  # single precision is used as GPU are computing 
                                             # most of the time with float having 32 bits

dest = np.zeros_like(a)
print(dest[:5])

multiply_them( drv.Out(dest)
              , drv.In(a)
              , drv.In(b)
              , block=(400, 1, 1)
              , grid=(1, 1) )

print(a[:5])
print(b[:5])
print(dest[:5])

[0. 0. 0. 0. 0.]
[-0.08468395  0.18728326  1.7134982  -1.051074   -0.35594532]
[-0.8216796  -0.07424407  0.8293001  -0.8038413  -1.1592786 ]
[ 0.06958307 -0.01390467  1.4210043   0.84489673  0.4126398 ]

print( sum(a*b - dest) )  # checking that result is good

0.0

Main­te­nant que nous avons véri­fié que le programme fonc­tion­nait correc­te­ment, nous pouvons calcu­ler son temps d’exé­cu­tion :

%timeit multiply_them( drv.Out(dest), drv.In(a), drv.In(b), block=(400, 1, 1), grid=(1, 1) )

164 μs ± 2.29 μs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)

164/0.422  # 388 times slower for the GPU :(

388.6255924170616

Misère !

Mais c’est encore pire ! Avec une fonc­tion écrite en langage C – et il est diffi­cile de faire plus concis, le GPU est cette fois plus de 388 fois plus lent que le CPU.

Mais que se passe-t-il avec cette carte ?

Afin d’ob­te­nir de véri­tables perfor­mances en faveur du GPU, il convient de comprendre un peu plus l’ar­chi­tec­ture d’un GPU.

Nous vous propo­sons de plon­ger plus en profon­deur dans l’ar­chi­tec­ture maté­rielle de votre carte avec la troi­sième partie de cette série d’ar­ticles.

En savoir plus

• Article Calcu­­lez sur GPU avec Python – Partie 1/3
• Article Calcu­­lez sur GPU avec Python – Partie 3/3
• Décou­vrez toutes nos forma­tions IA / Data Science
• Notre exper­tise IA / Data Science