Участник:Artyom-amg/Бидиагонализация матриц вращениями Гивенса
Зотов А.Э.
1 Свойства и структура алгоритма
1.1 Общее описание алгоритма
Алгоритм бидиагонализации матрицы — это численный метод в линейной алгебре, предназначенный для получения бидиагональной факторизации матрицы. Этот алгоритм, например, может быть использован для последующего решения проблемы сингулярного разложения, то есть отыскания всех сингуряных чисел и векторов матрицы.
Метод бидиагонализации используется для разложения матрицы в виде [math]A = VBU^*[/math] ( [math]V, U[/math] - унитарные матрицы; [math]B[/math] - бидиагональная матрица, т.е. матрица, все ненулевые элементы которой находятся на главной диагонали и на одной из под- или наддиагонали). При этом матрицы [math]V, U[/math] хранятся и используется не в явном виде, а в виде произведения матриц вращения. Каждая из матриц вращения (Гивенса)
[math]G_{ij} = \begin{bmatrix} 1 & \cdots & 0\quad 0\quad 0 & \cdots & 0\quad 0\quad 0 & \cdots & 0 \\ \vdots & \ddots & \vdots & \vdots & \vdots & \vdots \\ 0 & \cdots & 1\quad 0\quad 0 & \cdots & 0\quad 0\quad 0 & \cdots & 0 \\ 0 & \cdots & 0\quad c\quad 0 & \cdots & 0\ -s\ 0 & \cdots & 0 \\ 0 & \cdots & 0\quad 0\quad 1 & \cdots & 0\quad 0\quad 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots \\ 0 & \cdots & 0\quad 0\quad 0 & \cdots & 1\quad 0\quad 0 & \cdots & 0 \\ 0 & \cdots & 0\quad s\quad 0 & \cdots & 0\quad c\quad 0 & \cdots & 0 \\ 0 & \cdots & 0\quad 0\quad 0 & \cdots & 0\quad 0\quad 1 & \cdots & 0 \\ \vdots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & \cdots & 0\quad 0\quad 0 & \cdots & 0\quad 0\quad 0 & \cdots & 1 \\ \end{bmatrix}[/math]
может быть определена парой индексов и одним параметром. Это позволяет результаты разложения на месте матрицы [math]A[/math] без использования дополнительных массивов. На каждом шаге алгоритма две строки или два столбца преобразуемой матрицы подвергаются преобразованию вращения. При этом параметр такого преобразования подбирается так, чтобы один из поддиагональных или наддиагольныхх элементов преобразуемой матрицы стал нулевым. Сначала в нуль последовательно обращаются элементы 1-го столбца и 1-ой строки, затем 2-го слолбца и 2-ой строки, и т.д. до n-1-ого столбца и n-2-ой строки, после чего получившаяся матрица - это и есть матрица [math]B[/math].
Сам шаг алгоритма распадается на две части: подбор параметра вращения и выполнение вращения. Поскольку слева от "рабочего" столбца элементы вращаемых строк равны 0, то там выполнять вращение не нужно. Аналогично ситуация со строками: поскольку выше "рабочей" строк эдементы вращаемых столбцов равны 0, то и там выполнять вращения не нужно. Вращение элементов выполняется одновременно с подбором параметра вращения. Таким образом, вторая часть этапа заключается в выполнении вращения двумерных векторов, составленных из элементов вращаемых строк или столбцов.
Обычно для хранения матрицы вращения используется тангенс половинного угла [math]t[/math], с которым простыми формулами связаны косинус [math]c[/math] и синус [math]s[/math] самого угла: [math]c = (1 - t^2)/(1 + t^2), s = 2t/(1 + t^2)[/math].
1.2 Математическое описание алгоритм
Для получения [math]A = VBU^*[/math] разложения квадратной матрицы [math]A[/math] последнюю приводят к бидиагональной форме последовательностью умножений её слева и справа соответственно на [math] X_{1 2}, X_{1 3}, ..., X_{1 n}, Y_{2 3}^*, Y_{2 4}^*, ..., Y_{2 n}^*, ..., X_{n-2 n}, X_{n-1 n}, Y_{n-1 n}^*[/math], где [math]X_{i j}, Y_{i j}[/math] - матрицы Гивенса.
Далее рассмотрим операции, производимые над строками матрицы, т.е. уможения справа на матрицы вида [math] X_{i j} [/math]. Для столбцов все будет аналогично.
Каждая из матриц [math]G_{i j}[/math] задаёт вращение в 2-мерном подпространстве, связанном с [math]i[/math]-й и [math]j[/math]-й компонентами столбцов, остальные компоненты не меняются. При этом вращение подбирается так, чтобы элемент новой матрицы в [math]j[/math]-й строке и [math]i[/math]-м столбце стал нулевым. Поскольку вращения и тождественные преобразования нулевых векторов оставляют их нулевыми, то последующие вращения сохраняют полученные ранее нули слева и сверху от текущего обнуляемого элемента.
В конце процесса имеем [math]B=X_{n-1 n}X_{n-2 n}X_{n-2 n-1}...X_{1 3}X_{1 2}AY_{2 3}^*Y_{2 4}^*...Y_{2 n}^*...Y_{n-2 n-1}^*Y_{n-2 n}^*Y_{n-1 n}^*[/math]
Так как матрицы вращения унитарны, то естественно получается [math]V=(X_{n-1 n}X_{n-2 n}X_{n-2 n-1}...X_{1 3}X_{1 2})^* =X_{1 2}^* X_{1 3}^* ...X_{1 n}^* X_{2 3}^* X_{2 4}^* ...X_{2 n}^* ...X_{n-2 n}^* X_{n-1 n}^*[/math]
[math]U=(Y_{2 3}^*Y_{2 4}^*...Y_{2 n}^*...Y_{n-2 n-1}^*Y_{n-2 n}^*Y_{n-1 n}^*)^* =Y_{n-1 n}Y_{n-2 n}...Y_{2 n}...Y_{2 4}Y_{2 3}[/math]
Получаем разложение: [math]A = VBU^*[/math].
2 Литература
1. Тыртышников Е.Е. Методы численного анализа
