Получите операции с помощью шаблона на графическом процессоре

Этот пример использует "Игру Конуэя Жизни", чтобы продемонстрировать, как операции шаблона могут быть выполнены с помощью графического процессора.

Много операций над массивами могут быть выражены как "операция шаблона", где каждый элемент выходного массива зависит от небольшой области входного массива. Примеры включают конечные разности, свертку, среднюю фильтрацию и методы конечных элементов. Этот пример использует "Игру Конуэя Жизни", чтобы продемонстрировать два способа запустить операцию шаблона на графическом процессоре, начинающем с кода в электронной книге Клева Молера Эксперименты в MATLAB.

"Игра Жизни" следует нескольким простым правилам:

Ячейки располагаются в 2D сетке
На каждом шаге судьба каждой ячейки определяется живучестью ее восьми самых близких соседей
Любая ячейка точно с тремя живыми соседями оживает на следующем шаге
Живая ячейка точно с двумя живыми соседями остается живой на следующем шаге
Все другие ячейки (включая тех больше чем с тремя соседями) умирают на следующем шаге или остаются пустыми

"Шаблон" в этом случае поэтому 3x3 область вокруг каждого элемента. Вот некоторые примеры того, как обновляется ячейка:

Этим примером является функция, чтобы позволить использование вложенных функций:

function paralleldemo_gpu_stencil()

Сгенерируйте случайную начальную генеральную совокупность

Начальная генеральная совокупность ячеек создается на 2D сетке примерно с 25% живых местоположений.

gridSize = 500;
numGenerations = 100;
initialGrid = (rand(gridSize,gridSize) > .75);
gpu = gpuDevice();

% Draw the initial grid
hold off
imagesc(initialGrid);
colormap([1 1 1;0 0.5 0]);
title('Initial Grid');

Играя в игру жизни

Электронная книга Эксперименты в MATLAB обеспечивает начальное внедрение, которое может использоваться для сравнения. Эта версия полностью векторизована, обновив все ячейки в сетке в одной передаче на генерацию.

    function X = updateGrid(X, N)
        p = [1 1:N-1];
        q = [2:N N];
        % Count how many of the eight neighbors are alive.
        neighbors = X(:,p) + X(:,q) + X(p,:) + X(q,:) + ...
            X(p,p) + X(q,q) + X(p,q) + X(q,p);
        % A live cell with two live neighbors, or any cell with
        % three live neighbors, is alive at the next step.
        X = (X & (neighbors == 2)) | (neighbors == 3);
    end

grid = initialGrid;
% Loop through each generation updating the grid and displaying it
for generation = 1:numGenerations
    grid = updateGrid(grid, gridSize);

    imagesc(grid);
    title(num2str(generation));
    drawnow;
end

Теперь повторно выполните игру и измерьтесь, сколько времени она берет для каждой генерации.

grid = initialGrid;
timer = tic();

for generation = 1:numGenerations
    grid = updateGrid(grid, gridSize);
end

cpuTime = toc(timer);
fprintf('Average time on the CPU: %2.3fms per generation.\n', ...
    1000*cpuTime/numGenerations);

Average time on the CPU: 11.323ms per generation.

Сохраните этот результат проверить правильность каждой версии ниже.

expectedResult = grid;

Преобразование Игры Жизни, чтобы работать на графическом процессоре

Чтобы запустить Игру Жизни на графическом процессоре, начальная генеральная совокупность отправляется в графический процессор с помощью gpuArray. Алгоритм остается неизменным. Обратите внимание на то, что wait(gpu) используется, чтобы гарантировать, что графический процессор закончил вычислять, прежде чем таймер останавливается. Это требуется только для точной синхронизации.

grid = gpuArray(initialGrid);
timer = tic();

for generation = 1:numGenerations
    grid = updateGrid(grid, gridSize);
end

wait(gpu); % Only needed to ensure accurate timing
gpuSimpleTime = toc(timer);

% Print out the average computation time and check the result is unchanged.
fprintf(['Average time on the GPU: %2.3fms per generation ', ...
    '(%1.1fx faster).\n'], ...
    1000*gpuSimpleTime/numGenerations, cpuTime/gpuSimpleTime);
assert(isequal(grid, expectedResult));

Average time on the GPU: 1.655ms per generation (6.8x faster).

Создание поэлементной версии для графического процессора

Смотря на вычисления в функции updateGrid, очевидно, что те же операции применяются в каждом местоположении сетки независимо. Это предполагает, что arrayfun мог использоваться, чтобы сделать оценку. Однако каждая ячейка должна знать о своих восьми соседях, повреждая поэлементную независимость. Каждый элемент должен смочь получить доступ к полной сетке, также работая независимо.

Решение состоит в том, чтобы использовать вложенную функцию. Вложенные функции, даже используемые с arrayfun, могут переменные доступа, объявленные в их родительской функции. Это означает, что каждая ячейка может считать целую сетку с предыдущего такта и индекса в него.

grid = gpuArray(initialGrid);

    function X = updateParentGrid(row, col, N)
        % Take account of boundary effects
        rowU = max(1,row-1);  rowD = min(N,row+1);
        colL = max(1,col-1);  colR = min(N,col+1);
        % Count neighbors
        neighbors ...
            = grid(rowU,colL) + grid(row,colL) + grid(rowD,colL) ...
            + grid(rowU,col)                   + grid(rowD,col) ...
            + grid(rowU,colR) + grid(row,colR) + grid(rowD,colR);
        % A live cell with two live neighbors, or any cell with
        % three live neighbors, is alive at the next step.
        X = (grid(row,col) & (neighbors == 2)) | (neighbors == 3);
    end


timer = tic();

rows = gpuArray.colon(1, gridSize)';
cols = gpuArray.colon(1, gridSize);
for generation = 1:numGenerations
    grid = arrayfun(@updateParentGrid, rows, cols, gridSize);
end

wait(gpu); % Only needed to ensure accurate timing
gpuArrayfunTime = toc(timer);

% Print out the average computation time and check the result is unchanged.
fprintf(['Average time using GPU arrayfun: %2.3fms per generation ', ...
    '(%1.1fx faster).\n'], ...
    1000*gpuArrayfunTime/numGenerations, cpuTime/gpuArrayfunTime);
assert(isequal(grid, expectedResult));

Average time using GPU arrayfun: 0.795ms per generation (14.2x faster).

Обратите внимание на то, что мы также использовали другую новую функцию arrayfun здесь: расширение размерности. Мы должны были передать только векторы строки и столбца, и они были автоматически расширены в полную сетку. Эффект состоит в том, как будто мы вызвали:

[cols,rows] = meshgrid(cols,rows);

как часть вызова arrayfun. Это сохраняет нас обоих некоторое вычисление и некоторая передача данных между памятью памяти ЦП и графического процессора.

Заключение

В этом примере простая операция шаблона, "Игра Конуэя Жизни", была реализована на графическом процессоре с помощью arrayfun и переменных, объявленных в родительской функции. Этот метод может использоваться, чтобы реализовать область значений операций шаблона включая алгоритмы конечного элемента, свертки и фильтры. Это может также использоваться к элементам доступа в интерполяционной таблице, заданной в родительской функции.

fprintf('CPU:          %2.3fms per generation.\n', ...
    1000*cpuTime/numGenerations);
fprintf('Simple GPU:   %2.3fms per generation (%1.1fx faster).\n', ...
    1000*gpuSimpleTime/numGenerations, cpuTime/gpuSimpleTime);
fprintf('Arrayfun GPU: %2.3fms per generation (%1.1fx faster).\n', ...
    1000*gpuArrayfunTime/numGenerations, cpuTime/gpuArrayfunTime);

CPU:          11.323ms per generation.
Simple GPU:   1.655ms per generation (6.8x faster).
Arrayfun GPU: 0.795ms per generation (14.2x faster).

end

Документация