Участник:Малков Кирилл/Алгоритм бидиагонализации матрицы методом отражений

Свойства и структура алгоритмов

Общее описание алгоритма

Матрица [math]\displaystyle{ B = [b_{ij}] }[/math] называется (верхней) двухдиагональной, или бидиагональной, если [math]\displaystyle{ b_{ij} = 0 }[/math], [math]\displaystyle{ i \gt j }[/math] или [math]\displaystyle{ i + 1 \lt j }[/math]. Алгоритм бидиагонализации матрицы методом отражений приводит произвольную [math]\displaystyle{ n \times n }[/math] матрицу A к бидиагональному виду [math]\displaystyle{ B = PAQ }[/math], где [math]\displaystyle{ P }[/math] и [math]\displaystyle{ Q }[/math] - произведения конечного числа матриц отражения(или вращения).

Математическое описание алгоритма

Пусть дана матрица A размера [math]\displaystyle{ n \times n }[/math]. Вычисляется двухдиагональная матрица D размера [math]\displaystyle{ n \times n }[/math], такая что [math]\displaystyle{ A = P \times D \times Q }[/math] для некоторых матриц P и Q, являющихся произведением некоторого числа матриц отражения.Сначала мы умножаем A слева на матрицу отражения, аннулирующую все поддиагональные элементы вервого столбца. Затем умножаем результат справа на матрицу отражения, аннулирующую элементы первой строки в позициях с 3-й по n-ю. Далее умножением слева аннулируем все поддиагональные элементы второго столбца, затем умножением справа аннулируем получаем нули во второй строке в позициях с 4-й по n-ю и т.д. После каждого умножения на матрицу отражения все ранее полученные нули остаются. Формулы метода: Пусть [math]\displaystyle{ A = [a_1, a_2, ..., a_n] }[/math] [math]\displaystyle{ D_1 = A, P_1 = I, Q_1 = I. }[/math] [math]\displaystyle{ Если i четное, то D_i = D_{i-1} * (I - 2 * u_{i-1} * u_{i-1}^*) }[/math] [math]\displaystyle{ u_i = frac{[0, 0, ..., 0, d_{i+1} - ||d_{i+1}||_2, d_{i+2}, ..., d_n]}{||[0, 0, ..., 0, d_{i+1} - ||d_{i+1}||_2, d_{i+2}, ..., d_n]||_2}, Q_i = (I - 2 * u_{i} * u_{i}^*) * Q_{i-1}, P_i = P_{i-1}. }[/math] [math]\displaystyle{ Если i нечетное, то D_i = (I - 2 * u_{i-1} * u_{i-1}^*) * D_{i-1}, }[/math] [math]\displaystyle{ u_i = frac{[0, 0, ..., 0, d_i - ||d_i||_2, d_{i+1}, ..., d_n]}{||[0, 0, ..., 0, d_i - ||d_i||_2, d_{i+1}, ..., d_n]||_2}, Q_i = Q_{i-1}, P_i = P_{i-1} * (I - 2 * u_{i} * u_{i}^*). }[/math]

Вычислительное ядро алгоритма

Вычислительное ядро алгоритма составляет нахождение на каждом этапе матрицы отражения(всего [math]\displaystyle{ 2 * n - 2 }[/math].

Макроструктура алгоритма

Операции скалярного произведения, умножения матрицы на вектор

Схема реализации последовательного алгоритма

[math]\displaystyle{ D_1 = A, P_1 = I, Q_1 = I. }[/math] [math]\displaystyle{ 1) для четного i : D_i = D_{i-1} * (I - 2 * u_{i-1} * u_{i-1}^*), Q_i = (I - 2 * u_{i} * u_{i}^*) * Q_{i-1}, P_i = P_{i-1}; }[/math] [math]\displaystyle{ для нечетного i : D_i = (I - 2 * u_{i-1} * u_{i-1}^*) * D_{i-1}, Q_i = Q_{i-1}, P_i = P_{i-1} * (I - 2 * u_{i} * u_{i}^*). }[/math] [math]\displaystyle{ 2) для четного i : u_i = \frac{[0, 0, ..., 0, d_{i+1} - ||d_{i+1}||_2, d_{i+2}, ..., d_n]}{||[0, 0, ..., 0, d_{i+1} - ||d_{i+1}||_2, d_{i+2}, ..., d_n]||_2}, }[/math] [math]\displaystyle{ для нечетного i : u_i = \frac{[0, 0, ..., 0, d_i - ||d_i||_2, d_{i+1}, ..., d_n]}{||[0, 0, ..., 0, d_i - ||d_i||_2, d_{i+1}, ..., d_n]||_2}. }[/math]

Последовательная сложность

Для двухдиагонализации матрицы порядка n методом отражений при последовательном варианте требуется: <math>1) 2 * n - 2 вычислений квадратного корня, 2) 4 * n - 4 делений, 3) (4 * n - 4)(1 + n + n^2) умножений, 4) (4 * n - 4)(1 + n + n^2) сложений. Сложность последовательного алгоритма есть <math>O(n^3)<\math>.

Информационный граф

=== При параллельном вычислении двухдиагональной матрицы методом отражений требуется последовательно выполнить следующие ярусы: <math> -- k ярусов с вычислением 4 * (n/k) - 4 делений, -- k ярусов с вычислением (4 * (n/k) - 4)(1 + n + n^2) умножений, -- k ярусов с вычислением (4 * (n/k) - 4)(1 + n + n^2) сложений. Получим, что параллельная сложность алгоритма есть <math>O(n^2)<\math>

Входные и выходные данные алгоритма

На вход подается матрица A и ее размер n. Выходные данные -- матрицы P, Q, D, такие что <math>A = P * D * Q <\math>

Свойства алгоритма

Соотношение последовательной и параллельной сложности в случае неограниченных ресурсов является линейным. Вычислительная мощность алгоритма, как отношение числа операций к суммарному объему входных и выходных данных, линейна.

Литература

Е. Е. Тыртышников, "Методы численного анализа".