Cvičenie čerpá z knihy Hands-On GPU Programming with Python and CUDA: Explore high-performance parallel computing with CUDA (Dr. Brian Tuomanen).
V minulom cvičení sme sa stretli s konkurenciou v rámci jedného kernelu, kedy bol tento kernel spúšťaný vo viacerých vláknach. Rozhranie CUDA tiež poskytuje ďalšiu úroveň konkurencie a to beh viacerých kernelov a pamäťových operácií naraz. Keďže niektoré z tých operácií môžu byť na sebe závislé (napríklad kopírovanie vstupnej matice do pamäte GPU a následné spustenie kernelu, ktorý s ňou pracuje), musia byť použité určité synchronizačné mechanizmy.
CUDA prúd (angl. stream) je sekvencia operácií, ktorá je spúšťaná na GPU nezávisle od iných prúdov.
CUDA udalosť (angl. event) je objekt, ktorý umožňuje oznámiť hostovi, že sa v kerneli stala nejaká udalosť (napr. dosiahla sa určitá časť kódu). Taktiež umožňuje presné časovanie kernelov.
Pri programovaní CUDA pomocou modulu numba
sme sa zatiaľ nestretli so synchronizáciou zariadenia, keďže bola robená automaticky za nás po každej operácii na grafickej karte. Pozrime sa preto na kód napísaný v CUDA C, aby sme videli, kde všade sa synchronizácia používa:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 |
// Skopiruj pole desatinnych cisel z hosta na GPU cudaMemcpy(device_array, host_array, size_of_array*sizeof(float), cudaMemcpyHostToDevice); // Zablokuj vykonavanie az kym nie je pamatova operacia hotova cudaDeviceSynchronize(); // Spusti CUDA kernel Some_CUDA_Kernel <<< block_size, grid_size >>> (device_array, size_of_array); // Zablokuj vykonavanie az kym GPU kernel nedokonci pracu cudaDeviceSynchronize(); // Skopiruj vystup z kernulu do hosta cudaMemcpy(host_array, device_array, size_of_array*sizeof(float), cudaMemcpyDeviceToHost); // Zablokuj vykonavanie az kym nie je hotovy pamatovy prenos cudaDeviceSynchronize(); |
V kóde je možné vidieť, že po každej GPU operácii robíme synchronizáciu. Ak by sme ale mali kernely a dáta, ktoré navzájom nesúvisia, mohli by sme ich spúšťať nezávisle na seba a nemuseli by sme čakať, kým predchádzajúca operácia skončí. V príkladoch, ktoré nasledujú sa na to pozrieme.
Prvý program vytvorí niekoľko polí s náhodným obsahom, potom každé pole spracuje jednoduchým kernelom a na záver skopíruje výsledok späť na hosta. Najprv si ukážeme verziu, v ktorej idú operácie za sebou v jednom prúde a následne program upravíme tak, aby jednotlivé kernely bežali vo viacerých prúdoch naraz. Zároveň odmeriame získané vylepšenie vo výkone.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 |
from numba import cuda import numpy as np from time import perf_counter NUM_ARRAYS = 200 ARRAY_LEN = 1024**2 @cuda.jit def kernel(array): thd = cuda.grid(1) num_iters = array.size / cuda.blockDim.x for j in range(num_iters): i = j * cuda.blockDim.x + thd for k in range(50): array[i] *= 2 array[i] /= 2 data = [] data_gpu = [] gpu_out = [] for _ in range(NUM_ARRAYS): data.append(np.random.randn(ARRAY_LEN).astype('float32')) t_start = perf_counter() for k in range(NUM_ARRAYS): data_gpu.append(cuda.to_device(data[k])) for k in range(NUM_ARRAYS): kernel[1, 64](data_gpu[k]) for k in range(NUM_ARRAYS): gpu_out.append(data_gpu[k].copy_to_host()) t_end = perf_counter() for k in range(NUM_ARRAYS): assert(np.allclose(gpu_out[k], data[k])) print(f'Total time: {t_end - t_start:.2f} s') |
1 |
Total time: 7.71 s |
Pre každú maticu vytvoríme jeden prúd. Kopírovanie dát a vykonávanie kernelov bude asynchrónne vzhľadom na iné prúdy.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 |
from numba import cuda import numpy as np from time import perf_counter NUM_ARRAYS = 200 ARRAY_LEN = 1024**2 @cuda.jit def kernel(array): thd = cuda.grid(1) num_iters = array.size / cuda.blockDim.x for j in range(num_iters): i = j * cuda.blockDim.x + thd for k in range(50): array[i] *= 2 array[i] /= 2 data = [] data_gpu = [] gpu_out = [] streams = [] for _ in range(NUM_ARRAYS): streams.append(cuda.stream()) data.append(np.random.randn(ARRAY_LEN).astype('float32')) t_start = perf_counter() for k in range(NUM_ARRAYS): data_gpu.append(cuda.to_device(data[k], stream=streams[k])) for k in range(NUM_ARRAYS): kernel[1, 64, streams[k]](data_gpu[k]) for k in range(NUM_ARRAYS): gpu_out.append(data_gpu[k].copy_to_host(stream=streams[k])) t_end = perf_counter() for k in range(NUM_ARRAYS): assert(np.allclose(gpu_out[k], data[k])) print(f'Total time: {t_end - t_start:.2f} s') |
1 |
Total time: 2.72 s |
Použitím prúdov sme dosiahli skoro trojnásobné zlepšenie výkonu. Je to vďaka tomu, že jednotlivé operácie (presun polí na GPU, spustenie kernelu a presun dát späť do hosta) bežali v jednotlivých prúdoch asynchrónne.
Teraz ukážeme použitie prúdov na zložitejšom príklade. Bude ním hra života od amerického matematika Johna Conwaya.
Nasleduje implementácia pre jeden bežiaci kernel.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 |
# Conway's game of life in Python / CUDA C # written by Brian Tuomanen for "Hands on GPU Programming with Python and CUDA" from numba import cuda import numpy as np import matplotlib.pyplot as plt import matplotlib.animation as animation @cuda.jit(device=True) def nbrs(x, y, matrix): return (matrix[x - 1, y + 1] + matrix[x - 1, y] + matrix[x - 1, y - 1] + matrix[x, y + 1] + matrix[x, y - 1] + matrix[x + 1, y - 1] + matrix[x + 1, y] + matrix[x + 1, y + 1]) @cuda.jit def kernel(lattice_out, lattice): x, y = cuda.grid(2) n = nbrs(x, y, lattice) if lattice[x, y] == 1: if n in (2, 3): lattice_out[x, y] = 1 else: lattice_out[x, y] = 0 elif lattice[x, y] == 0: if n == 3: lattice_out[x, y] = 1 else: lattice_out[x, y] = 0 def update_gpu(frameNum, img, new_lattice_gpu, lattice_gpu, N): blockdim = (N // 32, N // 32) griddim = (32, 32) kernel[griddim, blockdim](new_lattice_gpu, lattice_gpu) img.set_data(new_lattice_gpu.copy_to_host()) lattice_gpu[:] = new_lattice_gpu[:] return img if __name__ == '__main__': # set lattice size N = 128 lattice = np.int32( np.random.choice([1, 0], N * N, p=[0.25, 0.75]).reshape(N, N)) lattice_gpu = cuda.to_device(lattice) new_lattice_gpu = cuda.device_array_like(lattice_gpu) fig, ax = plt.subplots() img = ax.imshow(lattice_gpu.copy_to_host(), interpolation='nearest') ani = animation.FuncAnimation(fig, update_gpu, fargs=( img, new_lattice_gpu, lattice_gpu, N, ), interval=33, frames=1000) plt.show() |
Pokračujeme verziou, v ktorej budú konkurentne prebiehať štyri hry. Na to vytvoríme štyri prúdy a taktiež vytvoríme štyri dvojice matíc.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 |
# Conway's game of life in Python / CUDA C # written by Brian Tuomanen for "Hands on GPU Programming with Python and CUDA" from numba import cuda import numpy as np import matplotlib.pyplot as plt import matplotlib.animation as animation @cuda.jit(device=True) def nbrs(x, y, matrix): return (matrix[x - 1, y + 1] + matrix[x - 1, y] + matrix[x - 1, y - 1] + matrix[x, y + 1] + matrix[x, y - 1] + matrix[x + 1, y - 1] + matrix[x + 1, y] + matrix[x + 1, y + 1]) @cuda.jit def kernel(lattice_out, lattice): x, y = cuda.grid(2) n = nbrs(x, y, lattice) if lattice[x, y] == 1: if n in (2, 3): lattice_out[x, y] = 1 else: lattice_out[x, y] = 0 elif lattice[x, y] == 0: if n == 3: lattice_out[x, y] = 1 else: lattice_out[x, y] = 0 def update_gpu(frameNum, imgs, new_lattices_gpu, lattices_gpu, N, streams, num_concurrent): blockdim = (N // 32, N // 32) griddim = (32, 32) for k in range(num_concurrent): kernel[griddim, blockdim, streams[k]](new_lattices_gpu[k], lattices_gpu[k]) imgs[k].set_data(new_lattices_gpu[k].copy_to_host(stream=streams[k])) lattices_gpu[k].copy_to_device(new_lattices_gpu[k], stream=streams[k]) return imgs if __name__ == '__main__': # set lattice size N = 128 num_concurrent = 4 streams = [] lattices_gpu = [] new_lattices_gpu = [] for k in range(num_concurrent): streams.append(cuda.stream()) lattice = np.int32( np.random.choice([1, 0], N * N, p=[0.25, 0.75]).reshape(N, N)) lattices_gpu.append(cuda.to_device(lattice)) new_lattices_gpu.append(cuda.device_array_like(lattices_gpu[k])) fig, ax = plt.subplots(nrows=1, ncols=num_concurrent) imgs = [] for k in range(num_concurrent): imgs.append(ax[k].imshow( lattices_gpu[k].copy_to_host(stream=streams[k]), interpolation='nearest')) ani = animation.FuncAnimation(fig, update_gpu, fargs=(imgs, new_lattices_gpu, lattices_gpu, N, streams, num_concurrent), interval=33, frames=1000) plt.show() |
V poslednej časti cvičenia sa pozrieme na využitie udalostí. Udalosti sú objekty, ktorými môžeme zistiť, v akom stave je spracovanie prúdu. Zavolaním metódy record()
vložíme do aktuálneho alebo zvoleného prúdu míľnik. Táto metóda neblokuje, ale po spracovaní všetkých operácií v prúde bude udalosť nastavená ako hotová. Ak v programe nastavíme dve udalosti, môžeme zistiť časový rozdiel medzi nimi.
V tomto príklade spustíme náš jednoduchý kernel z predchádzajúcich príkladov na jednej matici s veľkým počtom prvkov (skoro 82 miliónov). Vytvoríme dve udalosti, prvú umiestnime pred začiatkom výpočtu kernelu a druhú po jeho skončení.
Pozrieme sa na celkový čas výpočtu kernelu, k čomu bude potrebná synchronizácia druhej udalosti.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 |
from numba import cuda import numpy as np array_len = 78 * 1024**2 @cuda.jit def kernel(array): thd = cuda.grid(1) num_iters = array.size / cuda.blockDim.x for j in range(num_iters): i = j * cuda.blockDim.x + thd for k in range(50): array[i] *= 2.0 array[i] /= 2.0 data = np.random.randn(array_len).astype('float32') data_gpu = cuda.to_device(data) start_event = cuda.event() end_event = cuda.event() start_event.record() kernel[1, 64](data_gpu) end_event.record() # pocka, kym nebude `end_event` oznaceny za hotovy # end_event.synchronize() # vypocita kolko trvalo spustenie kernelu. print('Kernel execution time in milliseconds: %f ' % cuda.event.elapsed_time(start_event, end_event)) |
Ak necháme zakomentovaný riadok v príklade, môžeme vidieť asynchrónnosť udalostí. Prvá udalosť je označená za hotovú, ale druhá ešte nie, keďže sa jednalo len o umiestnenie udalosti do prúdu a grafická karta nestihla potvrdiť dokončenie udalosti.
1 |
Kernel execution time in milliseconds: 88.773141 |
Po odkomentovaní synchronizačného riadku, v ktorom počkáme na ukončenie udalosti end_event
môžeme vypočítať celkový čas uplynutý medzi ukončením oboch udalostí – v tomto prípade sa jedná o čas behu kernelu.
1 |
Kernel execution time in milliseconds: 2180.442871 |
V tejto ukážke použijeme udalosti spoločne s prúdmi. Jedna udalosť patrí jednému prúdu, preto musíme vytvoriť toľko udalostí, koľko máme prúdov. Pre 200 prúdov budeme mať 400 udalostí, keďže chceme časovať dĺžku výpočtu kernelu.
Prvú udalosť daného prúdu umiestnime hneď pred začiatok kernelu. Zvyšok udalostí umiestnime v samostatnom cykle po tom, čo sú všetky kernely spustené. Je to práve z toho dôvodu, aby sme nezdržiavali spustenia kernelov.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 |
from numba import cuda import numpy as np from time import perf_counter num_arrays = 200 array_len = 1024**2 @cuda.jit def kernel(array): thd = cuda.grid(1) num_iters = array.size / cuda.blockDim.x for j in range(num_iters): i = j * cuda.blockDim.x + thd for k in range(50): array[i] *= 2 array[i] /= 2 data = [] data_gpu = [] gpu_out = [] streams = [] start_events = [] end_events = [] # generate random arrays. for _ in range(num_arrays): streams.append(cuda.stream()) start_events.append(cuda.event()) end_events.append(cuda.event()) data.append(np.random.randn(array_len).astype('float32')) t_start = perf_counter() # copy arrays to GPU. for k in range(num_arrays): data_gpu.append(cuda.to_device(data[k], stream=streams[k])) # process arrays. for k in range(num_arrays): start_events[k].record(streams[k]) kernel[1, 32, streams[k]](data_gpu[k]) for k in range(num_arrays): end_events[k].record(streams[k]) # copy arrays from GPU. for k in range(num_arrays): gpu_out.append(data_gpu[k].copy_to_host(stream=streams[k])) t_end = perf_counter() for k in range(num_arrays): assert (np.allclose(gpu_out[k], data[k])) kernel_times = [] for k in range(num_arrays): kernel_times.append(cuda.event_elapsed_time(start_events[k], end_events[k])) print('Total time: %f' % (t_end - t_start)) print('Mean kernel duration (milliseconds): %f' % np.mean(kernel_times)) print('Mean kernel standard deviation (milliseconds): %f' % np.std(kernel_times)) |
1 2 3 |
Total time: 4.563637 Mean kernel duration (milliseconds): 2173.889211 Mean kernel standard deviation (milliseconds): 181.448503 |
Vo výpise vidíme, že zatiaľ čo celkový čas výpočtu bol 4.6 sekundy, beh jedného kernelu trval okolo 2 sekúnd. Štandardná odchýlka s hodnotou 0.2 sekundy je oproti dĺžke behu kernelu malá, čo hovorí o dobrom využití zdrojov GPU. Ak by bola štandardná výchylka väčšia, mali by sme upraviť parametre výpočtu pre lepšie využitie grafickej karty.
Program z minulého týždňa upravte tak, aby 1) ste optimálne využili zariadenie a 2) v rámci riešenia použite prúdy. Na meranie času výpočty použite udalosti. Vypracujte sprievodnú dokumentáciu: akú úlohu riešite, ako ste hľadali optimálne parametre výpočtu (použitie profilerov, kalkulátor obsadenosti GPU, dynamické zisťovanie obsadenosti zariadenia pomocou API), aké zrýchlenie ste dosiahli a čím bolo zaručené.