Например, Бобцов

АНАЛИЗ ЭФФЕКТИВНОСТИ ПРИМЕНЕНИЯ ТЕХНОЛОГИИ CUDA ДЛЯ РЕШЕНИЯ СИСТЕМ ЛИНЕЙНЫХ УРАВНЕНИЙ С ТРЕХДИАГОНАЛЬНЫМИ МАТРИЦАМИ В ЗАДАЧАХ РАСЧЕТА ЦЕН ОПЦИОНОВ

20 Д. Н. Шинкарук, Ю. А. Шполянский, М. С. Косяков
УДК 519.688
Д. Н. ШИНКАРУК, Ю. А. ШПОЛЯНСКИЙ, М. С. КОСЯКОВ
АНАЛИЗ ЭФФЕКТИВНОСТИ ПРИМЕНЕНИЯ ТЕХНОЛОГИИ CUDA ДЛЯ РЕШЕНИЯ СИСТЕМ ЛИНЕЙНЫХ УРАВНЕНИЙ С ТРЕХДИАГОНАЛЬНЫМИ МАТРИЦАМИ В ЗАДАЧАХ РАСЧЕТА ЦЕН ОПЦИОНОВ
Исследованы возможности графического процессора при решении линейных систем с трехдиагональными матрицами. Показано, что целесообразно формировать матрицы непосредственно в глобальной памяти видеокарты, при этом достигается более чем 20-кратное ускорение вычислений по сравнению с однопоточным расчетом. С учетом обмена данными между оперативной памятью и памятью видеокарты достигнуто ускорение 5—8 раз при использовании механизма отображаемой памяти.
Ключевые слова: CUDA, GPU общего назначения, система линейных уравнений, параллельная циклическая редукция, метод прогонки.
Введение. В последние годы все более широко для высокопроизводительной обработки данных применяются дополнительные ускорители, принимающие на себя часть вычислительной нагрузки. В частности, к таким ускорителям относятся видеокарты с графическими процессорами (ГП), изначально предназначенные для использования в качестве устройств обработки графики. Производительность этих устройств растет с каждым годом, а реализации алгоритмов обработки данных на ГП позволяют значительно повысить скорость вычислений по сравнению с использованием центрального процессора (ЦП) [1].
Результаты многочисленных исследований технологии CUDA (Compute Unified Device Architecture), позволяющей производить массово-параллельные вычисления с использованием ГП компании NVIDIA, показали значительное увеличение скорости обработки данных в не связанных с графическими приложениями областях [2—5].
В настоящей работе рассматривается задача расчета теоретических цен опционов как численного решения дифференциального уравнения в частных производных (ДУЧП), которое часто сводят к последовательному решению систем линейных алгебраических уравнений (СЛАУ) с трехдиагональными матрицами [6]. Целью работы является анализ эффективности применения технологии CUDA для ускорения решения СЛАУ с трехдиагональными матрицами на ГП и выявление условий, при которых достигается наибольшее преимущество ГП перед ЦП. При сравнении производительности таких вычислений использованы современные ЦП Intel и ГП компании NVIDIA. Детально изучены особенности применения различных методов обмена данными между оперативной памятью компьютера (ОЗУ) и памятью видеокарты, сформулированы соответствующие рекомендации.
Особенности архитектуры CUDA. Технология CUDA основана на следующей концепции: ГП выступает в роли массово-параллельного сопроцессора к ЦП и обладает своей иерархией памяти [7]. Наиболее востребованы разделяемая память видеокарты, располагающаяся на кристалле ГП и обладающая малым временем доступа, но малым объемом, и сравнительно „медленная“ глобальная память, имеющая наибольший объем.
Программа, написанная с использованием технологии CUDA, состоит из двух частей: последовательного кода, выполняющегося на ЦП, и параллельного кода, который выполняется на ГП. Параллельная часть программы, называемая ядром (kernel), исполняется одновременно большим количеством потоков (threads), которые организуются в блоки. С использованием последовательного кода подготавливаются и пересылаются данные из ОЗУ компью-
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2012. Т. 55, № 10

Анализ эффективности применения технологии CUDA

21

тера в память видеокарты, после чего осуществляется запуск ядра. По окончании работы ядра результаты передаются обратно в оперативную память.
Такие особенности архитектуры и программной организации работы ГП приводят к тому, что алгоритмы, предназначенные для работы на ЦП, не могут быть успешно применены для графического ускорителя. ГП требует применения специализированных параллельных алгоритмов обработки данных.
Методы решения СЛАУ. При организации вычислений на центральном процессоре СЛАУ обычно решают методом прогонки или LU-декомпозиции [8]. Алгоритм решения строго последователен и не может быть эффективно распараллелен. Поэтому при решении СЛАУ на ГП необходимы другие методы. В настоящей работе применен метод параллельной циклической редукции (ПЦР, Parallel Cyclic Reduction), который характеризуется большим объемом вычислений, но обладает высокой степенью параллелизма [9].
На рис. 1 представлена схема метода ПЦР. СЛАУ обрабатывается на ГП одним блоком, при этом расчетом каждого неизвестного xi занимается отдельный поток Pi. Метод выполняется пошагово. На шаге j производится вычисление новых коэффициентов Pij для уравнений системы с использованием коэффициентов, рассчитанных на предыдущем шаге j – 1. Взаимодействие потоков осуществляется через разделяемую память. На последнем этапе вычислений (на рис. 1 обозначен ромбами) производится итоговый расчет неизвестных xi. Работа алгоритма требует log2(N) шагов, где N — количество уравнений в СЛАУ.
x1 x2 x3 x4 x5 x6 x7

x1 x2 x3 x4 x5 x6 x7

P13 P23 P33 P43 P53 P63 P73

P12 P22 P32 P42 P52 P62 P72

P11 P21 P31 P41 P51 P61 P71
P1 P2 P3 P4 P5 P6 P7 Рис. 1
Вычислительная сложность ПЦР пропорциональна Nlog2(N), в то время как сложность метода прогонки линейно зависит от количества уравнений в системе. Но при массовопараллельном исполнении каждый отдельный поток выполняет лишь небольшую часть общего объема вычислений, пропорциональную log2(N).
Решение одной системы. Расчеты цен опционов и измерения времени вычислений проводились на персональном компьютере (ПК) со следующими параметрами: ЦП Intel Core i5 3,8 ГГц; 4 ГБ ОЗУ; ГП NVIDIA GTX 580; ОС Ubuntu 11.04 (64 бит); версия драйвера CUDA 270.41.19 (64 бит).
При организации расчетов на ГП в ОЗУ компьютера формировались 4 массива данных, соответствующие трем диагоналям матрицы и правому столбцу СЛАУ. Для матрицы СЛАУ
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2012. Т. 55, № 10

22 Д. Н. Шинкарук, Ю. А. Шполянский, М. С. Косяков

обеспечивалось требование диагонального преобладания. Массивы передавались в память

видеокарты, после чего запускалось ядро для обработки данных на ГП. Вектор результатов

копировался из памяти устройства обратно в ОЗУ. Измеряемыми характеристиками были

полное время Tf (full) от начала передачи данных на видеокарту до момента возвращения результатов в ОЗУ и чистое время вычислений на видеокарте, т.е. время выполнения ядра Tk (kernel). Согласно идеологии ПЦР, число запускаемых потоков выполнения совпадало с чис-

лом уравнений N. Для ускорения расчетов все потоки запускались в одном блоке, а требуемые

данные располагались в разделяемой памяти.

Время решения одной СЛАУ

Таблица 1

N 2 4 8 16 32 64 128 256 512 1024

Tk, мс

4 4 5 6 7 7 8 9 12 23

Tf, мс

117 113 119 117 119 119 119 121 125 263

TCPU, мс

1

1 2 2 2 3 5 11 19 39

Как видно из табл. 1, время работы ядра Tk на ГП сопоставимо со временем решения

системы уравнений методом прогонки на ЦП TCPU, однако для СЛАУ с N ≤ 128 превышает

его. Связано это с тем, что при массово-параллельных вычислениях на ГП существуют затра-

ты на создание и инициализацию потоков выполнения. Для CUDA время инициализации

практически не зависит от количества создаваемых потоков, поэтому наблюдается более мед-

ленный рост Tk по сравнению с TCPU. Важно отметить, что значение TCPU линейно зависит от числа уравнений N в решаемой СЛАУ.

Значение полного времени Tf слабо изменяется с ростом размерности СЛАУ N, но в десятки раз превышает Tk и TCPU (см. табл. 1). На рис. 2 представлена временная диаграмма работы видеокарты: передача данных из ОЗУ в память видеокарты (1) и обратно (5); работа яд-

ра (3); простои (2, 4). Время работы устройства составило 7,2 % от общего времени. Длитель-

ные простои (белые области) связаны с отработкой драйвера устройства, планированием про-

цессов в системе, синхронизацией и пр.

12

3 45

0 116 Тf, мс Рис. 2
Достичь преимущества ГП перед ЦП возможно, увеличив занятость устройства и эффективно используя его вычислительные ресурсы. Одним из способов достижения этой цели является параллельное решение большого числа СЛАУ.
При анализе погрешностей вычислений рассматривалась невязка решения СЛАУ. Наибольшее наблюдаемое значение невязки для метода ПЦР составило ∆ = 8,5⋅10–5, в то время как для метода прогонки (МП) ∆ = 7,8⋅10–6 при N = 300. Причина такого различия заключается в разном количестве арифметических операций, производимых над элементом данных. В реализации метода ПЦР осуществляется несколько больше операций умножения по сравнению с МП. Важно отметить, что для ПЦР количество арифметических операций определяется размерностью системы N: при увеличении числа уравнений в СЛАУ возрастает количество шагов решения. В наших экспериментах для N = 8 невязка составила ∆ = 8⋅10–8, а для N = 64 уже ∆ = 1,2⋅10–6. Эти значения не зависят того, где запускаются методы расчета: на ЦП или ГП, следовательно, использование ГП не приводит к дополнительным погрешностям.
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2012. Т. 55, № 10

Анализ эффективности применения технологии CUDA

23

Параллельное решение многих систем. На рис. 3, а представлена зависимость време-
ни Tk и Tf расчетов на ГП от количества M одновременно обрабатываемых СЛАУ при N = 256. Как видно, при небольшой нагрузке значения Tk (1) и Tf (2) мало меняются, что связано с недостаточной загрузкой устройства. При M = 64 достигается предел одновременно запущен-
ных потоков выполнения (16 384 потоков для видеокарты NVIDIA GTX 580), после чего вре-
мя вычислений растет линейно. Стоит отметить, что время расчета на ЦП TCPU (3) изменяется линейно во всей области рассматриваемой нагрузки.
На рис. 3, б представлена зависимость скорости решения СЛАУ V (определяемой как
отношение числа M решаемых СЛАУ ко времени их расчета T) от нагрузки на систему. Рост
скорости выполнения расчетов Vf (1) на ГП с увеличением M происходит за счет сокрытия задержек доступа к памяти при одновременной работе большого количества потоков выпол-
нения, при этом скорость ЦП (кривая 2) практически не изменяется.
а) б)

Т, мс

V, 105 сист./с

100
10 3
1 0,1 2 0,01 1 0,001
20 22 24 26 28 210 212 214 М

4 1
3
22 1
0 20 22 24 26 28 210 212 214 М

Рис. 3

При решении большого числа СЛАУ объем передаваемых данных между ОЗУ и памятью видеокарты значительно возрастает, что существенно сказывается на общем времени расчетов (табл. 2). Так, при увеличении числа решаемых систем M возможно сокращение чистого времени вычислений на ГП перед ЦП в 25 раз. Однако при пересылке данных между ОЗУ и памятью видеокарты преимущество в скорости расчетов с использованием ГП составляет 2,6 раза, 90 % общего времени работы устройства отводится на передачу данных. Один из возможных выходов — такая организация вычислительного процесса, при которой матрицы СЛАУ формируется непосредственно в глобальной памяти видеокарты.

Таблица 2 Время решения различного количества систем уравнений

M
1 16 256 1024 4096

Tk, мс
0,027 0,026 0,065 0,231 0,893

Tf, мс
0,167 0,171 0,701 2,418 8,763

TCPU, мс
0,005 0,08 1,263 5,255 23,077

TCPU / Tk
0,18 3,01 19,41 22,70 25,84

TCPU / Tf
0,03 0,47 1,80 2,17 2,63

Видеокарты последних поколений позволяют осуществлять асинхронную передачу данных из ОЗУ в память видеокарты и обратно одновременно с исполнением кода ядра [10]. Таким образом достигается эффективное функционирование приложения при наличии большого объема входных данных. Одно из ограничений при этом — необходимость обозначения области ОЗУ, в которой хранятся передаваемые данные, как невыгружаемой области памяти. При реализации этого метода входные данные „вручную“ разбиваются на фрагменты таким образом, чтобы время передачи одного фрагмента приближалось ко времени его обработки.

ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2012. Т. 55, № 10

24 Д. Н. Шинкарук, Ю. А. Шполянский, М. С. Косяков

Таблица 3 Полное время решения СЛАУ при синхронном и асинхронном обмене данными

M
1 16 256 1024 4096

Tfs, мс
0,167 0,171 0,701 2,418 8,763

Tfa, мс
0,041 0,060 0,389 1,374 5,141

TCPU, мс
0,005 0,08 1,263 5,255 23,077

TCPU / Tfa
0,12 1,32 3,25 3,82 4,49

Применение асинхронной передачи данных позволило уменьшить общее время расчетов Tfa более чем в 1,5 раза по сравнению со случаем синхронного обмена Tfs (см. табл. 3). Для небольшого числа решаемых систем M наблюдается уменьшение времени расчетов в 3—4
раза, поскольку при асинхронной передаче данных исчезает значительная часть простоев (см.
рис. 2). Однако необходимость „ручного“ разбиения входных данных на части весьма услож-
няет применение данного метода.
Нами исследована эффективность применения механизма отображаемой (mapped) памя-
ти для прямого доступа к данным в ОЗУ из исполняемого на ГП ядра. В этом случае соответ-
ствующая область ОЗУ также должна быть обозначена как невыгружаемая. Результаты сравнения полного времени выполнения для синхронного Tfs, асинхронного Tfa обмена данными и Tfm при использовании отображаемой памяти приведены в табл. 4.
Таблица 4
Сравнение времени расчета

M
1 16 256 1024 4096

Tfs, мс
0,167 0,171 0,701 2,418 8,763

Tfa, мс
0,041 0,060 0,389 1,374 5,141

Tfm, мс
0,026 0,038 0,269 1,001 3,101

TCPU, мс
0,005 0,08 1,263 5,255 23,077

TCPU / Tfm
0,19 2,10 4,69 5,25 7,44

Из таблицы видно, что использование отображаемой памяти дает значительное преимущество перед другими методами обмена данными. В этом случае достигается наилучшее перекрытие процессов передачи данных и решения СЛАУ во времени. Полное время выполнения на ГП по сравнению с ЦП может сокращаться в 5—8 раз для большого числа СЛАУ (N > 1000).
Выводы. В работе проанализирована эффективность применения технологии CUDA, позволяющей производить массово-параллельные вычисления с использованием ГП компании NVIDIA, для ускорения решения СЛАУ с трехдиагональными матрицами в задачах расчета цен опционов как численного решения ДУЧП. При организации вычислений на ЦП СЛАУ решались методом прогонки, который является строго последовательным, поэтому при расчетах на видеокарте использовался метод ПЦР.
Результаты показали допустимую для практических задач погрешность математических вычислений вещественных чисел с плавающей точкой для метода ПЦР. Использование ГП не вносит дополнительных погрешностей.
Сравнение времени выполнения расчетов показало неэффективность применения ГП для решения одной СЛАУ ввиду малой загрузки устройства. Для увеличения доли занятости ГП авторы рекомендуют решать одновременно большое число СЛАУ, что актуально для современных систем алгоритмической торговли. Матрицы СЛАУ целесообразно формировать непосредственно в глобальной памяти видеокарты.
Наша реализация алгоритма ПЦР на видеокарте NVIDIA GTX 580 позволила достичь более чем 20-кратного ускорения вычислений по сравнению с однопоточным расчетом на ЦП

ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2012. Т. 55, № 10

Анализ эффективности применения технологии CUDA

25

Intel Core i5 3,8 ГГц при хранении данных в памяти видеокарты и одновременном решении более 1000 СЛАУ. Для задач, в которых формирование матриц СЛАУ непосредственно в памяти видеокарты невозможно, был проведен детальный анализ различных способов обмена данными между ОЗУ и памятью видеокарты. В этом случае для ГП удалось добиться ускорения в 5—8 раз при использовании механизма отображаемой памяти.

СПИСОК ЛИТЕРАТУРЫ

1. Simoes B. General-purpose computing on the GPU (GPGPU) [Электронный ресурс]: .

2. Zhang Y., Cohen J., Owens J. D. Fast Tridiagonal Solvers on the GPU // Proc. of the 15th ACM SIGPLAN Symp. on Principles and Practice of Parallel Programming (PPoPP 2010). 2010. January. P. 127—136.

3. Egloff D. High Performance Finite Difference PDE Solvers on GPUs. Technical report. QuantAlea GmbH. February. 2010.

4. Lyuu Y. D., Wen K. W., Wu Y. C. Performance of GPU for Pricing Financial Derivatives Convertible Bonds // To appear in J. of Information Science and Engineering [Электронный ресурс]: .

5. Gaikwad A., Muni Toke I. Parallel Iterative Linear Solvers on GPU: A Financial Engineering Case // Proc. of the 18th Intern. Conf. on Parallel, Distributed and Network-Based Computing. 2010. February. P. 607—614.

6. Wilmott P. Paul Wilmott on Quantitative Finance. Wiley, 2006. 1380 p.

7. Боресков А. В., Харламов А. А. Основы работы с технологией CUDA. М.: ДМК Пресс, 2010. 232 с.

8. Бахвалов Н. С., Жидков Н. П., Кобельков Г. М. Численные методы. М.: Бином. Лаборатория знаний, 2003. 632 с.

9. Sweet R. A parallel and vector variant of the cyclic reduction algorithm // SIAM J. Sci. Statits. Computing. 1988. Vol. 9, N 4. P. 761—765.

10. NVIDIA CUDA C Programming Guide. Version 4.2 2011 [Электронный ресурс]: .

Дмитрий Николаевич Шинкарук Юрий Александрович Шполянский Михаил Сергеевич Косяков

Сведения об авторах — аспирант; Санкт-Петербургский национальный исследовательский
университет информационных технологий, механики и оптики, кафедра вычислительной техники; E-mail: dimashink@gmail.com — д-р физ.-мат. наук, доцент; Санкт-Петербургский национальный исследовательский университет информационных технологий, механики и оптики, кафедра фотоники и оптоинформатики; Тбрикс АБ; ведущий математик; E-mail: shpolyan@mail.ru — канд. техн. наук; Санкт-Петербургский национальный исследовательский университет информационных технологий, механики и оптики, кафедра вычислительной техники; E-mail: mkosyakov@gmail.com

Рекомендована кафедрой вычислительной техники

Поступила в редакцию 08.02.12 г.

ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2012. Т. 55, № 10