Прямая подстановка (вещественный вариант)
Прямая подстановка для нижней треугольной матрицы с единичной диагональю | |
Последовательный алгоритм | |
Последовательная сложность | [math]O(n^2)[/math] |
Объём входных данных | [math]O(n^2)[/math] |
Объём выходных данных | [math]n[/math] |
Параллельный алгоритм | |
Высота ярусно-параллельной формы | [math]O(n)[/math] |
Ширина ярусно-параллельной формы | [math]O(n)[/math] |
Основные авторы описания: А.В.Фролов.
Содержание
- 1 Свойства и структура алгоритма
- 1.1 Общее описание алгоритма
- 1.2 Математическое описание алгоритма
- 1.3 Вычислительное ядро алгоритма
- 1.4 Макроструктура алгоритма
- 1.5 Схема реализации последовательного алгоритма
- 1.6 Последовательная сложность алгоритма
- 1.7 Информационный граф
- 1.8 Ресурс параллелизма алгоритма
- 1.9 Входные и выходные данные алгоритма
- 1.10 Свойства алгоритма
- 2 Программная реализация алгоритма
- 3 Литература
1 Свойства и структура алгоритма
1.1 Общее описание алгоритма
Прямая подстановка - решение СЛАУ [math]Lx = y[/math] с нижней треугольной матрицей [math]L[/math]. Матрица [math]L[/math] - одна из составляющих матрицы [math]A[/math] и получается либо из [math]LU[/math]-разложения последней каким-либо из многочисленных способов (например, простое разложение Гаусса, разложение Гаусса с выбором ведущего элемента, компактная схема Гаусса), либо из других разложений. В силу треугольности [math]L[/math] решение СЛАУ является одной из модификаций общего метода подстановки и записывается простыми формулами.
В[1] метод решения СЛАУ с левой треугольной матрицей назван методом обратной подстановки. Там же отмечено, что в литературе иногда под обратной подстановкой имеют в виду, как и здесь, только решения СЛАУ с правой треугольной матрицей, а решение левых треугольных систем называют прямой подстановкой. Такой же системы названий будем придерживаться и мы, во избежание одноимённого названия разных алгоритмов. Кроме того, обратная подстановка, представленная в этой энциклопедии алгоритмов, одновременно является частью метода Гаусса для решения СЛАУ, а именно - его обратным ходом, чего нельзя сказать про представленную здесь прямую подстановку.
Общая структура прямой подстановки с неособенной нижней треугольной матрицей, тем не менее, практически полностью совпадает со структурой обратной подстановки. Поэтому здесь мы рассмотрим случай, когда матрица [math]L[/math], как полученная из разложения типа Гаусса, имеет единичные диагональные элементы.
1.2 Математическое описание алгоритма
Исходные данные: нижняя треугольная матрица [math]L[/math] (элементы [math]l_{ij}[/math]), вектор правой части [math]b[/math] (элементы [math]b_{i}[/math]).
Вычисляемые данные: вектор решения [math]y[/math] (элементы [math]y_{i}[/math]).
Формулы метода:
- [math] \begin{align} y_{1} & = b_{1} \\ y_{i} & = b_{i} - \sum_{j = 1}^{i-1} l_{ij} y_{j}, \quad i \in [2, n]. \end{align} [/math]
Существует также блочная версия метода, однако в данном описании разобран только точечный метод.
1.3 Вычислительное ядро алгоритма
Вычислительное ядро прямой подстановки можно составить из множественных (всего их [math]n-1[/math]) вычислений скалярных произведений строк матрицы [math]L[/math] на уже вычисленную часть вектора [math]y[/math]:
- [math] \sum_{j = 1}^{i-1} l_{ij} y_{j} [/math]
в режиме накопления или без него, в зависимости от требований задачи, с их последующим вычитанием из компоненты вектора [math]b[/math] и деления на диагональный элемент матрицы [math]L[/math]. В отечественных реализациях, даже в последовательных, упомянутый способ представления не используется. Дело в том, что даже в этих реализациях метода вычисление сумм типа
- [math] b_{i} - \sum_{j = 1}^{i-1} l_{ij} y_{j} [/math]
в которых и встречаются скалярные произведения, ведутся не в порядке «вычислили скалярное произведение, а потом вычли его из элемента», а путём вычитания из элемента покомпонентных произведений, являющихся частями скалярных произведений. Поэтому следует считать вычислительным ядром метода прямой подстановки не вычисления скалярных произведений, а вычисления выражений
- [math] b_{i} - \sum_{j = 1}^{i-1} l_{ij} y_{j} [/math]
в режиме накопления или без него, в зависимости от требований задачи.
1.4 Макроструктура алгоритма
Как уже записано в описании ядра алгоритма, основную часть метода прямой подстановки составляют множественные (всего [math]n-1[/math]) вычисления сумм
- [math] b_{i} - \sum_{j = 1}^{i-1} l_{ij} y_{j} [/math]
в режиме накопления или без него.
1.5 Схема реализации последовательного алгоритма
Чтобы понять последовательность исполнения, перепишем формулы метода так:
1. [math]y_{1} = b_{1}[/math]
Далее для всех [math]i[/math] от [math]2[/math] до [math]n[/math] по возрастанию выполняются
2. [math]y_{i} = b_{i} - \sum_{j = 1}^{i-1} l_{ij} y_{j}[/math]
Особо отметим, что вычисления сумм вида [math]b_{i} - \sum_{j = 1}^{i-1} l_{ij} y_{j}[/math] производят в режиме накопления вычитанием из [math]b_{i}[/math] произведений [math]l_{ij} y_{j}[/math] для [math]j[/math] от [math]1[/math] до [math]i-1[/math], c возрастанием [math]j[/math]. Другие порядки выполнения суммирования приводят к резкому ухудшению параллельных свойств алгоритма.
1.6 Последовательная сложность алгоритма
Для прямой подстановки порядка n в последовательном (наиболее быстром) варианте требуется:
- по [math]\frac{n^2-n}{2}[/math] умножений и сложений (вычитаний).
При этом использование режима накопления требует совершения умножений и вычитаний в режиме двойной точности (или использования функции вроде DPROD в Фортране), что ещё больше увеличивает затраты во времени, требуемом для выполнения прямой подстановки.
При классификации по последовательной сложности, таким образом, метод обратной подстановки относится к алгоритмам с квадратической сложностью.
1.7 Информационный граф
Опишем граф алгоритма как аналитически, так и в виде рисунка.
Граф алгоритма прямой подстановки состоит из одной группы вершин, расположенной в целочисленных узлах двумерной области, соответствующая ей операция [math]a-bc[/math].
Естественно введённые координаты области таковы:
- [math]i[/math] — меняется в диапазоне от [math]2[/math] до [math]n[/math], принимая все целочисленные значения;
- [math]j[/math] — меняется в диапазоне от [math]1[/math] до [math]i-1[/math], принимая все целочисленные значения.
Аргументы операции следующие:
- [math]a[/math]:
- при [math]j = 1[/math] элемент входных данных [math]b_{i}[/math];
- при [math]j \gt 1[/math] — результат срабатывания операции, соответствующей вершине с координатами [math]i, j-1[/math];
- [math]b[/math] — элемент входных данных, а именно [math]l_{ij}[/math];
- [math]c[/math]:
- при [math]j = 1[/math] элемент входных данных [math]b_{1}[/math];
- при [math]j \gt 1[/math] — результат срабатывания операции, соответствующей вершине с координатой [math]j, j-1[/math];
Результат срабатывания операции является:
- при [math]j \lt i-1[/math] - промежуточным данным алгоритма;
- при [math]j = i-1[/math] - выходным данным.
Описанный граф можно посмотреть на рисунке, выполненном для случая [math]n = 5[/math]. Здесь вершины обозначены зелёным цветом. Подача входных данных из вектора [math]b[/math], идущая в вершины левого столбца, и матрицы [math]L[/math], идущая во все вершины, на рисунке не представлена.
1.8 Ресурс параллелизма алгоритма
Для прямой подстановки порядка n в параллельном варианте требуется последовательно выполнить следующие ярусы:
- по [math]n - 1[/math] ярусов умножений и сложений/вычитаний (в каждом из ярусов — линейное количество операций, от [math]1[/math] до [math]n-1[/math].
При этом использование режима накопления требует совершения умножений и вычитаний в режиме двойной точности, а в параллельном варианте это означает, что практически все промежуточные вычисления для выполнения алгоритма в режиме накопления должны быть двойной точности. В отличие от последовательного варианта это означает некоторое увеличение требуемой памяти.
При классификации по высоте ЯПФ, таким образом, метод прямой подстановки относится к алгоритмам с линейной сложностью. При классификации по ширине ЯПФ его сложность также будет линейной.
1.9 Входные и выходные данные алгоритма
Входные данные: нижняя треугольная матрица [math]L[/math] (элементы [math]l_{ij}[/math]), вектор правой части [math]b[/math] (элементы [math]b_{i}[/math]).
Объём входных данных: [math]\frac{n (n + 1)}{2}[/math] (в силу треугольности и единичности диагональных элементов достаточно хранить только поддиагональные элементы матрицы [math]L[/math]).
Выходные данные: вектор решения [math]y[/math] (элементы [math]y_{i}[/math]).
Объём выходных данных: [math]n~.[/math]
1.10 Свойства алгоритма
Соотношение последовательной и параллельной сложности в случае неограниченных ресурсов, как хорошо видно, является линейным (отношение квадратической к линейной).
При этом вычислительная мощность алгоритма прямой подстановки, как отношение числа операций к суммарному объему входных и выходных данных – всего лишь константа.
При этом алгоритм прямой подстановки полностью детерминирован. Использование другого порядка выполнения ассоциативных операций в данной версии нами не рассматривается, поскольку в корне меняет структуру алгоритма и делает параллельную сложность квадратической.
2 Программная реализация алгоритма
2.1 Особенности реализации последовательного алгоритма
В простейшем варианте метод прямой подстановки на Фортране можно записать так:
Y(1) = B(1)
DO I = 2, N-1
S = B(I)
DO J = 1, I-1
S = S - DPROD(L(I,J), Y(J))
END DO
END DO
При этом для реализации режима накопления переменная [math]S[/math] должна быть двойной точности.
2.2 Возможные способы и особенности параллельной реализации алгоритма
2.3 Результаты прогонов
2.4 Выводы для классов архитектур
3 Литература
- ↑ В.В.Воеводин, Ю.А.Кузнецов. Матрицы и вычисления. М.: Наука, 1984, стр. 182.