ИССЛЕДОВАНИЕ ЭФФЕКТИВНОСТИ РЕАЛИЗАЦИИ ВОЛНОВОГО УРАВНЕНИЯ НА ГРАФИЧЕСКОМ ПРОЦЕССОРЕ

АННОТАЦИЯ

Цель работы реализация на векторном процессоре разностного решения волнового уравнения. Разработана программная реализация разностного решения волнового уравнения в виде удаленного сервиса. Посредством вычислительных экспериментов было установлено ускорение до 5,5 раз при сравнении времени выполнения программы на центральном процессоре и на графическом процессоре.

Ключевые слова: волновое уравнение, OpenACC, алгоритм, векторный процессор.

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

Уравнение любой волны является решением дифференциального уравнения, которое называется волновым [1].

                                                                (1)

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

Рассмотрим в n-мерном пространстве волновое уравнение следующего вида:

                                        (2) 

Для решения задачи требуется определить начальные и краевые условия:

                                                                       (3)

Составим несложную, но эффективную разностную схему для численного решения задачи , . В качестве шаблона разностной схемы возьмем схему “крест”, вид которого для случая двух измерений показан на рисунке Рисунок 1. При произвольном числе измерений разностная схема “крест” может быть записана в виде:

                                              (4)

где:

;

;

.

Рисунок 1. Разностная схема крест для двумерного случая

Трёхслойная схема является явной. Вычисления с помощью такой схемы одинаково просты для любого числа измерений. Схема имеет порядок аппроксимации .

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

В итоге для множителя роста  получим квадратное уравнение:

                              (5)

При анализе корней можно выяснить, что разностной схема является устойчивой, если

                                                    (6)

Изучим схему “крест” на примере численного решения двумерной задачи , . Пусть . Тогда уравнение примет [2] вид:

                                         (7)

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

                       (8)

Нетрудно проверить, что задача , в области G имеет аналитическое решение  

Определим сетку по пространству:    где шаги по направлениям  и  определяются согласно формулам , . Запишем разностную схему уравнения :

                                               (9)

Для проведения вычислительных экспериментов запишем полученный алгоритм разностной схемы на языке Fortran с применением директив технологии OpenACC.

Технология OpenACC (Open Accelerators) является современным программным стандартом для параллельного программирования. При разработке его первой версии принимали участие такие компании, как NVidia, PGI, CAPS, Cray, а сейчас список уже гораздо больше. Главной целью было поставлено упрощение разработки параллельных программ. OpenACC описывает набор директив для компилятора, с помощью которых тот производит автоматическое распараллеливание программы нужным нам образом для вычислений на графическом процессоре [3].

  1   !$acc enter data copyin(E1(1:Nx,1:Ny), E2(1:Nx,1:Ny))

  2   do t=1,Nt,2        ! переход по временным слоям

  3   !$acc kernels            

  4      !$acc loop independent

  5      do i=2,Nx-1

  6         !$acc loop independent

  7         do j=2,Ny-1

  8            E1(i,j) = 2.0*E2(i,j)-E1(i,j) + ( E2(i-1,j) + E2(i+1,j) + E2(i,j-1) + E2(i,j+1) - 4.0*E2(i,j) ) * c2;

  9         end do  

10      end do

11      E1(Izl,2) = sin(c5*float(t));     ! "жесткое" излучающее условие

12      !$acc loop independent

13      do i=2,Nx-1

14         !$acc loop independent

15         do j=2,Ny-1

16            E2(i,j) = 2.0*E1(i,j)-E2(i,j) + ( E1(i-1,j) + E1(i+1,j) + E1(i,j-1) + E1(i,j+1) - 4.0*E1(i,j) ) * c2;

17         end do  

18      end do

19      E2(Izl,2) = sin(c5*float(t+1));    ! "жесткое" излучающее условие

20   !$acc end kernels

21   end do

22   !$acc exit data copyout(E1(1:Nx,1:Ny), E2(1:Nx,1:Ny)) delete(E1(1:Nx,1:Ny), E2(1:Nx,1:Ny))

Произведём серию вычислительных экспериментов для выявления ускорения работы алгоритма с использованием директив относительно однопоточного выполнения на центральном процессоре.

Цели экспериментов:

1. Нахождение максимального ускорения для фиксированного числа узлов сеточной области по времени на длину волны .
2. Нахождение максимального ускорения для фиксированного числа узлов сеточной области по пространству на длину волны .

Все вычисления были произведены на системе со следующими характеристиками:

• процессор Intel Core i5-3570K с тактовой частотой 3,4 ГГц;
• графический процессор NVIDIA GeForce GTX 770 – 1536 CUDA ядра;
• объём видеопамяти 4 Гб;
• объём оперативной памяти 8 Гб;
• операционная система Windows 10.

Рисунок 2. График зависимости ускорения от Q при Q_t = 1000

Рисунок 3. График зависимости ускорения от Q_t при Q = 50

Результаты первой серии вычислений представлены на рисунке 2. Здесь можно увидеть, что при малых значениях Q ускорение минимально, однако наблюдается его довольно быстрый рост и начиная с 40 сильно замедляется. В первом случае это объясняется тем, что вспомогательные расходы, такие как инициализация видеопроцессора и копирование сетки, не покрываются достаточным образом ускорением от параллелизма. Резкое падение ускорения при Q > 40 объясняется исчерпыванием памяти видеопроцессора. Вторая серия экспериментов представлена на рисунке 3. Ускорение работы начинается с 3 и плавно увеличивается. Это происходит из-за того, что общий размер пространственной сетки неизменен, но доля операций по вычислению на видеокарте к операциям по передаче данных постепенно увеличивается.

Список литературы:

1. Кузнецов С.И. Колебания и волны. Геометрическая и волновая оптика: учебное пособие. 2-е изд., перераб., дополн. - Томск: Изд-во Томского политехнического университета, 2007. - 170 с.
2. Воротникова, Д.Г. Разностное решение волнового уравнения на графических процессорах с повторным использованием попарных сумм дифференциального шаблона / Д.Г. Воротникова, Д.Л. Головашкин // Компьютерная оптика. – 2017. – Т. 41, № 1. – С. 134-138.
3. OpenACC 2.5 API Reference Card [Электронный ресурс] // Сайт программного стандарта OpenACC. – Электрон. дан. – [Б. м.], 2015. – URL: https://www.openacc.org/sites/default/files/inline-files/OpenACC_2.5_ref_guide_1.pdf (дата обращения: 10.06.2019).