Exponenta Banner
exponenta event banner

Сравнительный анализ A\b на GPU

Открыть сценарий

В этом примере показано, как мы можем сопоставить решение линейной системы на GPU. Код MATLAB ® для решения x в A*x = b очень просто. Чаще всего мы используем матричное левое деление, также известное как mldivide или оператор обратной косой черты (\) для вычисления x (то есть x = A\b).

Связанные примеры:

Код, показанный в этом примере, можно найти в следующей функции:

function results = paralleldemo_gpu_backslash(maxMemory)

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

if nargin == 0
    g = gpuDevice;
    maxMemory = 0.4*g.AvailableMemory/1024^3;
end

Функция бенчмаркинга

Мы хотим сопоставить матрицу с левым разделением (\), а не стоимость передачи данных между CPU и GPU, время, необходимое для создания матрицы, или другие параметры. Поэтому мы отделяем генерацию данных от решения линейной системы и измеряем только время, необходимое для выполнения последней.

function [A, b] = getData(n, clz)
    fprintf('Creating a matrix of size %d-by-%d.\n', n, n);
    A = rand(n, n, clz) + 100*eye(n, n, clz);
    b = rand(n, 1, clz);
end

function time = timeSolve(A, b, waitFcn)
    tic;
    x = A\b; %#ok<NASGU> We don't need the value of x.
    waitFcn(); % Wait for operation to complete.
    time = toc;
end

Выбор размера проблемы

Как и в случае с большим количеством других параллельных алгоритмов, производительность решения линейной системы параллельно в значительной степени зависит от размера матрицы. Как видно из других примеров, таких как сравнительный анализ A\b, мы сравниваем производительность алгоритма для различных размеров матрицы.

% Declare the matrix sizes to be a multiple of 1024.
maxSizeSingle = floor(sqrt(maxMemory*1024^3/4));
maxSizeDouble = floor(sqrt(maxMemory*1024^3/8));
step = 1024;
if maxSizeDouble/step >= 10
    step = step*floor(maxSizeDouble/(5*step));
end
sizeSingle = 1024:step:maxSizeSingle;
sizeDouble = 1024:step:maxSizeDouble;

Сравнение производительности: гигафлопс

Мы используем количество операций с плавающей запятой в секунду в качестве показателя производительности, потому что это позволяет сравнивать производительность алгоритма для различных размеров матрицы.

Учитывая размер матрицы, функция бенчмаркинга создает матрицу A и правая сторона b один раз, а затем решает A\b несколько раз, чтобы получить точную меру времени, которое требуется. Мы используем число операций с плавающей запятой HPC Challenge, так что для матрицы n-by-n мы считаем операции с плавающей запятой как 2/3*n^3 + 3/2*n^2.

Функция передается в дескрипторе функции «wait». На CPU эта функция ничего не делает. На графическом процессоре эта функция ожидает завершения всех ожидающих операций. Таким образом, ожидание обеспечивает точную синхронизацию.

function gflops = benchFcn(A, b, waitFcn)
    numReps = 3;
    time = inf;
    % We solve the linear system a few times and calculate the Gigaflops
    % based on the best time.
    for itr = 1:numReps
        tcurr = timeSolve(A, b, waitFcn);
        time = min(tcurr, time);
    end

    % Measure the overhead introduced by calling the wait function.
    tover = inf;
    for itr = 1:numReps
        tic;
        waitFcn();
        tcurr = toc;
        tover = min(tcurr, tover);
    end
    % Remove the overhead from the measured time. Don't allow the time to
    % become negative.
    time = max(time - tover, 0);
    n = size(A, 1);
    flop = 2/3*n^3 + 3/2*n^2;
    gflops = flop/time/1e9;
end

% The CPU doesn't need to wait: this function handle is a placeholder.
function waitForCpu()
end

% On the GPU, to ensure accurate timing, we need to wait for the device
% to finish all pending operations.
function waitForGpu(theDevice)
    wait(theDevice);
end

Выполнение тестов

Выполнив все настройки, легко выполнить тесты. Однако выполнение вычислений может занять много времени, поэтому мы печатаем некоторую промежуточную информацию о состоянии по мере выполнения сравнительного анализа для каждого размера матрицы. Мы также инкапсулируем цикл по всем размерам матрицы в функцию, чтобы проверить как вычисления с одной, так и двойной точностью.

function [gflopsCPU, gflopsGPU] = executeBenchmarks(clz, sizes)
    fprintf(['Starting benchmarks with %d different %s-precision ' ...
         'matrices of sizes\nranging from %d-by-%d to %d-by-%d.\n'], ...
            length(sizes), clz, sizes(1), sizes(1), sizes(end), ...
            sizes(end));
    gflopsGPU = zeros(size(sizes));
    gflopsCPU = zeros(size(sizes));
    gd = gpuDevice;
    for i = 1:length(sizes)
        n = sizes(i);
        [A, b] = getData(n, clz);
        gflopsCPU(i) = benchFcn(A, b, @waitForCpu);
        fprintf('Gigaflops on CPU: %f\n', gflopsCPU(i));
        A = gpuArray(A);
        b = gpuArray(b);
        gflopsGPU(i) = benchFcn(A, b, @() waitForGpu(gd));
        fprintf('Gigaflops on GPU: %f\n', gflopsGPU(i));
    end
end

Затем мы выполняем тесты с одинарной и двойной точностью.

[cpu, gpu] = executeBenchmarks('single', sizeSingle);
results.sizeSingle = sizeSingle;
results.gflopsSingleCPU = cpu;
results.gflopsSingleGPU = gpu;
[cpu, gpu] = executeBenchmarks('double', sizeDouble);
results.sizeDouble = sizeDouble;
results.gflopsDoubleCPU = cpu;
results.gflopsDoubleGPU = gpu;

Starting benchmarks with 7 different single-precision matrices of sizes
ranging from 1024-by-1024 to 19456-by-19456.
Creating a matrix of size 1024-by-1024.
Gigaflops on CPU: 43.805496
Gigaflops on GPU: 78.474002
Creating a matrix of size 4096-by-4096.
Gigaflops on CPU: 96.459635
Gigaflops on GPU: 573.278854
Creating a matrix of size 7168-by-7168.
Gigaflops on CPU: 184.997657
Gigaflops on GPU: 862.755636
Creating a matrix of size 10240-by-10240.
Gigaflops on CPU: 204.404384
Gigaflops on GPU: 978.362901
Creating a matrix of size 13312-by-13312.
Gigaflops on CPU: 218.773070
Gigaflops on GPU: 1107.983667
Creating a matrix of size 16384-by-16384.
Gigaflops on CPU: 233.529176
Gigaflops on GPU: 1186.423754
Creating a matrix of size 19456-by-19456.
Gigaflops on CPU: 241.482550
Gigaflops on GPU: 1199.151846
Starting benchmarks with 5 different double-precision matrices of sizes
ranging from 1024-by-1024 to 13312-by-13312.
Creating a matrix of size 1024-by-1024.
Gigaflops on CPU: 34.902918
Gigaflops on GPU: 72.191488
Creating a matrix of size 4096-by-4096.
Gigaflops on CPU: 74.458136
Gigaflops on GPU: 365.339897
Creating a matrix of size 7168-by-7168.
Gigaflops on CPU: 93.313782
Gigaflops on GPU: 522.514165
Creating a matrix of size 10240-by-10240.
Gigaflops on CPU: 104.219804
Gigaflops on GPU: 628.301313
Creating a matrix of size 13312-by-13312.
Gigaflops on CPU: 108.826886
Gigaflops on GPU: 681.881032

Печать производительности

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

Во-первых, мы рассматриваем производительность оператора обратной косой черты в одной точности.

fig = figure;
ax = axes('parent', fig);
plot(ax, results.sizeSingle, results.gflopsSingleGPU, '-x', ...
     results.sizeSingle, results.gflopsSingleCPU, '-o')
grid on;
legend('GPU', 'CPU', 'Location', 'NorthWest');
title(ax, 'Single-precision performance')
ylabel(ax, 'Gigaflops');
xlabel(ax, 'Matrix size');
drawnow;

Теперь рассмотрим производительность оператора обратной косой черты с двойной точностью.

fig = figure;
ax = axes('parent', fig);
plot(ax, results.sizeDouble, results.gflopsDoubleGPU, '-x', ...
     results.sizeDouble, results.gflopsDoubleCPU, '-o')
legend('GPU', 'CPU', 'Location', 'NorthWest');
grid on;
title(ax, 'Double-precision performance')
ylabel(ax, 'Gigaflops');
xlabel(ax, 'Matrix size');
drawnow;

Наконец, мы рассмотрим ускорение оператора обратной косой черты при сравнении графического процессора с CPU.

speedupDouble = results.gflopsDoubleGPU./results.gflopsDoubleCPU;
speedupSingle = results.gflopsSingleGPU./results.gflopsSingleCPU;
fig = figure;
ax = axes('parent', fig);
plot(ax, results.sizeSingle, speedupSingle, '-v', ...
     results.sizeDouble, speedupDouble, '-*')
grid on;
legend('Single-precision', 'Double-precision', 'Location', 'SouthEast');
title(ax, 'Speedup of computations on GPU compared to CPU');
ylabel(ax, 'Speedup');
xlabel(ax, 'Matrix size');
drawnow;


end

ans = 

         sizeSingle: [1024 4096 7168 10240 13312 16384 19456]
    gflopsSingleCPU: [1x7 double]
    gflopsSingleGPU: [1x7 double]
         sizeDouble: [1024 4096 7168 10240 13312]
    gflopsDoubleCPU: [34.9029 74.4581 93.3138 104.2198 108.8269]
    gflopsDoubleGPU: [72.1915 365.3399 522.5142 628.3013 681.8810]

Документация по инструментам для параллельных вычислений

Поддержка