Другие журналы
|
научное издание МГТУ им. Н.Э. БауманаНАУКА и ОБРАЗОВАНИЕИздатель ФГБОУ ВПО "МГТУ им. Н.Э. Баумана". Эл № ФС 77 - 48211. ISSN 1994-0408![]()
Паралеллный алгоритм решения систем линейных алгебраических уравнений с многодиагональной матрицей коэффициентов
# 07, июль 2013 DOI: 10.7463/0713.0590785
Файл статьи:
![]() УДК 004.942, 519.876.5 Россия, МГТУ им. Н.Э. Баумана
Введение При решении задач математической физики численными методами, такими как метод конечных элементов [1], метод конечных разностей [2], часто возникает необходимость нахождения решения систем линейных алгебраический уравнений (СЛАУ) [3] высоких порядков, матрица коэффициентов которых имеет выраженный диагональный (ленточный) характер. С распространением многопроцессорных вычислительных систем появилась возможность адаптации точных алгоритмов решения СЛАУ к данной архитектуре. Работы в этом направлении ведутся продолжительное время [4], разработано большое количество алгоритмов для вычислительных систем различной архитектуры [5-7]. Однако, предложенные алгоритмы рассматривают диагональную ленту матрицы коэффициентов как полностью заполненную ненулевыми элементами, хотя в практических задачах эта лента, как правило, содержит большое количество ненулевых диагоналей. Поэтому актуальной является задача разработки алгоритма, учитывающего внутреннюю разреженность диагональной ленты матрицы коэффициентов, что должно повысить его экономичность. В работе предложен алгоритм параллельного решения СЛАУ с многодиагональной (многоленточной) матрицей коэффициентов на многопроцессорной вычислительной системе с общей памятью. Данный алгоритм реализует метод Гаусса [3] и состоит из трёх основных этапов: трансформация задачи, прямой ход алгоритма, обратный ход алгоритма.
1. Постановка задачи Задача нахождения решения СЛАУ может быть представлена в виде
где A - квадратная ленточная матрица коэффициентов с шириной ленты L, x - вектор-столбец неизвестных, с – вектор-столбец свободных членов. В общем случае A, x и c имеют вид
Под решением задачи понимают такую совокупность чисел Следует отметить следущую особенность решения СЛАУ с ленточной матрицей коэффициентов. Исходная матрица в практических задачах является сильно разреженной, однако, в ходе решения прямыми методами (таким, например, как метод Гаусса) лента полностью заполняется новыми ненулевыми элементами, что резко замедляет ход как последовательного, так и параллельного решения.
2. Трансформация задачи На данном этапе производится изменение задачи (1) с целью уменьшения количества новых ненулевых элементов матрицы коэффициентов и упрощения процедуры распределения нагрузки между процессорами вычислительной системы. Пусть N - число диагоналей, содержащих отличные от нуля элементы матрицы A(ненулевые диагонали), Введем смещение Расширим вектор xзадачи (1), добавив в его начало Sпеременных, а также введём в систему где
где A` - матрица коэффициентов, x` - вектор-столбец неизвестных и с` - вектор-столбец свободных членов. Размерность матрицы Далее к строкам матрицы Рассмотрим пример трансформации исходной задачи (1) размерности 5x5 со смещением
Отметим, что в данной задаче D= (d1, d2, d3). Модифицируем задачу Так как S = 2, добавим две дополнительные переменные и два уравнения в задачу (4):
Прибавим последние две строки матрицы
Отметим, что в модифицированной задаче (5) векторы ненулевых диагоналей имеют вид На основе теоремы Кронекера – Капелли легко показать существование и единственность решения задачи (3) при условии существования и единственности решения задачи (1). Необходимо отметить, что проведённая трансформация системы уравнений привела её матрицу коэффициентов к верхнетреугольному виду с нижним окаймлением. Это означает, что в прямом ходе метода Гаусса все подлежащие исключению ненулевые элементы матрицы (в том числе и вновь возникающие) располагаются только в её S последних строках.
3. Представление данных В задаче (3), как правило, матрица Матрица · n – размерность одной строки матрицы A; · S – смещение матрицы A; · L– ширина ленты матрицы A; · N – число ненулевых диагоналей матрицы A; · · p– вектор размерности Вектора cи xсохраняют свое представление, указанное в (2). Для матрицы A` дополнительно введем векторы При таком подходе к представлению данных система из 10000 уравнений с
4. Прямой ход алгоритма Первоначально первые S уравнений задачи (3) вычитаются из последних S уравнений системы по следующему правилу: из уравнения с номером
Тогда векторы Далее начинается итерационный процесс, состоящий из
В выражении (6) индекс i- номер итерации, На втором шаге итерации выполняется циклический сдвиг векторов Из выражения (6) легко видеть возможность независимого вычисления компонентов векторов После выполнения
При этом векторы Перед выполнением следующего этапа требуется синхронизировать работу программы. Для последних S уравнений задачи (3) выполняется прямой ход метода Гаусса. Для примера (4) имеем
Тогда вектора
5. Обратный ход алгоритма Обратный ход алгоритма полностью совпадает с аналогичным этапом канонического метода Гаусса.
6. Схема метода для многопроцессорной вычислительной системы Рассматриваем вычислительную систему с общей памятью, включающую p процессоров. Схема метода решения задачи на вычислительной системе данной архитектуры выглядит следующим образом. 1) Какой-либо из процессоров назначается главным (master-процессор), а все остальные - вспомогательными (slave-процессоры). На master-процессоре запускается основной процесс алгоритма. 2) В оперативную память системы загружаются исходные данные решаемой задачи. 3) На master–процессоре выполняем модификацию задачи (1) согласно изложенному выше алгоритму. 4) Запускается 5) Каждый поток управления, используя (6), производит модификацию векторов 6) На master-процессоре выполняется прямой ход стандартного метода Гаусса для последних S уравнений модифицированной задачи. 7) На master-процессоре выполняется обратный ход метода Гаусса.
7. Аналитическая оценка эффективности алгоритма Для начала оценим время выполнения итерации алгоритма на одном процессоре. На каждой итерации требуется выполнить не более чем S операций деления, Пусть среднее время выполнения одной операции деления над числами с плавающей точкой двойной точности задаётся константой
где Оценим время выполнения данного алгоритма на многопроцессорной системе с общей оперативной памятью, используя рассмотренную выше схему метода:
В (8) время
В формуле (9) величина
Отсюда получаем выражение для ускорения вычислительного процесса
Ниже приведены результаты аналитического исследования алгоритма. Пример 1. Отобразим зависимость
Рисунок 1. Зависимость
На рисунке 2 дан график зависимости ускорения от смещения S ленты матрицы коэффициентов при фиксированном количестве процессов p=8.
Рисунок 2. Зависимость
Пример 2. Здесь приводятся зависимости ускорения
Рисунок 3. Зависимость
Рисунок 4. Зависимость
Пример 3. Здесь даны зависимости ускорения
Рисунок 5. Зависимость
Рисунок 6. Зависимость
Пример 4. Здесь исследуется зависимость ускорения
Рисунок 7. Зависимость
Из рассмотренных выше примеров очевидно, что ускорение слабо зависит от размерности задачи, а зависит лишь от смещения. Ступеньки на рисунках 1, 3, 5, 7 появились из-за несбалансированности загрузки системы, полученной в результате некратности чисел S и p. Так при увеличении p не всегда стоит ожидать кратного увеличения ускорения. Из рисунков 2, 4, 6 видно, как с ростом смещения падает ускорение, это связано с преобладанием слагаемого
8. Программная реализация Алгоритм реализован на языке высокого уровня C++ [8] с использованием стандартных библиотек многопоточного программирования операционной системы UNIX [9]. Ниже приведены результаты численных экспериментов для некоторых СЛАУ.
n = 100000, S = 400, N = 5.
n = 100000, S = 700, N =5.
n = 3000000, S = 400, N = 5.
n = 3000000, S = 700, N = 5.
Далее для n = 3000000, N = 5 приведена зависимость ускорения
Данные результаты полностью согласуются с аналитической оценкой эффективности алгоритма. Так с ростом ширины ленты эффективность алгоритма падает, а при одинаковой ширине ленты ускорение практически не зависит от размерности исходной задачи n.
Заключение
В статье предложен новый параллельный алгоритм решения систем линейных алгебраических уравнений с диагонально-ленточной матрицей коэффициентов, учитывающий наличие нулевых элементов в отдельных диагоналях внутри ленты. Основными достоинствами алгоритма являются следующие. - Простота программной реализации. - Высокие показатели ускорения для задач с малым смещением ленты матрицы коэффициентов. - Экономия памяти вычислительной системы. - Возможность простой адаптации под различные вычислительные системы с общей памятью. - Отсутствие множественного доступа по записи к исходным данным. К недостаткам алгоритма относятся следующие. - Падение ускорения с увеличением смещения S ленты матрицы A. - Необходимость выделения дополнительной памяти на этапе модификации исходной задачи. - Возможная несбалансированность нагрузки процессоров вычислительной системы.
Список литературы 1. Ильин В.А., Позняк Э.Г. Линейная алгебра. М.: Физматлит, 1999. 280 с. 2. Ильин В.А., Ким Г.Д. Линейная алгебра и аналитическая геометрия. М.: Проспект, 2007. 400 с. 3. Ортега Дж. Введение в параллельные и векторные методы решения линейных систем: пер. с англ. М.: Мир, 1991. 367 с. [Ortega J.M. Introduction to Parallel and Vector Solution of Linear Systems. NY.: Plenum Press, 1988. 305 p.] 4. Tewarson R.P. Solution of linear equations with coefficient matrix in band form // BIT Numerical Mathematics. 1968. Vol. 8, iss. 1. P. 53-58. DOI: 10.1007/BF01939979 5. Saad Y., Schultz H. Parallel direct methods for solving banded linear systems // Linear Algebra Appl. 1987. Vol. 88-89. P. 623-650. http://dx.doi.org/10.1016/0024-3795(87)90128-5 6. Polizzi E., Sameh A. SPIKE: A parallel environment for solving banded linear systems // Computers & Fluids. 2007. Vol. 36, no. 1. P. 113-141. http://dx.doi.org/10.1016/j.compfluid.2005.07.005 7. Krüger J., Westermann R. GPU Gems 2. Chapter 44. A GPU Framework for Solving Systems of Linear Equations. Режим доступа: http://http.developer.nvidia.com/GPUGems2/gpugems2_chapter44.html (дата обращения 19.06.2013). 8. Эккель Б. Философия C++. Введение в стандартный C++. СПб.: Питер, 2004. 572 с. [Eckel B. Thinking in C++: Introduction to Standard C++. NY.: PrenticeHall, 2000. 814 p.] 9. Чан Т. Системное программирование на C++ для Unix. М.: BHV, 1997. 490 с. [Chan T. UNIX System Programming using C++. Prentice Hall, 1997. 598 p.] Публикации с ключевыми словами: параллельный алгоритм, система линейных алгебраических уравнений, многодиагональная матрица коэффициентов, многоленточная матрица коэффициентов, многопроцессорная вычислительная система Публикации со словами: параллельный алгоритм, система линейных алгебраических уравнений, многодиагональная матрица коэффициентов, многоленточная матрица коэффициентов, многопроцессорная вычислительная система Смотри также: Тематические рубрики: Поделиться:
|
|
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|