Метод двухуровневого распараллеливания прогонки для решения систем алгебраических трехточечных уравнений на гибридных ЭВМ с многоядерными сопроцессорами

Abstract

В работе содержится описание метода двухуровневого распараллеливания прогонки (на общей памяти средствами OpenMP, на распределенной памяти средствами MPI) для решения систем линейных алгебраических трехточечных уравнений, возникающих при моделировании двумерных и трехмерных физических процессов. Здесь также описываются особенности реализации метода как на ЭВМ с универсальными процессорами, так и на гибридных ЭВМ с многоядерными сопроцессорами Intel Xeon Phi. Оценивается арифметическая сложность реализованного метода, и приводятся результаты численных экспериментов по исследованию масштабируемости метода. Ключевые слова: система линейных алгебраических трехточечных уравнений, трехдиагональные матрицы, метод прогонки, распараллеливание прогонки, параллельноконвейерный метод, метод Яненко, параллельные ЭВМ.

-

Российский федеральный ядерный центр - Всероссийский научно

исследовательский институт экспериментальной физики неизвестных, равным числу процессоров. Редуцированная система, состоящая из граничнопроцессорных точек, решается методом прогонки. Чтобы распараллелить решения этой системы, можно применить такие методы, как метод встречной прогонки и метод параллельноциклической редукции. Также известен параллельно-конвейерный метод [ 21, 22 ], предназначенный для решения множества трехдиагональных СЛАУ, возникающих при решении двумерных и трехмерных задач.

В работе [ 23 ] показано, как метод Яненко удалось применить при работе на ЭВМ с графическими ускорителями. В работе [ 24 ] было показано сравнение параллельно-конвейерного метода с фиксированным числом порций и метода Яненко, на втором этапе которого работал метод правой прогонки, и приведены результаты эффективности распараллеливания не более чем на 20 универсальных процессорах.

В данной работе приводятся описания метода Яненко, параллельно-конвейерного метода и комбинированного метода, представляющего собой метод Яненко, на втором этапе которого работает параллельно-конвейерный метод встречной прогонки с автоматическим подбором числа порций. Также приводятся оценки арифметической сложности методов, особенности их реализации для работы на ЭВМ с многоядерными сопроцессорами, численные исследования эффективности распараллеливания комбинированного метода и сравнение его с параллельноконвейерным методом. 2. Параллельно-конвейерный метод

Основная идея метода состоит в том, что процессор, выполнив прямой или обратный ход прогонки по части точек, находящихся на нем, передает данные соседнему процессору, чтобы тот мог продолжить прогонку, а сам в это время начинает обрабатывать следующую порцию точек. Если исключить время разгона и торможения конвейера, то все остальное время в счете задействованы все процессоры. Поясним вышесказанное с помощью рисунка 1. На рисунке конвейер состоит из 3 процессоров, строки на рисунке (сверху вниз) соответствуют последовательным моментам времени, стрелки вправо – протягиванию прямого хода прогонки в одном направлении, а влево – в другом. Светло-серым цветом отмечены процессоры, протянувшие прогонку в одном направлении, темно-серым – в обоих.

Рис. 1. Схема параллельно-конвейерного алгоритм Следует отметить, что число порций, определяющее число тактов конвейера, подбирается автоматически во время решения задачи на ЭВМ, исходя из времени счета на предыдущем временном шаге задачи. В случае, когда число порций становится равно единице, ПКМ вырождается в метод правой прогонки, когда в один момент времени работает только один процессор. Подробнее этот метод описан в [ 22 ]. 3. Метод Яненко

В отличие от параллельно-конвейерного метода Яненко можно применить не только для решения множества систем трехточечных уравнений, но и для решения одной системы. Рассмотрим систему линейных уравнений с трехдиагональной матрицей следующего вида: ai yi1  bi yi  ci yi1  di ,bi  0,i  1, 2, , N 1, b0 y0  c0 y1  d0 , aN yN1  bN yN  dN . (1) Для простоты пусть на каждом процессоре будет одинаковое количество точек m=K/M, где K – число неизвестных (в нашем случае K = N+1), M – число процессоров, индексация будет глобальная. Таким образом, на отдельном процессоре с номером j будет находиться лишь часть уравнений системы (1) с номерами от (j-1)*m+1 до j*m, где j – номер процессора. Обозначим y j*m через z j , j=0,..,M и будем искать решения системы (1) в виде

y(j1)*mi  z j1ui  z jvi  wi , i  1, ,m 1, j  1, , M , где u, v, w есть решения следующих систем уравнений: (2) (3) (4) aiui1  biui  ciui1  0, u j1*m  1,

u j*m  0; aivi1  bivi  civi1  0, v j1*m  0, v j*m  1; aiwi1  biwi  ciwi1  di , w j1*m  0, wj*m  0; i   j 1* m+1, ,j * m 1, j  1, , M .

Решения этих трех систем можно найти методом прогонки, причем независимо на каждом процессоре. Будем называть их предрешениями, а этот этап решения задачи – этапом нахождения предрешений.

В уравнения с номерами j*m из системы (1): a j*m y j*m1  bj*m y j*m  c j*m y j*m1  d j*m , j  0, , M , подставляем вместо «y» комбинации (2). Таким образом, получаем систему трехточечных уравнений для нахождения z j , имеющую следующий вид:

Aj z j1  Bj z j  C j z j1  Dj , j  1, , M 1,

B0 z0  C0 z1  D0 , AM zM 1  BM zM  DM с коэффициентами:

B0  b0  c0u1, C0  c0v1, D0  d0  c0w1, Aj  a j*mu j*m1, Bj  a j*mvj*m1  bj*m  cj*mu j*m1, Cj  cj*mvj*m1, Dj  d j*m  aj*mwj*m1  cj*mwj*m1, j  1, , M1,

AM  aM*muM*m1, BM  bM*m  aM*mvM*m1, DM  dM*m  aM*mwM*m1.

Назовем этот этап – этапом нахождения гранично-процессорных решений. Размерность этой системы уравнений равна количеству процессоров, что существенно меньше количества точек задачи.

На последнем этапе восстанавливаем окончательное решение по формуле (2). Итак, метод Яненко содержит три этапа: 1) нахождение предрешений, 2) нахождения гранично-процессорных решений, 3) восстановление решения.

Первый и третий этапы выполняются независимо на каждом устройстве (универсальном процессоре, ускорителе или сопроцессоре). А второй этап требует коммуникаций между MPIпроцессами. Таким образом, эффективность распараллеливания метода Яненко напрямую зависит от эффективности распараллеливания этого этапа.

Возникает вопрос, каким методом решать систему (4), состоящую из граничнопроцессорных точек? В качестве возможных методов рассматривались метод правой прогонки, метод встречной прогонки [ 7 ] и метод параллельно-циклической редукции [6]. Рассмотрим эти методы подробнее. Формулы для метода правой прогонки выглядят следующим образом:  0  c0 / b0 ,

0  d0 / b0 , i  di  aii1  / bi  ai i1 , i  1,..., M 1,

 i  ci / bi  ai i1 , yM  dM  aM  M 1  / bM  aM M 1 , yi  i yi1  i , i  M 1,...,0. Для метода встречной прогонки формулы выглядят так: 1,0  c0 / b0 ,

1,0  d0 / b0 , 1,i  di  ai1,i1  / bi  ai1,i1 , i  1,..., k,  2,M  aM / bM ,  2,M  dM / bM ,

 2,i  ai / bi  ci 2,i1 ,  2,i  di  ci 2,i1  / bi  ci 2,i1 , i  M 1,..., k 1,

1,i  ci / bi  ai1,i1 , yk  1,k 1,k 2,k1  / 11,k 2,k1 , yi  2,i yi1   2,i , i  k 1,..., M .

yi 1,i yi1  1,i , i  k 1,...,0, (5) (6) 4. Вычисление арифметической сложности методов

Оценим арифметическую сложность описанных выше методов при решении одной системы линейных уравнений.

Поскольку параллельно-конвейерный метод применяется при решении множества систем методом правой прогонки по формулам (5), то число арифметических операций в алгоритме при решении одной системы совпадает с числом арифметических операций в методе правой прогонки, а именно равно 8K-7, где K – число неизвестных в системе (1) [ 25 ].

Для того чтобы оценить, насколько хорошо метод Яненко будет работать в параллельном режиме, вычислим количество операций, которые нужно сделать при нахождении решений системы (1) методом Яненко, а при нахождении решений системы (4) на втором этапе метода Яненко – методом правой прогонки.

Число арифметических операций в методе складывается как сумма операций на всех этапах метода, описанных ранее. Искомое число есть S  s1  s2,1  s2,2  s3 , где s1 - число операций на этапе нахождения предрешений, s2,1 - число операций в расчете коэффициентов матрицы системы для нахождения гранично-процессорных решений, s2,2 - число операций в методе прогонки на этом этапе, s3 - число операций при восстановлении решения исходной системы.

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

Сначала заметим, что если решать системы (3) методом прогонки, не учитывая особенностей этих систем, то мы получим число арифметических операций s1  3* (8K 7)  24 K 21.

Покажем, как можно упростить нахождение предрешений, которые будем искать методом правой прогонки по формулам (5). Обозначим прогоночные коэффициенты первой системы (относительно u) через 1,i , 1,i , второй системы (относительно v) через  2,i , 2,i , третьей системы (относительно w) через  3,i ,  3,i .

Заметим, что все коэффициенты 1,i , 2,i , 3,i равны между собой, так как Также заметим, что  2,i  0 . Покажем это методом индукции. Сначала вычислим коэффичим 1,i  2,i  3,i  i . циент на левой границе задачи: 2, j1*m1  d j1*m1  ai2, j1*m  b j1*m1  ai j1*m

0  aiv j1*m b j1*m1  ai j1*m 

ai  0 b j1*m1  ai j1*m  0 . Пусть теперь  2,i  0 , тогда найдем значение  2,i1 : 2,i1  di1  ai12,i  0  ai  0  0 , i   j 1*m+1, ,j *m  2 .

bi1  ai1i bi1  ai1i Используя полученные выше утверждения, формулы прогоночных коэффициентов для нахождения предрешений будут выглядеть следующим образом: i   ci , 1,i  di  ai1,i1  ai1,i1 ,  2,i  0, 3,i  di  ai3,i1 .

bi  aii1 bi  aii1 bi  aii1 bi  aii1 При этом формулы нахождения предрешений будут выглядеть следующим образом ui iui1  1,i , vi ivi1, wi iwi1  3,i , i   j 1 * m+1, ,j * m 1 , j  0, , M . Таким образом, для прямого хода прогонки необходимо число операций равное s1,1  (2 1 2  0  3)  (m 1)  8(m 1), а для обратного хода прогонки – s1,2  (2 1 2)  (m 1)  5(m 1).

Для всего этапа требуется s1  s1,1  s1,2  8(m 1)  5(m 1) 13(m 1), что гораздо меньше, чем s1  24 * (m 1)  21  24m  3 . То есть, учитывая структуру систем уравнений, нам удалось уменьшить число арифметических операций почти в 2 раза, по сравнению с трехкратным применением метода правой прогонки для решения трех систем уравнений.

В итоге, число операций на всех этапах в методе Яненко (при использовании правой прогонки на втором этапе) выглядит следующим образом (Из (4) число неизвестных равно M+1): s1  13M (m 1), s2,1  5 10(M 1)  5  10M , s2,2  8(M 1)  7  8M 1, s3  4M (m 1), S  13M (m 1) 18M 1  4M (m 1)  17M (m 1) 18M 1. В случае, когда число процессоров равно 1 (M=1, m=K), то S 17(K 1) 18 1 17K  2 для системы (1). Видно, что при работе на одном процессоре метод Яненко будет проигрывать методу правой прогонки, в котором S  8K  7 .

Теперь попробуем оценить арифметическую сложность алгоритмов, которые мы планируем применить на втором этапе методе Яненко (либо встречная прогонка, либо параллельноциклическая редукция). Получаемые нами в дальнейшем оценки будут составлять слагаемое s2,2 в вышеприведенной формуле арифметической сложности метода Яненко. Из (4) число уравнений в системе на этом этапе равно M + 1.

Из (6) видно, что число арифметических операций для метода встречной прогонки равно 2  6k  2  6(M  k 1)  5  2(k)  2(M  k)  8M  3 , что немногим больше числа операций в методе правой прогонки (8M + 1).

Теперь рассмотрим метод параллельно-циклической редукции. Из (7) видно, что на этапе преобразования коэффициентов матрицы для гранично-процессорных точек число операций будет выражаться величиной S2,1  6  (2 1 4 1 4)(M 1)  6  6 12M 12  6 12M . После таких преобразований мы получаем две независимые системы трехточечных уравнений, которые можно решать одновременно на разных процессорах (с четными и нечетными номерами, соответственно). Формулы преобразований будут теми же, но изменится индексация. Для системы, в которой на четных местах стоят ненулевые коэффициенты диапазон изменения индекса следующий: i  0, 2,..., 2 *M / 2, а для системы, в которой на нечетных местах стоят ненулевые коэффициенты – i 1,3,...,2*(M 1) / 2 1 . В этом случае число арифметических операций будет равно

S2,2 12M / 2 12(M 1) / 2 12M / 2  M / 2  (M 1) / 2  (M 1) / 2   12M / 2  (M 1) / 2  M / 2  (M 1) / 2 12M  0.5  0.5 12M 1, где M  – целая часть числа M, а M – дробная часть числа M.

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

S2,3  12M / 4 12(M 1) / 4 12(M  2) / 4 12(M  3) / 4   12 M / 4  (M 1) / 4  (M  2) / 4  (M  3) / 4  12M / 4  (M 1) / 2  (M  2) / 4  (M  3) / 4   12 M  6 / 4 123 / 4  2 / 4 1 / 4  0  12 M  3.

Этот процесс продолжается до тех пор, пока в каждой из получившихся систем не останется по одному уравнению с одним неизвестным, из которого мы и найдем решение нашей системы. Число таких итераций будет H  log2  M 1  d , где d=0, если M является степенью двойки и d=1, если M не является степенью двойки. Тогда общее число арифметических операций в этом методе будет равно

S2  S2,1  S2,2  ...  S2,H  12(M  M 1 M  3  ...  M  (2H  3))   12HM 12(1 3  ...  (2H  3))  12HM 12(1 2H  3) / 2 * (H 1)  ли M является степенью двойки, d=1, если M не является степенью двойки. Очевидно, что это гораздо больше, чем в методе встречной прогонки, но здесь практически всегда все процессоры работают независимо (за исключением тех случаев, когда на процессоре уже найдено решение системы уравнений).

Сведем полученные результаты для методов распараллеливания прогонки в один список: 1. Метод правой прогонки: 1.1. Число арифметических операций: 8M+1; 1.2. Число независимо работающих MPI-процессов: 1. 2. Метод встречной прогонки: 2.1. Число арифметических операций: 8M+3; 2.2. Число независимо работающих MPI-процессов: 2 (почти всегда). 3. Метод параллельно-циклической редукции: 3.1. Число арифметических операций: 2 12M log2 M 1 12log2 M 1  d 1 12dM , d=0, если M является степенью двойки, d=1, если M не является степенью двойки;

3.2. Число независимо работающих MPI-процессов: M-1 (почти всегда).

При трехмерной процессорной декомпозиции M – число процессоров вдоль направления прогонки.

Из полученных оценок следует, что неизвестно, какой метод прогонки будет эффективнее работать на втором этапе метода Яненко: метод встречной прогонки (число операций меньше, но и число независимо работающих MPI-процессов мало) либо метод параллельноциклической редукции (число операций существенно больше, но и число независимо работающих MPI-процессов почти максимально). При работе на графических ускорителях метод параллельно-циклической редукции на втором этапе метода Яненко оказался лучше метода встречной прогонки примерно на 20%. 5.Реализация двухуровневого распараллеливания метода прогонки Реализация на общей памяти была выполнена с использованием стандарта OpenMP, а реализация на распределенной памяти – с использованием стандарта MPI.

Как уже было отмечено, при решении двумерных или трехмерных задач, использующих метод расщепления, возникает множество независимых систем трехточечных уравнений. При решении задач на параллельных ЭВМ как с общей, так и распределенной памятью на каждом MPI-процессе лежит часть уравнений таких систем, и, чтобы найти их решения, мы применяем методом Яненко, на втором этапе которого в силу простоты реализации был применен метод встречной прогонки. При этом каждый OpenMP-поток выполняет вычисления для нескольких систем уравнений. Для повышения эффективности распараллеливания при решении множества систем трехточечных уравнений решено было использовать параллельно-конвейерный метод на втором этапе метода Яненко, то есть множество систем (4) разбивается на порции, в каждой из которых OpenMP-потоки выполняли арифметическую работу для нескольких систем. Следует отметить, что здесь, как и в [ 22 ] был реализован автоматический подбор числа порций конвейера в зависимости от времени работы. А именно:

1. На первом шаге задаем начальное число порций точек, которые необходимо обработать за один такт конвейера, равным 1, то есть все линии точек на процессоре будут обработаны за 1 такт конвейера.

2. На каждом последующем шаге число порций увеличивается на 1.

3. Если на текущем шаге время работы метода прогонок стало больше, чем на предыдущем, то число порций уменьшается на 1.

4. На каждом последующем шаге это число будет уменьшаться до тех пор, пока время прогонок на текущем шаге снова не станет больше, чем на предыдущем. Тогда это число начнет увеличиваться. 6. Численные исследования эффективности распараллеливания В данном разделе приводятся результаты численных экспериментов для определения следующих характеристик рассмотренных методов:

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

Тестирование проведено на трехмерной задаче газовой динамики с неявной схемой аппроксимации дифференциальных уравнений на кубической трехмерной сетке, в каждом направлении было взято 200 ячеек. Таким образом, общее число точек задачи равно 200х200х200=8 000 000 ячеек.

Задача была просчитана в различных сочетаниях количества OpenMP-потоков и MPIпроцессов. Было получено, что физические результаты расчета не зависят от соотношения числа OpenMP-потоков и MPI-процессов.

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

Тестирование проводилось на вычислительном комплексе Центра компетенций и обучения (ЦКО) [ 26 ], который занимает 39-ю строчку списка «50 наиболее мощных компьютеров СНГ» (в редакции от 28.09.2915) [ 27 ]. В каждом вычислительном узле находилось 2 CPU Intel Xeon E5-2680 (10 вычислительных ядер в каждом) и 2 Intel Xeon Phi coprocessor 7120P (61 ядро в каждом). Распараллеливание внутри устройств выполнено на общей памяти с использованием OpenMP, а между устройствами – с помощью MPI (при использовании библиотеки S-MPI [ 28 ]).

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

Масштабирование проводилось методом деления – при увеличении числа MPI-процессов в 2 раза число точек задачи не менялось (сильная масштабируемость).

На 64-х универсальных узлах (процессорная декомпозиция задачи – 8х4х4) второй этап с методом встречной прогонки выполнился в 1,28 раза быстрее, чем с методом правой прогонки. А параллельно-конвейерный метод встречной прогонки выполнился в 1,46 раза быстрее, чем метод встречной прогонки без конвейерности.

Сравним масштабируемость метода Яненко с параллельно-конвейерным методом встречной прогонки на втором этапе и параллельно-конвейерного метода на универсальной архитектуре. Введем следующие обозначения: - ПКМ – время работы параллельно-конвейерного метода [сек]; - Яненко – время работы метода Яненко [сек];

V p  t1 / ( ptp ) 100% – эффективность распараллеливания при сильной масштабируемости [%], где - t1 – время работы программы на одном вычислительном узле; - t p – время работы программы на p вычислительных узлах.

На рисунке 2 приведена эффективность распараллеливания методов прогонки с помощью различных методов в зависимости от числа универсальных узлов.

1 2 4

8 Число узлов 16 32 64 Рис. 2. Эффективность распараллеливания методов прогонки в зависимости от числа универсальных узлов (сильная масштабируемость) 100,00 90,00 80,00 ] 70,00 [%60,00 p V 50,00 40,00 30,00 Теперь будем масштабировать задачу методом умножения – при увеличении числа MPIпроцессов в 2 раза число точек задачи также увеличивается в 2 раза (слабая масштабируемость).

На рисунке 3 приведена эффективность распараллеливания методов прогонки с помощью различных методов в зависимости от числа универсальных узлов.

Введем следующие обозначения:

V p  t1 / tp  100% – эффективность распараллеливания при слабой масштабируемости [%], где - t1 – время работы программы на одном вычислительном узле; - t p – время работы программы на p вычислительных узлах.

1 2 4 Числ8о узлов 16 32 64 Рис. 3. Эффективность распараллеливания методов прогонки в зависимости от числа универсальных узлов (слабая масштабируемость) На рисунках 2 и 3 можно заметить, что на универсальной архитектуре метод Яненко работает быстрее, чем ПКМ (более чем в 1,8 раза на 64-х узлах), то есть метод Яненко на данной задаче показал себя лучше, чем параллельно-конвейерный метод. При этом эффективность распараллеливания метода Яненко составила не ниже 70% на 64-х универсальных узлах.

Теперь посмотрим, как ведет себя реализованный метод Яненко на ЭВМ с использованием сопроцессоров. В таблице 1 приведена эффективность распараллеливания между сопроцессорами в native-режиме (слабая масштабируемость). Из этой таблицы следует, что эффективность распараллеливания на 4-х сопроцессорах (92,83) выше, чем на 4-х универсальных узлах на рисунке 2 (ниже 70%).

Таблица 1. Эффективность распараллеливания на сопроцессорах в native-режиме Число сопроцессоров

Vp, [%] 1 2

4 В таблице 2 приведены значения эффективности распараллеливания (Vp) и ускорения (слабая масштабируемость) от использования сопроцессоров в симметричном режиме для реализованного метода.

Таблица 2. Эффективность распараллеливания и ускорение от использования сопроцессоров в симметричном режиме Число гибридных узлов

Vp, [%] Ускорение от использования сопроцессоров 1

2 7. Заключение

Выполнена реализация метода распараллеливания прогонки для решения множества систем трехточечных уравнений, возникающих при моделировании двумерных и трехмерных физических процессов, на гибридных ЭВМ с сопроцессорами Intel Xeon Phi, имеющими как общую, так и распределенную память. Метод заключается в том, что для распараллеливания на распределенной памяти применяется метод Яненко с встречной прогонкой на втором этапе. Для распараллеливания множества встречных прогонок реализован параллельно-конвейерный метод с автоматическим подбором количества порций. При распараллеливании на общей памяти каждый OpenMP-поток решает несколько систем трехточечных уравнений. На тестовой задаче показана независимость решения от числа MPI-процессов и OpenMP-потоков. По эффективности распараллеливания метод Яненко выигрывает у параллельно-конвейерного метода на универсальной архитектуре, при этом эффективность распараллеливания на универсальной ЭВМ составила около 70% на 64-х узлах. Получено ускорение счета в 2.16 раза от использования четырех сопроцессоров на двух вычислительных узлах.

В дальнейшем планируется попробовать применить метод параллельно-циклической редукции на втором этапе метода Яненко, чтобы исследовать возможность повышения эффективности распараллеливания, и использовать векторизацию, чтобы повысить ускорение от использования сопроцессоров. Литература 2. Jeffers J., Reinders J. Iтеуд Xeon Phi Coprocessor High-Performance Programming. Morgan

Kaufmann Publishers Inc. San Francisco. 2013. 409 p. 10. Goddeke D., Strzodka R. Cyclic reduction tridiagonal solvers on GPUs applied to mixed-precision multigrid // IEEE Transactions on Parallel and Distributed Systems. 2011. Vol. 22. P. 22–32. 13. Sakharnykh N. Efficient tridiagonal solvers for ADI methods and fluid simulation. NVIDIA GPU

Technology Conference. September 2010. 15. Stone H. S. An efficient parallel algorithm for the solution of a tridiagonal linear system of equations // Journal of the ACM. 1973. Vol. 20. P. 27–38. 16. Chang L.-W., Stratton J. A., Kim H.-S., Hwu W.-M. A scalable, numerically stable, highperformance tridiagonal solver using GPUs. High Performance Computing, Networking, Storage and Analysis (SC), 2012 International Conference for. IEEE Computer Society Press, 2012. 11 p. 17. Polizzi E., Sameh A. A parallel hybrid banded system solver: The SPIKE algorithm // Parallel

Computing. 2006. Vol. 32, No. 2. P. 177–194. 18. Polizzi E., Sameh A. SPIKE: A parallel environment for solving banded linear systems // Computers and Fluids. 2007. Vol. 36, No. 1. P. 113– 120. 21. Povitsky A. Parallelization of the pipelined Thomas algorithm. ICASE Report No. 98-48. Hampton. NASA Langley Research Center. 1998. 22 p. 22. Сапронов И.С., Быков А.Н. Параллельно-конвейерный алгоритм // Атом. 2009. № 44. С. 2425. 23. Быков А. Н., Ерофеев А. М., Сизов Е. А., Федоров А. А. Метод распараллеливания прогонки на гибридных ЭВМ // Вычислительные методы и программирование. 2013. Т. 14. Раздел 2. С. 43-47. The method of two-level parallelization for solving tridiagonal systems on hybrid computers with multicore coprocessors А.А. Fedorov, А.N. Bykov Russian federal nuclear center All-russian research institute of experimental physics This paper describes a method of two-level parallelization for solving tridiagonal systems (thread-level parallelism using OpenMP and process-level parallelism using MPI) arising in modeling of two-dimensional and three-dimensional physical processes. It also describes the features of implementation for parallel multiprocessor systems and for hybrid multiprocessor systems with multicore coprocessors. The authors also evaluate arithmetic complexity of this method and show the results of numerical experiments to investigate the scalability of the method. 2. Jeffers J., Reinders J. Iтеуд Xeon Phi Coprocessor High-Performance Programming. Morgan

Kaufmann Publishers Inc. San Francisco. 2013. 409 p.

Davidson A., Zhang Y., D. Owens J. An auto-tuned method for solving large tridiagonal systems on the GPU. URL: http:// idav.ucdavis.edu/func/return_pdf?pub_id=1052 (accessed: 07.11.2015). Kim H.-S., Wu S., Chang L.-W., Hwu W.-M. A scalable tridiagonal solver for gpus. In Parallel Processing (ICPP), 2011 International Conference on. P. 444 –453. 10. Goddeke D., Strzodka R. Cyclic reduction tridiagonal solvers on GPUs applied to mixed-precision multigrid // IEEE Transactions on Parallel and Distributed Systems. 2011. Vol. 22. P. 22–32. 13. Sakharnykh N. Efficient tridiagonal solvers for ADI methods and fluid simulation. NVIDIA GPU

Technology Conference. September 2010.

References

1. Bykov A.N. , Veselov V.A. , Voronin B.L. , Erofeev A.M. Metodika RAMZES-CP dlya rascheta prostranstvennykh dvizheniy mnogokomponentnykh teploprovodnykh sred v eylerolagranzhevykh koordinatakh [RAMZES-KP Technique for Computation of Multicomponent HeatConducting Space Motions in Euler-Lagrange Coordinates] // RFNC-VNIIEF Proceedings . 2008 . Vol 13. P. 50 - 57 .

NVIDIA. TESLA GPU ACCELERATORS FOR SERVERS . URL: http:// http://www.nvidia.com/object/tesla-servers. html (accessed: 07.11 . 2015 ).

Hilger. Bristol. 1986 . P. 274 - 280 .

7. Samarskiy A. A. Vvedenie v chislennye metody [Introduction to Numerical Methods] . Moscow. 1982 . P. 34 - 38

NVIDIA Corporation. CUDA CUSPARSE Library . URL: https://developer.nvidia.com/cusparse (accessed: 07 . 11 . 2015 ).

9. Sengupta S. , Harris M. , Zhang Y. , Owens J. D. Scan primitives for GPU computing . URL: http:// idav.ucdavis.edu/func/return_pdf?pub_id= 915 ( accessed : 07 . 11 . 2015 ).

11. Davidson A. , Owens J. D. Register packing for cyclic reduction: A case study . URL: http:// idav.ucdavis.edu/func/return_pdf?pub_id= 1052 ( accessed : 07 . 11 . 2015 ).

12. Egloff D. High performance finite difference PDE solvers on GPUs . URL: http:// www.researchgate.net/publication/228942149_High_ performance_finite_difference_PDE_solvers _on_GPUs (accessed : 07 . 11 . 2015 ).

14. Zhang Y., Cohen J. , et al. Fast tridiagonal solvers on the GPU . URL: http:// idav.ucdavis.edu/func/return_pdf?pub_id= 978 ( accessed : 07 . 11 . 2015 ).

15. Stone H. S. An efficient parallel algorithm for the solution of a tridiagonal linear system of equations // Journal of the ACM . 1973 . Vol. 20 . P. 27 - 38 .

16. Chang L .-W., Stratton J. A. , Kim H .-S., Hwu W.-M. A scalable, numerically stable, highperformance tridiagonal solver using GPUs. High Performance Computing, Networking, Storage and Analysis (SC), 2012 International Conference for . IEEE Computer Society Press, 2012 . 11 p.

17. Polizzi E. , Sameh A . A parallel hybrid banded system solver: The SPIKE algorithm // Parallel Computing . 2006 . Vol. 32 , No. 2. P. 177 - 194 .

18. Polizzi E. , Sameh A . SPIKE: A parallel environment for solving banded linear systems // Computers and Fluids. 2007 . Vol. 36 , No. 1. P. 113 - 120 .

19. Venetis I. , et al. A general tridiagonal solver for coprocessors: Adapting g-SPIKE for the Intel Xeon Phi . URL: http://www.researchgate.net/publication/282286515 (accessed: 07 . 11 . 2015 ).

20. Yanenko N.N. , Konovalov A.N. , Bugrov A.N. , Shustov G.V. Ob organizatsii parallel'nykh vychisleni i “rasparallelivanii” progonki [On organization of parallel computations and parallelization of the tridiagonal matrix algorithm] // Numerical Methods of Continuum Mechanics. 1978 . Vol. 9 . №7. P. 139 - 146 .

21. Povitsky A. Parallelization of the pipelined Thomas algorithm . ICASE Report No. 98-48. Hampton. NASA Langley Research Center . 1998 . 22 p.

22. Sapronov I.S. , Bykov A.N. Parallel'no-konveyernyy algoritm [Pipelined Thomas algorithm]// Atom. 2009 . № 44. P. 24 - 25 .

23. Bykov A. N. , Erofeev A. M. , Sizov E. A. , Fedorov A. A. Metod rasparallelivaniya progonki na gibridnykh EVM [A parallel sweep method for hybrid supercomputers] // Numerical methods and programming . 2013 . Vol. 14 . Section 2. P. 43 - 47 .

24. Il 'in S.A. Starchenko A.V. Rasparallelivanie skhemy pokomponentnogo rasshchepleniya dlya chislennogo resheniya uravneniya teploprovodnosti [Parallelization of scheme component splitting for the numerical solution of the heat equation]. Parallelnye vychislitelnye tekhnologii (PaVT' 2015 ): Trudy mezhdunarodnoj nauchnoj konferentsii (Ekaterinburg, 31 marta - 2 aprelya 2015) [Parallel Computational Technologies (PCT' 2015 ): Proceedings of the International Scientific Conference (Ekaterinburg, Russia, March , 31 - April, 2, 2015 )]. Chelyabinsk, Publishing of the South Ural State University, 2015 . P. 399 - 402 .

25. Verzhbitskiy V.M. Vychislitel 'naya lineynaya algebra [Numerical Linear Algebra] . Мoscow . 2013 . P. 296 - 299 .

26. Tsentr kompetentsiy i obucheniya. Predostavlenie vychislitel'nykh resursov TsKO [The center of competence and training. Providing computing resources CCT] . URL: http://compcenter.org/predostavlenie-vychislitelnyh-resur (accessed: 07.11 . 2015 ).

27. Top- 50 . 23 -ya redaktsiya ot 28.09.2015 [Top-50. 23th edition of 28.09 . 2015 ]. URL: http://top50.supercomputers.ru/?page= rating (accessed: 07.11 . 2015 ).

28. Voronov G.I. , Trushchin V.D. , Shumilin V.V. , Ezhov D .V. Programmnyy kompleks S-MPI dlya obespecheniya razrabotki, optimizatsii i vypolneniya vysokoparallel'nykh prilozheniy na superkomp'yuternykh klasternykh sistemakh [S-MPI software system for the development, optimization and execution of parallel applications on supercomputer clusters] // Voprosy atomnoy nauki i tekhniki . Mathematical modeling of physical processes . 2013 . Vol. 3 . P. 55 - 60 .