CC BY
55
ИЗБРАННЫЕ ПУБЛИКАЦИИ I СЕССИИ НАУЧНОЙ ШКОЛЫ «ТЕХНОЛОГИИ ВЫСОКОПРОИЗВОДИТЕЛЬНЫХ ВЫЧИСЛЕНИЙ
И СИСТЕМ»
УДК 681.3.069, 681.324
ОСОБЕННОСТИ АДАПТАЦИИ ВЫЧИСЛИТЕЛЬНЫХ АЛГОРИТМОВ ПОД ПАРАЛЛЕЛЬНУЮ АРХИТЕКТУРУ ГРАФИЧЕСКИХ АКСЕЛЕРАТОРОВ
С.В. Ковальчук, С.М. Вишняков
В работе обсуждаются вопросы отображения вычислительных алгоритмов на параллельную архитектуру GPU-акселератора. Описываются важнейшие характеристики рассматриваемой аппаратной архитектуры и их влияние на эффективность работы параллельных алгоритмов. Предлагается модель оценки производительности параллельного алгоритма, использующего вычислительные ресурсы графического акселератора. В качестве примера рассматривается задача оптимизации параметрической аппроксимации двумерной функции методом случайного поиска.
Ключевые слова: графический акселератор, GPGPU, параллельные алгоритмы.
Введение
В последние годы интенсивное развитие получила специфическая отрасль высокопроизводительных вычислений - расчеты на системах с параллельными акселераторами как дополнительными устройствами, принимающими на себя существенную часть вычислительной нагрузки работающего приложения [1]. К таким устройствам относятся, в частности, GPU (Graphic Processor Шк)-устройства, предназначенные для работы с графикой. Современные GPU, по сути, являются многопроцессорными системами SIMD-архитектуры с достаточно высокой (до 0,5 Тфлопс) пиковой производительностью. По сравнению с традиционными архитектурами (например, кластерами), они обладают несопоставимо низкой характеристикой «цена/производительность», что стимулирует интерес к использованию GPU не только для обработки графической информации, но и для решения произвольных вычислительных задач [2]. В частности, к таким задачам относятся сортировки больших объемов данных [3], решения систем линейных уравнений [4], уравнений математической физики [5-7] и пр.
Параллельное программирование для GPU традиционно использует графический программный инструментарий, например, на основе библиотеки OpenGL [8] или DirectX [9]. Для удобства работы с ними применяются пакеты более высокого уровня, например, Gg [10], HLSL [11] и OpenGL Shading Language [12]. Однако реализация вычислительных задач общего плана требует использования более универсальных средств, облегчающих отображение задач различного типа на архитектуру GPU-устройств и позволяющих абстрагироваться от особенностей работы графических алгоритмов. Например, компанией nVidia разрабатывается набор библиотек CUDA SDK [13], предназначенный для использования вычислительных мощностей видеокарт и, в перспективе, специализированных устройств для решения задач общего назначения.
Несмотря на высокую производительность, вычислительные акселераторы на базе GPU являются специфическим классом многопроцессорных систем, характеризуемым существенными ограничениями на масштабируемость, использование памяти и управление данными и задачами. Потому представляет интерес изучение способов отображения различных типов вычислительных алгоритмов на GPU-архитектуру с целью выявления ключевых особенностей этого процесса и факторов, влияющих на получаемую
2
производительность. В данной работе рассматриваются особенности этого процесса на примере алгоритма параметрической оптимизации методом случайного поиска.
Особенности архитектуры GPU-устройств и средства отображения вычислительных алгоритмов
Акселератор на базе GPU представляет собой набор вычислительных узлов (мультипроцессоров), состоящих из некоторого числа арифметико-логических устройств (АЛУ). В любой момент времени все АЛУ одного мультипроцессора выполняют одинаковую последовательность инструкций над разными наборами данных, расположенных в памяти GPU (т.е. мультипроцессоры имеют SIMD-архитектуру). В данной работе в качестве примера рассматривается акселератор nVidia GeForce 8800 GTX, который имеет 16 мультипроцессоров по 8 АЛУ в каждом.
Поскольку в SIMD-устройствах каждый мультипроцессор конфигурируется определенным образом для выполнения заданной последовательности команд на множестве входных данных, то перед разработчиком стоит задача формирования последовательности инструкций, решающей поставленную задачу конфигурации устройства и передачи в устройство потока входных данных. При использовании nVidia CUDA SDK эта процедура выглядит следующим образом: разработчик описывает ядро (процедуру, которая будет исполняться на GPU) и при запуске ядра на GPU задает ему конфигурацию на основе описанной ниже иерархии. Каждый поток, физически выполняющийся на АЛУ мультипроцессора, исполняет инструкции, описанные в ядре. При этом благодаря SIMD-архитектуре на каждом мультипроцессоре несколько потоков параллельно выполняют одну и ту же последовательность инструкций. Логически потоки объединяются в блоки, ограничивающие возможность обмена данными между потоками (обмен данными через общую память только внутри одного блока). Потоки из одного блока выполняются на одном и том же мультипроцессоре. Таким образом, перед запуском ядра разработчик задает его конфигурацию (grid) - размер блока и общее число блоков. Для более удобного структурирования конкретных задач в CUDA предусмотрена возможность конфигурирования блоков и потоков в двумерные сетки. Данная возможность обеспечивает удобство при обработке структурированных данных. Отметим еще один важный структурный элемент конфигурации ядра - основу (warp). Под основой (в терминологии nVidia) подразумевается набор потоков, инструкции которых выполняются на мультипроцессоре одновременно. Размер основы определяется архитектурой конкретного GPU-устройства, и размер блока, определяемый при конфигурировании ядра, должен быть пропорционален размеру основы.
Рис. 1. Схема отображения вычислительного алгоритма на ОРи-устройство
Ключевым моментом отображения вычислительного алгоритма на ОРИ-архитектуру является определение ядра, соответствующего решению поставленной задачи. При этом следует учитывать, что каждый из потоков, исполняющих код ядра, для получения входных данных должен обращаться к своей области памяти устройства. На
рис. 1 представлена общая схема преобразования последовательного циклического алгоритма в его GPU-реализацию.
При исполнении на CPU-системе алгоритм многократно последовательно вычисляет некоторую функцию f (X) на различных наборах данных Xt. При исполнении на GPU указанная функция преобразуется к виду f'(X) , удобному для исполнения в рамках одного потока, и затем вычисляется на GPU параллельно, на разных наборах данных. Разумеется, данный подход является эффективным лишь в том случае, когда итерации цикла, вычисляющего f (X), независимы. Заметим также, что при описании вычислительного ядра необходимо учитывать ряд важных особенностей работы GPU-устройств. Во-первых, SIMD-архитектура накладывает ограничения на использование ветвящихся конструкций (if/else). Это связано с тем, что потоки, исполняющиеся на одном и том же мультипроцессоре, должны оперировать одной и той же последовательностью инструкций, что может нарушаться при несовпадении условий ветвления у разных потоков. При нарушении этого ограничения может происходить сильное ухудшение производительности. Во-вторых, имеет место особенность работы GPU-устройств с памятью - определение ядра таким образом, чтобы одновременно выполняющиеся потоки читали данные из соседних участков памяти, существенно увеличивает производительность. Кроме того, необходимо учитывать трудоемкость обращения к памяти GPU-устройства. В-третьих, многое зависит от архитектуры конкретного GPU-устройства, на котором выполняется вычислительная задача. Более подробно эти и остальные особенности GPU-архитектуры рассмотрены в [13].
Отображение алгоритма параметрической оптимизации методом случайного поиска на архитектуру GPU-устройства
В качестве примера вычислительного алгоритма, отображаемого на GPU-архитектуру, в работе была выбрана задача параметрической аппроксимации двумерной функции с использованием метода случайного поиска [15].
В процессе работы данного алгоритма оптимизации производится случайный выбор направления в пространстве поиска, V е Rn. В работе используется вариация метода адаптивного случайного поиска с линейной тактикой:
Vn = Vn-1 +AVn $(1), (1)
где $(11 - единичный орт в случайном направлении j. В этом случае шаг AVn в случайно выбранном направлении повторяется до тех пор, пока он не перестанет приводить к уменьшению значения целевой функции J(V) :
fa$, если AJn 1 > 0,
AVn =\ ^ n(2)
n [AVn-1, если AJn-1 < 0.
Здесь a - начальный размер шага для нового направления. Для изменения длины шага в процессе работы алгоритма используется параметр релаксации. Решение об изменении шага принимается в зависимости от успешности предыдущей операции:
, , [у 1 |AV„_,I, в случае успеха,
AVn Ч IА V \ (3)
ly2 |AVn_11, в противном случае.
Для определения параметров релаксации используется следующее соотношение:
Yp Y21-p) = 1. (4)
Алгоритм имеет широкое применение, в частности, при решении задачи аппроксимации спектров морского волнения [16]. Одной из особенностей этой задачи является необходимость обработки больших массивов входных данных, что, с одной стороны,
приводит к большой трудоемкости вычислений, а с другой - позволяет эффективно использовать распараллеливание по данным. В результате эффективно используются вычислительные мощности GPU-устройств, поскольку для каждого набора входных данных решается одна и та же задача.
Анализ параллельной производительности
Анализ параллельной производительности проведен для алгоритма, основанного на декомпозиции по данным, при котором в каждом потоке выполняется подсчет значения целевой функции для своего набора входных данных. Предметом исследования являлось влияние конфигурации системы на получаемое ускорение, а также оценка составных частей времени работы алгоритма. Кроме того, был исследован эффект применения различных типов оптимизаций вычислительного ядра. Для оценки ускорения время работы с использованием GPU-устройства nVidia GeForce 8800 GTX сравнивалось со временем обработки тех же данных на CPU Intel Core 2 Duo 2.3 GHz.
На рис. 2 представлена зависимость получаемого ускорения от конфигурации GPU-устройства (число потоков в блоке, умноженное на число блоков) ядра.
64x4 64x8 64x16 64x24 64x32 64x48 64x64 конф игурация ядра
Рис. 2. Зависимость ускорения от конфигурации ядра
Рост ускорения при увеличении числа блоков от 4 до 16 объясняется тем, что в исследуемой системе 16 мультипроцессоров. Таким образом, при конфигурации в 4 и 8 блоков часть мультипроцессоров остается бездействующей. При дальнейшем увеличении числа блоков все мультипроцессоры загружены работой, тем не менее, ускорение продолжает расти. Этот факт объясняется особенностью обращения мультипроцессоров к памяти: если на одном мультипроцессоре, в соответствии с конфигурацией, выполняется одновременно несколько блоков, то, в то время как процессы одного блока ждут данные из памяти, на мультипроцессоре могут выполняться арифметические операции других блоков. Схема на рис. 3 поясняет данное утверждение: первоначально при увеличении числа блоков (а, значит, и количества входных данных) общее время работы практически не меняется.
1П I II II I I
В2 I I II II I I
В2 I ^^
~В4 I ^^ II II I ]
-Время-►
I I Ожидание данных из памяти I I Последовательность арифметических операций
Рис. 3. Одновременное выполнение блоков на мультипроцессоре ОРЫ-устройства
Аналогичный эффект имеет место при последовательности арифметических вычислений, в процессе которых результат одной операции используется в последую-
щей [13]. Задержки, возникающие в данной ситуации, компенсируются аналогичным образом. В тот момент, когда число блоков достигает некоторого критического значения, определяемого максимальным числом блоков, которое возможно запустить на одном мультипроцессоре при заданном размере блока и заданном ядре, время работы ядра резко увеличивается. Это число зависит от того, сколько ресурсов мультипроцессора используется одним потоком заданного ядра. Зависимость максимального числа блоков от числа регистров и размера разделяемой памяти, используемой ядром, можно найти в nVidia Occupancy Calculator [13]. На рис. 5 представлены результаты, полученные в процессе измерения составных частей времени работы алгоритма на GPU-устройстве:
Т Tmalloc+Tmemcpyin+Tkernel+Tmemcpyout+Tfree- (5)
Время Т состоит не только из времени выполнения ядра Tkernei, но и времени, затрачиваемого на выделение Tmalloc и освобождение Tfree памяти на устройстве и копирования
данных Т
Tm
memcpyin-, 1 memcpyout
Рис. 5. Составные части времени работы
Время исполнения ядра практически не меняется с увеличением числа блоков от 4 до 32 - сначала из-за неполной загруженности устройства, затем - из-за описанного выше эффекта, в то время как при увеличении числа блоков до 48 происходит резкий скачок. Время, затрачиваемое на работу с памятью, растет при этом линейно (так как с увеличением числа блоков увеличивается и количество входных и выходных данных). Заметим, что это время в данном случае составляет существенную часть от общего времени. Таким образом, можно сделать вывод о том, что в задачах, использующих повторяющиеся наборы данных, следует избегать повторного копирования и выделения памяти там, где это возможно.
Изучим более подробно влияние конфигурации ядра на производительность алгоритма и попробуем получить теоретическую модель производительности. Для этого упростим ситуацию, исключив из рассматриваемого выше алгоритма недетерминированность, сведя его к многократному подсчету интеграла методом прямоугольников, и проведем ряд измерений с различными параметрами конфигурации ядра.
На рис. 6, а, представлена зависимость производительности ядра (количества обработанных данных, отнесенного ко времени работы ядра) от числа блоков в конфигурации ядра для различного числа потоков на блок. Вид зависимости производительности обусловливается ступенчатым характером зависимости времени работы ядра от тех же параметров (рис. 6, б). Ступенчатый характер этой зависимости объясняется описанным выше эффектом: при увеличении общего числа блоков время работы ядра не
изменяется до определенного момента благодаря тому, что операции добавляемых блоков выполняются во время «простоя» мультипроцессоров.
500
400
2 300
I-
о ю
224 -А-160
96 32
80 112 144
число блоков
Рис. 6. Характеристики ядра
Попытаемся вывести зависимость, представленную на рис. 7, аналитически, учитывая тот факт, что «насыщение» мультипроцессоров происходит в первый раз в тот момент, когда общее число блоков, приходящееся на мультипроцессор, требует больше ресурсов, чем имеется на мультипроцессоре. Иначе говоря, в этот момент все блоки не могут работать на мультипроцессоре одновременно. Таким образом, становится верным неравенство
nblocks * (blockSize / ^атр^е)
M
> mWarps * occupancy,
(6)
где nblocks - общее число блоков в конфигурации ядра; blockSize - число потоков в блоке, warpSize - число потоков в основе (определяется архитектурой GPU-устройства); mWarps - максимальное число основ, которое может одновременно работать на одном мультипроцессоре GPU-устройства данной архитектуры; occupancy - отношение максимального числа основ, которое можно запустить на одном мультипроцессоре, исходя из использования данным ядром ресурсов и количества доступных ресурсов на мультипроцессоре, к mWarps. Параметр occupancy можно подсчитать при помощи nVidia Occupancy Calculator [13], зная количество регистров и разделяемой памяти, используемых ядром.
Исходя из приведенных выше соображений, можно предложить следующую аналитическую формулу для времени работы ядра:
Tkernel = 10(blockSize) ■
N
nblocks ■ blockSize
nblocks ■ blockSize
= 10 (blockSize) ■ workPerThread
M ■ mWarps ■ occupancy ■ warpSize warpsPerMP
(7)
mWarps ■ occupancy
где t0 - «эффективное» время работы ядра, которое вычисляется как время работы ядра при полной загрузке мультипроцессора перед насыщением (фактически - высота первой ступени графика), отнесенное к числу элементарных операций, выполняемых ядром; N - общее число элементарных операций, выполняемых ядром; workPerThread -число элементарных операций, приходящихся на один поток; warpsPerMP - количество основ, приходящихся на один мультипроцессор при данной конфигурации ядра.
Вычислив значения t0(blockSize) на основе данных, представленных на рис. 7, можно, используя формулу (7), построить теоретические графики зависимости времени работы ядра от его конфигурации для другого числа элементарных операций, выполняемых ядром. Сравнение теоретической зависимости с зависимостью, полученной на практике, представлено на рис. 7.
200
100
0
16
48
16 48 80 112 144 176
число блоков
Рис. 7. Сравнение теоретической и практической зависимостей времени работы ядра от конфигурации
Как видно, теоретическая зависимость довольно хорошо совпадает с экспериментальными данными. Расхождения наблюдаются лишь в началах «ступеней» графика в случаях, когда блок содержит четное число основ. Причина этого расхождения является открытым вопросом.
Рис. 8. Зависимость t0(blockSize) и ее теоретический вид (а), зависимость occupancy(blockSize) (б)
Экспериментальные данные показывают, что параметр t0 в формуле (7) зависит от размера блока. Соответствующая экспериментальная зависимость представлена на рис. 8, а. Заметим, что зависимость параметра occupancy от размера блока (рис. 8, б) носит похожий характер, согласно [13]. Этот факт позволяет предположить, что параметры t0 и occupancy связаны линейным соотношением: 10 = ta ■ occupancy + tb. Приняв это предположение и вычислив значения ta и tb из экспериментальных данных, можно сравнить теоретическую зависимость t0(blockSize) с полученной экспериментально.
Заключение
Результаты вычислительных экспериментов на задаче аппроксимации климатических спектрах морского волнения показывают, что адаптация алгоритма параметрической оптимизации к архитектуре GPU позволяет получить ускорение более 30 раз по сравнению с реализацией без использования GPU. Учитывая, что стоимость простейшей кластерной системы на основе стандартных комплектующих, обеспечивающих такое же ускорение, примерно в 30 раз превышает стоимость рассматриваемого GPU-устройства, это подтверждает целесообразность применения GPU-устройств при решении вычислительных задач подобного класса. В планах дальнейших исследований -разработка методов отображения других алгоритмов на архитектуру GPU-устройств.
Литература
1. Hasle G., Lie K.-A., Quak E. Geometric Modelling, Numerical Simulation, and Optimization: Applied Mathematics at SINTEF. - Springer, 2007. - 558 p.
2. General-Purpose Computation Using Graphics Hardware. - Режим доступа: http://gpgpu.org/, свободный.
3. Purcell T.J., Donner C., Commarano M., Jensen H.W., Hanrahan P. Photon mapping on programmable graphic hardware // Proceeding of the ACM SIGGRAPH/EUROGRAPHICS Conference on Graphics Hardware. - Eurographics Association. - 2003. - Р. 41-50.
4. Goddeke M., Strzodka R., Turek S. Accelerating Double Precision FEM Simulations with GPUs // Proceeding of ASIM 2005. - 18th Symposium on Simulation Technique.-2005. - Р. 139-144.
5. Hagen T.R., Henriksen M.O., Hjelmervik J.M., Lie K.-A. Using the graphic processor as a high-performance computational engine for solving system of hyperbolic conservation low // Geometric Modelling, Numerical Simulation, and Optimization: Applied Mathematics at SINTEF. - Springer, 2007. - Р. 211-264.
6. Hagen T.R., Hjelmervik J.M., Lie K.-A., Natvig J.R., Henriksen M.O. Visual simulation of shallow-water waves // Simulation Practice and Theory. Special Issue on Programmable Graphics Hardware. - 2005. - V.13. - № 9. - Р. 716-726.
7. Hagen T.R., Lie K.-A., Natvig J.R. Solving the Euler equation on graphical processing units // Computational Science - ICCS 2006: 6th International Conference, Reading, UK, May 28-31, 2006, Proceedings, Part IV, volume 3994 of Lecture Notes in Computational Science (LNCS). - Springer Verlag, 2006. - Р. 220-227.
8. OpenGL - The Industry Standard for High Performance Graphics. - Режим доступа http://www.opengl.org/, свободный.
9. Microsoft's DirectX site. - Режим доступа: http://www.microsoft.com/directx, свободный.
10. Fernando R., Kilgard M.J. The Cg Tutorial: The Definitive Guide to Programmable Real-Time Graphics. - Adisson-Wesley Longman Publishing Co., 2003.
11. GPUBench: How much does your GPU bench. - Режим доступа: http://graphics.stanford.edu/projects/gpubench/, свободный.
12. Rost R.J. OpenGL Shading Language. - Adisson-Wesley Longman Publishing Co, 2004.
13. NVIDIA CUDA Compute Unified Device Architecture Programming Guide. Ver 1.0. June 2007. - NVIDIA Corporation, 2007.
14. Адинец А.В., Сахарных Н.А. О программировании вычислений общего назначения на графических процессорах // Научный сервис в сети Интернет: многоядерный компьютерный мир. 15 лет РФФИ: Труды Всероссийской научной конференции.-М.: Издательство МГУ, 2007. - С. 249-256.
15. Растригин Л.А. Адаптация сложных систем - Рига: Зинатне, 1981. - 375 с.
16. Boukhanovsky A.V., Lopatoukhin L.J., Guedes Soares C. Spectral wave climate of the North Sea // Applied Ocean Research. - 2007. - Vol. 29. - Р. 146-154.
Ковальчук Сергей Валерьевич
Вишняков Сергей Михайлович
Санкт-Петербургский государственный университет информационных технологий, механики и оптики, м.н.с., kovalchuk@mail.ifmo.ru Санкт-Петербургский государственный университет информационных технологий, механики и оптики, студент, sergey.vishnyakov@gmail.com
