РОЗВ''ЯЗУВАННЯ ЗАДАЧ З ПОЧАТКОВИМИ УМОВАМИ ДЛЯ СИСТЕМ ЗВИЧАЙНИХ ДИФЕРЕНЦіАЛЬНИХ РіВНЯНЬ НА БАГАТОЯДЕРНОМУ КОМП’ЮТЕРі З ГРАФіЧНИМИ ПРИСКОРЮВАЧАМИ ІНПАРКОМ Текст научной статьи по специальности «Математика»

УДК 519.6

М.Ф. ЯКОВЛЕВ*, Т.О. ГЕРАСИМОВА*, В.М. БРУСНІКІН**

РОЗВ’ЯЗУВАННЯ ЗАДАЧ З ПОЧАТКОВИМИ УМОВАМИ ДЛЯ СИСТЕМ ЗВИЧАЙНИХ ДИФЕРЕНЦІАЛЬНИХ РІВНЯНЬ НА БАГАТОЯДЕРНОМУ КОМП’ЮТЕРІ З ГРАФІЧНИМИ ПРИСКОРЮВАЧАМИ ІНПАРКОМ

Інститут кібернетики імені В.М. Глушкова НАН України, Київ, Україна Інститут проблем математичних машин і систем НАН України, Київ, Україна

Анотація. Розроблено і досліджено гібридний алгоритм методу Рунге-Кутта 4-го порядку розв’язування задач з початковими умовами для систем звичайних диференціальних рівнянь (СЗДР). Розглянуто питання програмної реалізації алгоритму на комп’ютерах з графічними процесорами. Наведені результати апробації алгоритму на багатоядерному комп ’ютері з графічними прискорювачами Інпарком, коефіцієнти прискорення та ефективності використання запропонованого алгоритму.

Ключові слова: паралельні обчислення, CUDA, гібридний алгоритм, СЗДР, метод Рунге-Кутта 4го порядку, ефективність розпаралелювання.

Аннотация. Разработан и исследован гибридный алгоритм метода Рунге-Кутта 4-го порядка для решения задач с начальными условиями для систем обыкновенных дифференциальных уравнений (СОДУ). Рассмотрены вопросы программной реализации алгоритма на компьютерах с графическими процессорами. Приведены результаты апробации алгоритма на многоядерном компьютере с графическими ускорителями Инпарком, коэффициенты ускорения и эффективности использования предложенного алгоритма.

Ключевые слова: параллельные вычисления, CUDA, гибридный алгоритм, СОДУ, метод Рунге-Кутта 4-го порядка, эффективность распараллеливания.

Abstract. A hybrid algorithm of the Runge-Kutta 4-th order method intended for the solving of initial-value problems in systems of ordinary differential equations (SODE) has been developed and investigated. The paper deals with problems related to program implementation of algorithm on multi-core computers with graphic accelerators. The results gained during program implementation and testing of algorithm on multi-core computer with graphic accelerators Inparcom are presented; acceleration and efficiency coefficients for the proposed algorithm are determined as well.

Keywords: parallel computing, CUDA, hybrid algorithm, SODE, Runge-Kutta 4-th order method, parallelization efficiency.

1. Вступ

При моделюванні на комп’ютері будь-яких реальних процесів і явищ за допомогою систем звичайних диференціальних рівнянь часто виникає необхідність розв’язувати задачі з початковими умовами (задачі Коші). СЗДР можуть існувати як самостійний клас задач або можуть виникати на проміжних стадіях розв’язування більш складних математичних задач. Такі задачі виникають, наприклад, при описуванні за допомогою СЗДР рухів, процесів або явищ, що змінюються в часі. Це стосується, зокрема, процесів хімічної кінетики, процесів, що відбуваються в ядерних реакторах, а також динамічних задач аналізу міцності конструкцій при застосуванні методу скінченних елементів по просторових змінних та ін. При розв’язуванні деяких задач, наприклад, пов’язаних із рухом керованих об’єктів, виникає необхідність розв’язувати СЗДР швидше, ніж відбувається процес у реальному часі; більш того, їх розв’язування потребує багатоваріантних розрахунків та значних обчислювальних ресурсів. Дуже часто виникає необхідність розв’язувати задачі з початковими умовами у випадку, коли вихідні дані (тобто формули для обчислення правих частин СЗДР та початкові умови) задані наближено. Такі задачі можна ефективно розв’язувати на висо-

20

© Яковлев М.Ф., Герасимова Т.О., Бруснікін В.М., 2015 ISSN 1028-9763. Математичні машини і системи, 2015, № 2

копродуктивних комп’ютерах з паралельною організацією обчислень, зокрема, на комп’ютерах MIMD-архітектури [1].

Вимоги до високопродуктивної обчислювальної техніки випереджають можливості традиційних паралельних комп’ютерів. У даний час такі комп’ютери можуть налічувати сотні і навіть тисячі процесорів (ядер). Але збільшення числа процесорів у паралельних комп’ ютерах часто призводить до значного збільшення комунікаційних втрат і зниження їх ефективності. Одним із основних напрямів підвищення продуктивності комп’ютерів є використання графічних процесорів (GPU). Цей напрям є особливо ефективним у випадку, коли необхідно виконувати великі обсяги однорідних арифметичних операцій. Використання графічних процесорів дало можливість розробити комп’ютери гібридної архітектури, яка об'єднує багатоядерні процесори MIMD-архітектури та графічні прискорювачі SIMD-архітектури.

У статті розглядається питання організації обчислень при розв’язуванні задач з початковими умовами для СЗДР високого порядку на комп’ютерах гібридної архітектури на прикладі застосування гібридного алгоритму методу Рунге-Кутта 4-го порядку, який використовує обчислювальні можливості CPU та GPU.

2. Постановка задачі

Нехай дано систему n звичайних диференціальних рівнянь. Задачу з початковими умовами для СЗДР n -го порядку на інтервалі [t0, Т] розглядатимемо у вигляді

<dy = f (X УХ (1)

at

y(t о) = У ^ (2)

де У = (Уі,У2,...,Уп)Т - шуканий вектор, а права частина системи - n -вимірна неперервна

вектор-функція f (t, У) = (f1 (t, У), f2 (t, У),..., fn (t, У ))Т .

При моделюванні реальних процесів на комп'ютері за допомогою СЗДР виникає ряд труднощів, зокрема, потреба мати справу з задачами з наближеними вихідними даними. Наближений характер вихідних даних може бути обумовлений такими причинами:

- похибками в вихідних даних;

- похибками у значеннях правої частини;

- застосуванням чисельного (дискретного) методу інтегрування і округленням чисел при обчисленнях;

- дискретизацією динамічних задач по просторових змінних.

Тому на практиці, як правило, замість задачі (1), (2) маємо задачу з наближеними вихідними даними:

d = f (^ vX (3)

at

v(t 0) = V0, (4)

де ||у0 - v0|| <5 , j(t, w) = f (t, w)+A(t, w), ||A(t, w|<A для довільних значень w(t) .

Такі задачі розв’язуються на комп'ютерах за допомогою чисельних методів.

ISSN 1028-9763. Математичні машини і системи, 2015, № 2

21

3. Паралельний алгоритм методу Рунге-Кутта 4-го порядку для комп’ютерів гібридної архітектури

Зупинимося на питанні організації обчислень при розв’язуванні задач з початковими умовами для СЗДР високого порядку на комп’ютерах гібридної архітектури (CPU+GPU) методом Рунге-Кутта 4-го порядку.

3.1. Декомпозиція даних

Для розв’язування СЗДР чисельними методами необхідно задавати таку вихідну інформацію: порядок системи звичайних диференціальних рівнянь; початкову та кінцеву точки інтервалу інтегрування; визначальну точність отриманого розв’язку; похибку завдання початкових умов; похибку завдання вектор-функції; повний вектор розв’язку (оскільки в загальному випадку всі компоненти вектор-функції правих частин системи залежать від усього вектора розв’язку), на початку інтегрування цей вектор містить початкові значення розв’язку; функцію обчислення значень вектор-функції.

Реалізація методів розв’язування задач даного класу на комп'ютерах гібридної архітектури вимагає автоматичного розподілу вихідних даних та компонент векторів правих частин у локальній пам'яті ядер кожного процесора. Кожний процес одночасно і незалежно реалізує власну програму обчислень з використанням функцій CUBLAS на GPU [2]. Комунікаційні операції та синхронізація роботи процесів здійснюються засобами MPI [3].

Щоб розв’язати СЗДР n -го порядку на MIMD-комп'ютері, її заздалегідь необхідно розбити на p блоків (p — кількість ядер). В паралельних алгоритмах розв’язування СЗДР на комп’ ютерах гібридної архітектури використовується одновимірний блочний розподіл даних та обчислень. Блоки вибираються приблизно одного розміру в залежності від порядку СЗДР n та кількості процесів (CPU) p , що використовуються для розв’язування задачі. Тоді розмір блока обчислюється за формулою

5

q

\n / p| + 1, q<n-p [ n / p [n / p|, q3n-p[n /p|

де q = 0,1,...,p — 1 - логічний номер CPU, значення Гаї дорівнює цілій частині числа а. За цією схемою розподіляються компоненти векторів у та f (t, y) . Такий розподіл дозволяє виконувати автоматичний розподіл обчислень компонент вектор-функції СЗДР на p процесів.

3.2. Розпаралелювання обчислень вектор-функції правих частин СЗДР

При розв’язуванні задач з початковими умовами на MIMD-комп'ютерах СЗДР необхідно задавати спеціальним чином, щоб забезпечити максимально можливе розпаралелювання операцій на кожному кроці інтегрування, зокрема, при обчисленні вектор-функції правих частин, та, якщо необхідно, наближення до матриці Якобі. При автоматичному розпарале-люванні обчислень на GPU слід враховувати, що лише в рідкісних випадках правила обчислення значень компонент вектор-функції СЗДР істотно відрізняються одна від одної. В цих випадках системи рівнянь, як правило не мають великого порядку і тому розв’язувати їх на комп'ютерах гібридної архітектури недоцільно. В більшості ж випадків побудова вектор-функції правих частин підпорядковується деякому спеціальному правилу і такі системи можуть мати достатньо великий порядок.

Як правило, при розв’язуванні задач з початковими умовами для СЗДР значна частина арифметичних операцій припадає на обчислення вектор-функції правих частин системи. Тому у більшості чисельних методів в першу чергу розпаралелюється обчислення компо-

22

ISSN 1028-9763. Математичні машини і системи, 2015, № 2

нент вектор-функції правих частин системи рівнянь на вибраній кількості CPU багатопроцесорного комп’ютера з графічними прискорювачами [4]. Останні s = p(q +1)-n процесорів оброблятимуть блоки по q рівнянь, а перші p — s процесорів - блоки по (q +1) -му рівнянню; тут p — кількість процесорів, q = 0, 1,...,p — 1 - логічний номер CPU.

Опишемо один із можливих способів завдання вектор-функції правих частин на мові програмування Сі. В кожному процесорі праві частини системи обчислюються за програмою

diffun (n, l, т, t, у, f)

{.......

for ( i = l; i < m; i++)

{ -/=•••; }

}

У кожному процесорі обчислюється вектор-функція від l -ої до (т — 1) -ої компоненти. При цьому l = kq, т = (k +1)q — 1, де k — логічний номер процесора, k = 0,1,..., p — 1. Величини l і m для кожного процесора повинні бути обчислені перед зверненням до програми обчислення вектор-функції правих частин. Цей шаблон використовується для всіх випадків завдання правих частин. Тут введені такі позначення: n - порядок системи звичайних диференціальних рівнянь, t - точка інтегрування, у - значення вектор-функції в точці t . Такий спосіб завдання вектор-функції має переваги:

• по-перше, полегшення при програмуванні СЗДР великого порядку, оскільки такі системи, як правило, мають регулярну структуру і тому фактично необхідно написати заголовок циклу і рівняння, залежне від параметра циклу;

• по-друге, при прихованому паралелізмі обчислення вектор-функції правих частин автоматично розподіляється між вибраною кількістю процесорів.

Якщо ж система звичайних диференціальних рівнянь має нерегулярну структуру, то кожне рівняння доведеться позначати власною міткою типу «case k:» (k: - ціле число) і забезпечувати звернення до функції за описаним вище правилом, використовуючи у програмі оператор типу «switch (i)».

Слід зазначити, що кількість арифметичних операцій, необхідних для обчислення компонент вектор-функції правих частин СЗДР, істотно впливає на ефективність і прискорення обчислень при розпаралелюванні на комп'ютері з графічними прискорювачами. Слід також зазначити, що задача може бути ефективно розв’язана, коли порядок системи n є достатньо великим числом, а також достатньо великим числом є і кількість операцій, необхідна для обчислення компонент вектор-функції f (u).

3.3. Гібридний алгоритм методу Рунге-Кутта 4-го порядку

Одним із методів, що застосовується для чисельного інтегрування задач з початковими умовами для СЗДР, є метод Рунге-Кутта 4-го порядку. Цей метод є найбільш розповсюдженим методом чисельного розв’язування задач з початковими умовами для СЗДР.

Класичний метод Рунге-Кутта 4-го порядку реалізується за формулами:

У (i+1) = у(i) + k + 2k 2 + 2k3 + k4)/6, (5)

де

k1 = hj (ti, у (0X k2 = hif (ti + hi/2, у(!) + 0,5k1X k3 = hif (ti + hi /2, у(i) + 0,5k2),

ISSN 1028-9763. Математичні машини і системи, 2015, № 2

23

k4 = h,f + hi, У(i) + 0,5*3 X (i = 0,1,2,...) .

На кожному кроці інтегрування цей метод потребує чотириразового обчислення вектор-функції f (t, y). Але для оцінки головного члена похибки, який визначає вибір кроку інтегрування, доводиться тричі застосовувати метод Рунге-Кутта.

Наведемо обчислювальну схему паралельного алгоритму методу Рунге-Кутта 4-го порядку для знаходження розв’язку систем звичайних диференціальних рівнянь на комп'ютерах гібридної архітектури (CPU+GPU) [5]:

1. У пам'ять кожного з p CPU посилається вихідна інформація.

Для подальшого викладу позначимо через m кількість компонент вектор-функції, що будуть обчислюватись одним CPU.

2. На інтервалі від початкового значення незалежної змінної до координати точки

і -1

виводу при кожному значенні ti = 10 + ^ hj (i= 1, 2,...) проводяться обчислення з викорис-

j=0

танням CPU та графічних процесорів у такому порядку:

a) у CPU обчислюється прогнозована довжина кроку інтегрування hi і посилається

в GPU;

b) далі виконується інтегрування системи, кожний і -ий крок якого обчислюється за схемою:

- у GPU обчислюються m компонент вектора k1, обчислені компоненти вектора k1 посилаються в CPU;

- у кожному CPU обчислюються m компонент вектора y(і) + 0,5k1, проводиться му-льтизбирання цього вектора та пересилка його в GPU;

- у GPU обчислюються m компонент вектора k2, обчислені компоненти вектора k2 посилаються в CPU,

у кожному CPU обчислюються m компонент вектора y(і) + 0,5k2, проводиться му-льтизбирання цього вектора та пересилка його в GPU;

- у GPU обчислюються m компонент вектора k3, обчислені компоненти вектора k3 посилаються в CPU;

- у CPU обчислюються m компонент вектора y(і) + k3, проводиться мультизбиран-ня цього вектора та пересилка його в GPU;

у GPU обчислюються m компонент вектора k4 , а потім m компонент вектора

y(і+1) = y(i) + v/6 де v = k1 + 2k2 + 2k3 + k4 ; обчислені компоненти вектора y(і+1) посилаються в CPU;

- кожен CPU проводить мультизбирання всього вектора y(і+1) та посилає його в

GPU.

Цим завершується і -ий крок інтегрування системи;

c) за схемою, описаною в п. «Ь», двократним інтегруванням з довжиною кроку hj 2 обчислюється вектор y^і+і) та однократним інтегруванням з довжиною кроку hi обчислюється вектор y ^г+і); ці вектори посилаються в GPU;

d) у GPU обчислюється похибка апроксимації системи y(l+l) = max y(l+l) - y(і+1м; ця

1<j<n І 1 1 І

похибка апроксимації посилається в CPU;

e) при виконанні умов досягнення заданої точності в CPU обчислюється уточнена довжина кроку інтегрування;

24

ISSN 1028-9763. Математичні машини і системи, 2015, № 2

f) в CPU перевіряється умова tt + ht > t p; якщо умова виконується, то як довжина

кроку інтегрування вибирається ht = tp — ti де tp — координата точки виводу розв’язку, в

іншому випадку довжина кроку інтегрування не коригується, вибрана довжина кроку посилається для продовження обчислень у GPU;

g) для отримання розв’язку у(г+:) в точці ti+1 = tt + ht повторюються пункти а) - е) з обчисленою довжиною кроку інтегрування.

Процес розв’язування задачі продовжується, доки не буде досягнуто кінцевої точки інтервалу інтегрування T. Після обчислення розв’язку в точці виводу запам'ятовуються вектор розв’язку, константа Ліпшиця та оцінка похибки розв’язку.

4. Програмна реалізація та чисельний експеримент. Експериментальне дослідження ефективності гібридного алгоритму

Розрахунки проводились на вузлі кластера Інпарком-G з такими характеристиками:

- процесори: 2 Xeon 5606 (4 ядра з частотою 2,13 ГГц);

- графічні прискорювачі: 2 Tesla M2090 (6 Гб пам’яті);

- об’єм оперативної пам’яті: 24 Гб;

- комунікаційне середовище: InfiniBand 40 Гбіт/с (з підтримкою GPUDirect), Gigabit Ethernet.

Також на вузлах встановлена бібліотека MKL10.2.6 та CUDA, починаючи з версії 3.2.

Розв’язувалась така тестова задача.

На інтервалі [0,0; 0,4] розв’язати СЗДР n -го порядку:

dui

dt

n—1

u j — Ui + n(l +1) + 2 +1

j=0

з початковими умовами ui (0) = 1 i = 0,1,2,..., n — 1.

В обчислювальному вузлі використовується чотириядерний процесор, на одне ядро приходиться одна карточка GPU. На кожному вузлі як обчислювальні процесори використовуються дві ігрові відеокарти.

Зупинимося на програмній реалізації алгоритму при роботі з GPU. Основними блоками операцій, що виконуються з GPU, є виділення пам’яті для змінних, копіювання даних на GPU, запуск обчислень, копіювання результатів в оперативну пам’ять, звільнення пам’ яті GPU. Зауважимо, що кожний процес одночасно і незалежно реалізує власну програму обчислень з використанням функцій CUBLAS на GPU для матрично-векторних операцій [2]. Комунікаційні операції та синхронізація роботи процесів здійснюються засобами MPI [3].

Для даної тестової задачі вектор-функція правої частини може бути легко обчислена на графічних прискорювачах, оскільки система має регулярну структуру і порядок цієї системи n є достатньо великим числом. Якщо ж права частина дуже складна, то в описаному алгоритмі треба ігнорувати обчислення вектор-функції в GPU і залишити її обчислення в CPU.

Нижче наведено таблицю, що містить часи розв’язування цієї тестової задачі при n =20000 на експериментальному зразку паралельного комп’ютера гібридної архітектури Інпарком. Результати експериментів наведені для гібридної архітектури 1 CPU, 1 GPU. У таблиці порівнюються часи розв’язування задачі при обчисленні компонент вектор-функції правої частини на СPU (перший рядок у таблиці) та на GPU (другий рядок у табл. 1). Зауважимо, що час розв’язування задачі на CPU становить 220 секунд.

ISSN 1028-9763. Математичні машини і системи, 2015, № 2

25

Таблиця 1. Час розв’язування тестової задачі на CPU і GPU

np=2 np=4 np=8 np=16

1 114,96 69,48 42,46 22,44

2 11,73 6,84 3,9 2,04

Нижче наведено графіки, з яких видно залежність коефіцієнта прискорення та коефіцієнта ефективності від кількості процесорів.

З графіка видно, що коефіцієнт прискорення практично лінійно збільшується зі збільшенням кількості процесів, що використовуються.

На наступному графіку зображено коефіцієнт ефективності використання кількості процесів.

5. Висновки

Дослідження паралельних алгоритмів розв’язування систем звичайних диференціальних рівнянь показали, що багатоядерні комп’ютери з графічними процесорами дають можливість суттєво прискорити час розв’язування задач, особливо у випадку, коли праві частини СЗДР обчислюються на GPU. У випадку, коли права частина дуже складна і має нерегулярну структуру, то обчислення вектор-функції доведеться залишити в CPU, що суттєво погіршить час розв’язування задачі.

Запропоновано алгоритм, який забезпечує високу ефективність розпаралелювання на GPU, оптимізує використання пам’яті CPU. Показано ефективність гібридного алгоритму методу Рунге-Кутта 4-го порядку розв’язування задач з початковими умовами для систем звичайних диференціальних рівнянь (СЗДР). У подальшому науковий інтерес представляють алгоритми для архітектури з n CPU, m GPU зі спільною або розподіленою пам’яттю.

26

ISSN 1028-9763. Математичні машини і системи, 2015, № 2

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

1. Молчанов И.Н. Интеллектуальные параллельные компьютеры на графических процессорах для решения научно-технических задач / И.Н. Молчанов // Праці Міжнародної молодіжної математичної школи «Питання оптимізації обчислень (ПОО - XXXVII)». - К.: Інститут кібернетики ім. В.М. Глушкова НАН України, 2011. - С. 121 - 122.

2. http://developer.download.nvidia.eom/compute/cuda/1 0/CUBLAS Library.

3. http: //www .mpi .org.

4. Яковлев М.Ф. Особливості розв’язування систем нелінійних та диференціальних рівнянь на паралельних комп’ютерах / М.Ф. Яковлев, Т.О. Герасимова, А.Н. Нестеренко // Праці Міжнародного симпозіуму «Питання оптимізації обчислень (ПОО - XXXV)». - Київ: Інститут кібернетики ім.

В.М. Глушкова НАН України, 2009. - Т. 2. - С. 435 - 439.

5. Параллельные алгоритмы решения задач вычислительной математики / Химич А.Н., Молчанов И.Н., Попов А.В., Чистякова Т.В. [и др.]. - Киев: Наукова думка, 2008. - 248 с.

Стаття надійшла до редакції 28.01.2015

ISSN 1028-9763. Математичні машини і системи, 2015, № 2

27