Управление потоками данных в параллельных алгоритмах вычислительной линейной алгебры - (диплом)
p>Тогда вполне естественно приписать каждый блок отдельному процессору. При этом, если размерности блоков сильно отличаются, возникает проблема балансировки нагрузки. Данную проблему можно разрешить при помощи динамического распределения блоков между процессами.Пусть все Аii — матрицы размерности mґm, тогда в системе из p=s процессоров
Положим p=4, m=10 и построим график Sp(n) при различных значениях отношения Tsend/Tmul.
Видно, что рост производительности параллельной системы наблюдается как при больших значениях m, так и при больших значениях n.
Вычислим оптимальное число процессоров p’ и построим график Sp(n, m) при n=5000, m=10.
Значение p’ определяется из уравнения
, p і 1
Отсюда
Результаты эксперимента
m
t1(n, m), сек
tp(n, m), p=2, сек
Sp(n)
1000
2131, 162
1097, 971
1, 941
1100
2235, 178
1142, 729
1, 956
1200
3764, 320
1931, 411
1, 949
1300
4898, 564
2506, 942
1, 954
1. 4. Ленточные матрицы
Матрица А порядка n называется ленточной, если
aij = 0, i – j > b1, j – i > b2
Далее будем рассматривать лишь матрицы с симметричной лентой, для которых b1 = b2 = b. Таким образом, ненулевые элементы матрицы расположены только на главной диагонали и на 2b диагоналях, прилегающих к главной сверху и снизу.
Заметим, что если полуширина ленты для А равна a, а для В равна b, то полуширина ленты для АВ в общем случае будет равна a+b.
Пусть А — матрица размерности nґn с полушириной ленты a, В — матрица размерности nґn с полушириной ленты a. Причем А разбита на блоков
Число ненулевых элементов матрицы B не превышает (2a+1)n, в блоке Ai — (2a+1)m. В результирующей матрице будет не более (4a+1)n ненулевых элементов.
В системе из p процессоров построим вычисление по следующей схеме: первый процессор вычисляет А1В, второй — А2В и т. д. Соответственно, матрицу А выгодно хранить по строкам, а матрицу В — по столбцам. В результате умножения получим матрицу, которая будет храниться в памяти по строкам.
На вычисление одного ненулевого элемента матрицы АВ затрачивается время не более (2a + 1)Tmul. Тогда
Видно, что ускорение Sp зависит от полуширины ленты a и не зависит от n.
Теперь вычислим оптимальное значение числа процессоров p’ и построим график функции Sp(a) при a=500.
Значение p’ определяется из уравнения
Отсюда
Г Л А В А 2
ПРЯМЫЕ МЕТОДЫ РЕШЕНИЯ ЛИНЕЙНЫХ СИСТЕМ
Рассмотрим систему линейных уравнений
Ax = b
с невырожденной матрицей А размера nґn .
2. 1. LU-разложение
Постановка задачи
Построим разложение мтрицы A = LU, где L — нижнетреугольная матрица с единицами на главной диагонали, а U — верхнетреугольная матрица.
В качестве алгоритма LU-разложения выберем kij-форму с выбором главного элемента по строке (столбцу). Упрощенно алгоритм можно записать следующим образом
для k = 1 до n-1
для i = k+1 до n
lik = aik/akk
для j = k+1 до n
aij = aij - likakj
Распараллеливание и оценки производительности
Предположим вначале, что система состоит из p = n процессоров. Тогда один из возможных вариантов организации LU-разложения выглядит так. Пусть i-ая строка матрицы А хранится в процессоре i. На первом шаге первая строка рассылается всем процессорам, после чего вычисления
(*)
могут выполняться параллельно процессорами Р2, … , Pn. На втором шаге вторая строка приведенной матрицы рассылается из процессора P2 процессорам P3, … , Pn, и т. д.
На k-ом шаге процессоры P1, … , Pk-1 простаивают, Pk выполняет нормировку и рассылку строки, Pk+1, … , Pn — вычисляют по формуле (*). Рассылка строки выполняется за время (n-k)(n-k)Tsend, нормировка — за время (n-k)Tmul. Параллельная модификация строк выполняется за время (n-k)Tmul. Таким образом
Т. е. при больших значениях n параллельная версия может работать медленнее последовательной версии. Это объясняется значительным объемом обмена данными между каждыми двумя шагами и уменьшением числа активных процессоров на 1 на каждом шаге.
В ситуации, когда p
Для уменьшения простоя процессоров будем использовать опережающую рассылку: на k-м шаге (k+1)-я строка рассылается сразу после того, как окончится ее модификация.
Для обеспечения численной устойчивости алгоритма применим схему выбора главного элемента по строке.
В этом случае получим
При достаточно больших значениях n
Альтернативой для хранения элементов по строкам является столбцовая циклическая схема хранения. В этом случае по процессорам распределяются не строки, а столбцы, что позволяет уменьшить объем передаваемых данных. Главный элемент будем выбирать в активном столбце.
На k-м шаге процессор, хранящий k-й столбец, выбирает главный элемент, вычисляет значения (i > k) и пересылает их на остальные процессоры. Далее каждый процессор вычисляет по формуле
Передача данных происходит за время p(n-k)Tsend, вычисления — за время . Тогда
При достаточно больших значениях n
Зафиксируем p=4 и построим графики функции Sp(n) при различных значениях отношения Tsend/Tmul.
Вычислим оптимальное число процессоров p’ и построим график функции Sp(n) при n=5000.
Значение p’ найдем из уравнения
, p і 1
Отсюда получаем .
Результаты эксперимента
n
t1(n), сек
tp(n), p=2, сек
Sp(n)
1000
450, 975
309, 948
1, 455
1100
600, 691
408, 078
1, 472
1200
780, 679
554, 960
1, 408
1300
993, 420
689, 396
1, 441
2. 2. Решение треугольных систем
Постановка задачи
После выполнения LU-разложения нужно решать треугольные системы.
Ly = b, Ux = y
Процесс их решения назывантся прямой и обратной подстановками. Рассмотрим только решение нижнетреугольной системы (решение верхнетреугольной системы организуется аналогичным образом).
В качестве алгоритма решения выберем алгоритм скалярных произведений. Его запись на псевдокоде выглядит следующим образом
для i = 1 до n
для j = 1 до i-1
bi = bi – Lijyj
yi = bi/Lii
Распараллеливание
Пусть вычислительная система состоит из p процессоров, а размерность матрицы L равна n. Разобьем матрицу L на блоки в соответствии с рисунком
Все блоки (кроме, может быть, двух последних) состоят из строк. Также на p блоков разбиваются вектора y и b.
Далее, пусть блок S1 расположен в памяти процессора P1, блоки L2 и S2 — в памяти процессора P2 и т. д. Перед началом вычислений блок bk пересылается процессору Pk, k = 1...p.
Рассмотрим первый шаг алгоритма решения системы.
1. Процессор P1 решает систему S1y1=b1
2. Найденный вектор y1 пересылается процессорам P2, …, Pp
3. Процессор Pk (k=2...p) пересчитывает вектор bk по формуле
Остальные шаги алгоритма организуются аналогично.
Оценки производительности
Будем считать, что матрица L уже распределена по процессорам. В противном случае выгоднее решать систему последовательным способом.
На одном процессоре задача размерности n решается за время
Оценим время выполнения k-ой итерации в системе из p процессоров
Тогда
Зафиксируем p=4 и построим графики функции Sp(n) при различных значениях отношения Tsend/Tmul.
Вычислим оптимальное число процессоров p’ и построим график функции Sp(n) при n=1500.
Значение p’ найдем из уравнения
, p і 1
Отсюда получаем .
Результаты эксперимента
К сожалению, для данной задачи не удалось получить приемлимых результатов. Дело в том, что на момент проведения экспериментов в сети было недостаточное количество компьютеров, поэтому параллельный алгоритм выполнялся даже дольше, чем последовательный.
2. 3. QR-разложение
Постановка задачи
Рассмотрим ортогональное приведение матрицы A
A = QR
где Q —ортогональная, а R — верхнетреугольная матрицы. В качестве алгоритма решения задачи выберем преобразование Хаусхолдера. На последовательных компьютерах преобразование Хаусхолдера выполняется примерно в два раза медленне, чем LU-разложение. Поэтому данный метод редко используется для решения СЛУ, однако он широко используется при вычислении собственных значений.
Пусть A — матрица nґn. Рассмотрим запись алгоритма на псевдокоде
для k = 1 до n-1
,
uTk = (0, …, 0, akk – sk, ak+1, k, …, ank)
akk = sk
для j = k+1 до n
aj = gkukaj, aj = aj - ajuk
Распараллеливание
Первый шаг алгоритма выглядит следующим образом
Послать первый столбец матрицы А всем процессорам
Вычислить s1 в процессоре 1 и переслать его всем процессорам. Пересчитать матрицу А во всех процессорах в соответствии с алгоритмом Хаусхолдера.
Остальные шаги организуются аналогичным образом.
Оценки производительности
Пусть А — матрица размерности nґn. И пусть система состоит из p процессоров, причем n=pm.
Распределим матрицу А по процессорам в соответствии со столбцовой циклической схемой хранения.
Последовательный алгоритм выполняется за время
Пересылка начальных данных выполняется за время n2Tsend, пересылка результата — за время n2Tsend. Поэтому
Зафиксируем p=4 и построим графики функции Sp(n) при различных значениях отношения Tsend/Tmul.
Построим график функции Sp(n) при n=2000.
Результаты эксперимента
n
t1(n), сек
tp(n), p=2, сек
Sp(n)
1300
1787, 780
1406, 593
1, 271
1400
2242, 265
1719, 528
1, 304
1500
2762, 694
2126, 785
1, 299
1600
4085, 485
2969, 102
1, 376
Г Л А В А 3
ИТЕРАЦИОННЫЕ МЕТОДЫ РЕШЕНИЯ ЛИНЕЙНЫХ УРАВНЕНИЙ
3. 1. Метод Якоби
Постановка задачи
Пусть А — невырожденная матрица nґn и нужно решить систему
Ax = b,
где диагональные элементы aii матрицы А ненулевые. Тогда метод Якоби в матричной форме имеет вид
xk+1 = Hxk + d, k = 0, 1, …, H = D-1B, d = D-1b,
где D = diag(a11, …, ann), B = D – A и указано начальное значение х0.
Конечно, метод Якоби сходится не всегда. Приведем две стандартные теоремы сходимости.
Теорема. Если матрица А имеет строгое диагональное преобладание или является неприводимой с диагональным преобладанием, то итерации метода Якоби сходятся при любом начальном приближении х0.
Теорема. Если матрица A = D – B симметрична и положительно определена, то итерации метода Якоби сходятся при любом начальном приближении х0 тогда и только тогда, когда матрица D + B положительно определена.
Для проверки сходимости итерационных приближений будем использовать следующий критерий:
где .
Распараллеливание
Пусть система состоит из р процессоров и матрицы разбиты на блоки следующим образом
, ,
Припишем блоки H1, d1, xk1 процессору Р1, блоки H2, d2, xk2 — процессору Р2 и т. д. Тогда соответствующие компоненты вектора х можно вычислять параллельно по схеме:
Шаг 1. Переслать на i-й процессор Hi, di, x0.
Шаг 2. На к-й итерации i-й процессор выполняет следующие действия: хiк+1 = Hixik + di
переслать xik+1 остальным процессорам
получить недостающие компоненты вектора xk+1
проверить условие сходимости и, если условие выполнено, установить “флаг” Шаг 3. Если хотябы один флаг не установлен, то перейти к следующей итерации (Шаг 2).
Шаг 4. Передать результат.
Оценки производительности
Для простоты будем считать, что n = pm. Оценим время выполнения одной итерации метода Якоби.
Далее, пусть метод сошелся за k шагов. Тогда
Из последней формулы несложно установить, что если число итераций , то каким бы большим небыло значение n, нам все равно не удастся достичь ускорения Sp і 2.
Зафиксируем p=4, k=500 и построим график функции Sp(n) при различных значениях отношения Tsend/Tmul.
Кривая, которой соответствует значение Tsend/Tmul=1000 никогда не поднимется выше 1, т. к. число итераций k слишком мало.
Вычислим оптимальное число процессоров p’ и построим график функции Sp(n, k) при n=5000, k=500.
Значение p’ найдем из уравнения
, p і 1
Отсюда . Как и следовало ожидать, при больших значениях k значение p’ практически не зависит от k (уже при k=100 отношение ).
Результаты эксперимента
n
t1(n), сек
tp(n), p=2, сек
Sp(n)
1000
489, 634
333, 535
1, 468
1100
608, 457
404, 232
1, 505
1200
832, 035
514, 840
1, 616
1300
857, 569
559, 942
1, 531
1400
992, 017
661, 823
1, 499
ЗАКЛЮЧЕНИЕ
Результатом данной работы является программная среда FLOWer, предназначенная для конструирования параллельных программ. Основные идеи, включая язык DGL, на которых базируется система, были разработаны самостоятельно.
Для испытания возможностей системы были реализованы некоторые параллельные алгоритмы ВЛА, даны теоретические оценки их эффективности и проведены вычислительные эксперименты, которые показали состоятельность и эффективность системы.
Литература
Ортега Дж. Введение в параллельные и векторные методы решения линейных систем. - М. : Мир, 1991.
Водяхо А. И. , Горнец Н. Н. , Пузанков Д. В. Высокопроизводительные системы обработки данных. - М. : Высшая школа, 1997.
Бэбб Р. Программирование на параллельных вычислительных системах. - М. : Мир, 1991.
Белецкий В. Н. Многопроцессорные и параллельные структуры с организацией асинхронных вычислений. Киев: Наукова думка, 1988.
http: //parallel. ru. Информационно-аналитический центр по параллельным вычислениям.
