Улучшайте Производительность Небольших Матричных проблем на графическом процессоре с помощью PAGEFUN

Этот пример показывает, как использовать pagefun, чтобы улучшать производительность применения большого количества независимых вращений и переводов в объекты в 3-D среде. Это типично для области значений проблем, которые включают большой пакет вычислений на небольших массивах.

Графические процессоры являются самыми эффективными при выполнении вычислений на очень больших матрицах. В MATLAB® это обычно достигается путем векторизации кода, чтобы максимизировать работу, сделанную в каждой инструкции. Когда у вас есть большой набор данных, но ваши вычисления разделены на многие маленькие операции над матрицей, это может быть сложно, чтобы максимизировать производительность путем выполнения одновременно на многих сотнях ядер графического процессора.

arrayfun и функции bsxfun позволяют скалярным операциям быть выполненными параллельно на графическом процессоре. Функция pagefun добавляет возможность проведения операций над матрицей в пакете похожим способом. Функция pagefun доступна в Parallel Computing Toolbox™ для использования с gpuArrays.

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

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

function paralleldemo_gpu_pagefun

Настройте карту

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

numObjects = 1000;
roomDimensions = [50 50 5]; % Length * breadth * height in meters

Мы представляем положения и ориентации с помощью векторов 3 на 1 T и 3х3 матрицы вращения R. То, когда у нас есть N этих, преобразовывает, мы упаковываем переводы в 3 N матрицей и вращения в 3 3 N массивом. Следующая функция инициализирует N, преобразовывает со случайными значениями, обеспечивая структуру, как выведено:

    function Tform = randomTransforms(N)
        Tform.T = zeros(3, N);
        Tform.R = zeros(3, 3, N);
        for i = 1:N
            Tform.T(:,i) = rand(3, 1) .* roomDimensions';
            % To get a random orientation, we can extract an orthonormal
            % basis for a random 3-by-3 matrix.
            Tform.R(:,:,i) = orth(rand(3, 3));
        end
    end

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

Map = randomTransforms(numObjects);
Robot = randomTransforms(1);

Определите уравнения

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

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

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

Rrel(:,:,i) = Rbot' * Rmap(:,:,i)
Trel(:,i) = Rbot' * (Tmap(:,i) - Tbot)

Выполните много матриц, преобразовывает на центральном процессоре с помощью цикла for

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

    function Rel = loopingTransform(Robot, Map)
        Rel.R = zeros(size(Map.R), 'like', Map.R); % Initialize memory
        Rel.T = zeros(size(Map.T), 'like', Map.T); % Initialize memory
        for i = 1:numObjects
            Rel.R(:,:,i) = Robot.R' * Map.R(:,:,i);
            Rel.T(:,i) = Robot.R' * (Map.T(:,i) - Robot.T);
        end
    end

Ко времени вычисление мы используем функцию timeit, которая вызовет loopingTransform многократно, чтобы получить среднюю синхронизацию. Поскольку это требует функции без аргументов, мы используем синтаксис @(), чтобы создать анонимную функцию правильной формы.

cpuTime = timeit(@()loopingTransform(Robot, Map));
fprintf('It takes %3.4f seconds on the CPU to execute %d transforms.\n', ...
        cpuTime, numObjects);
It takes 0.0104 seconds on the CPU to execute 1000 transforms.

Попробуйте тот же код по графическому процессору

Запускать этот код по графическому процессору - просто вопрос копирования данных в gpuArray. Когда MATLAB столкнется с данными, хранимыми на графическом процессоре, это запустит любой код с помощью него на графическом процессоре, пока это поддерживается.

gMap.R = gpuArray(Map.R);
gMap.T = gpuArray(Map.T);
gRobot.R = gpuArray(Robot.R);
gRobot.T = gpuArray(Robot.T);

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

fprintf('Computing...\n');
gpuTime = gputimeit(@()loopingTransform(gRobot, gMap));

fprintf('It takes %3.4f seconds on the GPU to execute %d transforms.\n', ...
        gpuTime, numObjects);
fprintf(['Unvectorized GPU code is %3.2f times slower ',...
    'than the CPU version.\n'], gpuTime/cpuTime);
Computing...
It takes 0.5588 seconds on the GPU to execute 1000 transforms.
Unvectorized GPU code is 53.90 times slower than the CPU version.

Процесс пакетной обработки с помощью pagefun

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

    function Rel = pagefunTransform(Robot, Map)
        Rel.R = pagefun(@mtimes, Robot.R', Map.R);
        Rel.T = Robot.R' * bsxfun(@minus, Map.T, Robot.T);
    end

gpuPagefunTime = gputimeit(@()pagefunTransform(gRobot, gMap));
fprintf(['It takes %3.4f seconds on the GPU using pagefun ',...
    'to execute %d transforms.\n'], gpuPagefunTime, numObjects);
fprintf(['Vectorized GPU code is %3.2f times faster ',...
    'than the CPU version.\n'], cpuTime/gpuPagefunTime);
fprintf(['Vectorized GPU code is %3.2f times faster ',...
    'than the unvectorized GPU version.\n'], gpuTime/gpuPagefunTime);
It takes 0.0008 seconds on the GPU using pagefun to execute 1000 transforms.
Vectorized GPU code is 13.55 times faster than the CPU version.
Vectorized GPU code is 730.18 times faster than the unvectorized GPU version.

Объяснение

Первое вычисление было вычислением вращений. Это включило матрицу, умножаются, который переводит в функциональный mtimes (*). Мы передаем это pagefun наряду с двумя наборами вращений, которые будут умножены:

Rel.R = pagefun(@mtimes, Robot.R', Map.R);

Robot.R' является 3х3 матрицей, и Map.R является 3 3 N массивом. Функция pagefun совпадает с каждой независимой матрицей от карты до того же вращения робота и дает нам необходимые 3 3 N вывод.

Вычисление перевода также включает матрицу, умножаются, но нормальные правила умножения матриц позволяют этому прибывать вне цикла без любых изменений. Однако это также включает вычитание Robot.T от Map.T, которые являются различными размерами. Поскольку это поэлементно операция, мы можем использовать bsxfun, чтобы подойти размерности таким же образом, как pagefun сделал для вращений:

Rel.T = Robot.R' * bsxfun(@minus, Map.T, Robot.T);

На этот раз мы должны были использовать поэлементный оператор, который сопоставляет с функциональным minus (-).

Более усовершенствованная векторизация графического процессора - Решение "потерянного робота" проблема

Если наш робот был в неизвестной части карты, он может использовать алгоритм глобального поиска, чтобы определить местоположение себя. Алгоритм протестировал бы много возможных местоположений путем выполнения вышеупомянутого вычисления и поиска хорошего соответствия между объектами, замеченными датчиками робота и что это будет ожидать видеть в том положении.

Теперь у нас есть несколько роботов, а также несколько объектов. N объекты и роботы M должен дать нам, N*M преобразовывает. Чтобы отличить 'пробел робота' от 'объектного пространства', мы используем 4-ю размерность для вращений и 3-е для переводов. Это означает, что наши вращения робота будут 3 3 1 M, и переводы будут 3 1 M.

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

numRobots = 10;
Robot = randomTransforms(numRobots);
Robot.R = reshape(Robot.R, 3, 3, 1, []); % Spread along the 4th dimension
Robot.T = reshape(Robot.T, 3, 1, []); % Spread along the 3rd dimension
gRobot.R = gpuArray(Robot.R);
gRobot.T = gpuArray(Robot.T);

Наше новое цикличное выполнение преобразовывает функцию, требует двух вложенных циклов, чтобы циклично выполниться по роботам, а также по объектам.

    function Rel = loopingTransform2(Robot, Map)
        Rel.R = zeros(3, 3, numObjects, numRobots, 'like', Map.R);
        Rel.T = zeros(3, numObjects, numRobots, 'like', Map.T);
        for i = 1:numObjects
            for j = 1:numRobots
                Rel.R(:,:,i,j) = Robot.R(:,:,1,j)' * Map.R(:,:,i);
                Rel.T(:,i,j) = ...
                    Robot.R(:,:,1,j)' * (Map.T(:,i) - Robot.T(:,1,j));
            end
        end
    end

cpuTime = timeit(@()loopingTransform2(Robot, Map));
fprintf('It takes %3.4f seconds on the CPU to execute %d transforms.\n', ...
        cpuTime, numObjects*numRobots);
It takes 0.1493 seconds on the CPU to execute 10000 transforms.

Для наших синхронизаций графического процессора мы используем tic и toc на этот раз, потому что в противном случае вычисление заняло бы слишком много времени. Это будет достаточно точно в наших целях. Гарантировать любую стоимость, сопоставленную созданием выходных данных, включено, мы вызываем loopingTransform2 с одной выходной переменной, как timeit и gputimeit делают по умолчанию.

fprintf('Computing...\n');
tic;
gRel = loopingTransform2(gRobot, gMap); %#ok<NASGU> Suppress unused variable warning
gpuTime = toc;

fprintf('It takes %3.4f seconds on the GPU to execute %d transforms.\n', ...
        gpuTime, numObjects*numRobots);
fprintf(['Unvectorized GPU code is %3.2f times slower ',...
    'than the CPU version.\n'], gpuTime/cpuTime);
Computing...
It takes 7.0564 seconds on the GPU to execute 10000 transforms.
Unvectorized GPU code is 47.26 times slower than the CPU version.

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

Новая версия pagefun должна включить оператор transpose, а также mtimes в вызов pagefun. Нам также нужны к squeeze транспонированные ориентации робота, чтобы поместить распространение по роботам в 3-ю размерность, совпадать с переводами. Несмотря на это, получившийся код значительно более компактен.

    function Rel = pagefunTransform2(Robot, Map)
        Rt = pagefun(@transpose, Robot.R);
        Rel.R = pagefun(@mtimes, Rt, Map.R);
        Rel.T = pagefun(@mtimes, squeeze(Rt), ...
                        bsxfun(@minus, Map.T, Robot.T));
    end

Еще раз pagefun и bsxfun расширяют размерности соответственно. Таким образом, где мы умножаемся 3 на 3 на 1 M матричным Rt с 3 на 3 N на 1 матричный Map.R, мы получаем 3 3 N M матрицей.

gpuPagefunTime = gputimeit(@()pagefunTransform2(gRobot, gMap));
fprintf(['It takes %3.4f seconds on the GPU using pagefun ',...
    'to execute %d transforms.\n'], gpuPagefunTime, numObjects*numRobots);
fprintf(['Vectorized GPU code is %3.2f times faster ',...
    'than the CPU version.\n'], cpuTime/gpuPagefunTime);
fprintf(['Vectorized GPU code is %3.2f times faster ',...
    'than the unvectorized GPU version.\n'], gpuTime/gpuPagefunTime);
It takes 0.0025 seconds on the GPU using pagefun to execute 10000 transforms.
Vectorized GPU code is 59.97 times faster than the CPU version.
Vectorized GPU code is 2834.45 times faster than the unvectorized GPU version.

Заключение

Функция pagefun поддерживает много 2D операций, а также большинство скалярных операций, поддержанных arrayfun и bsxfun. Вместе, эти функции позволяют вам векторизовать область значений вычислений, включающих матричную алгебру и манипуляцию с массивами, устраняя циклы for потребности и делая огромное увеличение производительности.

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

end