Додавання двух масивів

Нехай перед нами стоїть задача по-елементного додавання двух масивів. Для початку, створимо два масиви розміром 1024, типу Float32 заповнені випадковими даними.

using CUDA

a = CUDA.rand(Float32, 1024^2)
b = CUDA.rand(Float32, 1024^2)

GPU-масиви підтримують більшість операцій що й звичайні масиви, тому для обчислення по-елементного додавання достатньо використати оператор .+.

Замітка: крапка перед оператором .+ дозволяє застосувати його до кожного елемента масивів та називається broadcast.

y = a .+ b

Тепер, перевіримо чи результат правильний порівнявши з обчисленням на процесорі. Для цього, перенесемо дані з GPU на процесор та обчислимо результат:

ah = Array(a)
bh = Array(b)
yh = ah .+ bh

Для порівняння, перенесемо результат з GPU на процесор та порівняємо два масиви:

res = isapprox(Array(y), yh)
@show res
@assert res

res = true

Julia також дозволяє використовувати Unicode символи, тому функуцію isapprox можна замінити на оператор ≈:

res = Array(y) ≈ yh
@show res
@assert res

res = true

---

Швидкість роботи на GPU та CPU

Для порівняння швидкості роботи використаємо пакет BenchmarkTools.jl, який виконує код декілька разів для збору статистики та більш точного заміру часу.

Виміряємо швидкість роботи на CPU:

using BenchmarkTools

@btime $ah .+ $bh

ArgumentError: Package BenchmarkTools not found in current path.
- Run `import Pkg; Pkg.add("BenchmarkTools")` to install the BenchmarkTools package.

Тепер виміряємо швидкість роботи GPU.

Замітка: весь код, що виконується на відеокарті, виконується асинхронно, тому для правильного вимірювання часу, нам наобхідно синхронізувати відеокарту після кожної операції за допомогою макро CUDA.@sync.

@btime CUDA.@sync $a .+ $b

LoadError: UndefVarError: `@btime` not defined
in expression starting at none:1

Бачимо, що відеокарта виконує обчислення в декілька разів швидше ніж CPU.

Також, для повноти картини, покажемо час у випадку якщо не синхронізовувати відеокарту. В такому випадку час буде не правильний, а ми будемо вимірювати швидкість створення команд для відеокарти, а не час їх виконання.

@btime $a .+ $b

LoadError: UndefVarError: `@btime` not defined
in expression starting at none:1

---

Перше ядро для відеокарти

Хоча для обчислення по-елементної суми можна використати оператор .+, важливо розуміти, що відбувається "під капотом" цієї операції.

Для багатьох операцій, Julia автоматично компілює функції (ядра), які обчислюють певні операції, як у випадку по-елементної суми.

Замітка: ядро - це функція, яка виконується на відеокарті.

Напишемо таке ядро вручну. Дякуючи можливостям Julia, нам не потрібно використовувати для цього С/C++ як це роблять у випадку з Python, адже ми можемо писати такі ядра одразу на Julia.

function vadd!(y, a, b)
    i = threadIdx().x + (blockIdx().x - 1) * blockDim().x
    if i ≤ length(a)
        y[i] = a[i] + b[i]
    end
    return
end

Як бачимо, ядро в Julia це звичайна функція. Ядра виконуються на решітці, розмір якої задає користувач. Решітка може бути 1D, 2D або 3D, та поділена на блоки, які також задаються користувачем.

В нашому випадку, решітка одновимірна (1D), розмір її дорівнює довжині нашого масиву даних (1024^2), а розмір блока задамо 256.

Розмір блока визначає кількість потоків які будуть виконуватись одним потоковим процесором (stream processor), але про це пізніше.

Для обчислення положення поточного потоку i на решітці, використаємо:

• його індекс в блоці (threadIdx().x) який в нашому випадку буде в межах від 1 до 256.

• поточний індекс блоку (blockIdx().x) який в нашому випадку буде від 1 до ceil(length(a) / 256).

• та власне розмір блоку (blockDim().x = 256).

Таким чином з кожним елементом масива асоційований унікальний потік з індексом i, який приймає значення від 1 до length(a), що ми й використали для обчислення суми.

Створимо не-ініціалізований масив, в який ми будемо записувати результат обчислення:

y = CuArray{Float32}(undef, length(a))

Для запуску ядра на відеокарті, необхідно викликати функцію vadd!, попередньо додавши макро @cuda та передавши агрументи в ядро:

threads = 256
blocks = cld(length(y), threads)
@cuda threads=threads blocks=blocks vadd!(y, a, b)

Після запуску необхідно синхронізувати відеокарту перші ніж зчитувати дані. Це можна зробити вручну викликом функції CUDA.synchronize(), проте при копіюванні даних з відеокарти на процесор, синхронізація відбувається автоматично. Тепер перевіримо результат, порівнявши з обчисленнями зробленими на процесорі:

res = Array(y) ≈ yh
@show res
@assert res

res = true

Таким чином ми успішно написали наше перше ядро для по-елементного додавання двух масивів!