QR-алгоритм

QR-алгоритм — это численный метод в линейной алгебре, предназначенный для решения полной проблемы собственных значений, то есть отыскания всех собственных чисел и собственных векторов матрицы. Был разработан в конце 1950-х годов независимо В. Н. Кублановской и Дж. Фрэнсисом^[англ.].

Алгоритм[править]

Пусть A — вещественная матрица, для которой мы хотим найти собственные числа и векторы. Положим A₀=A. На k-м шаге (начиная с k = 0) вычислим QR-разложение A_k=Q_kR_k, где Q_k — ортогональная матрица (то есть Q_k^T = Q_k⁻¹), а R_k — верхняя треугольная матрица. Затем мы определяем A_k+1 = R_kQ_k.

Заметим, что

A_{k + 1} = R_{k} Q_{k} = Q_{k}^{- 1} Q_{k} R_{k} Q_{k} = Q_{k}^{- 1} A_{k} Q_{k} = Q_{k}^{T} A_{k} Q_{k},

то есть все матрицы A_k являются подобными, то есть их собственные значения равны.

Пусть все диагональные миноры матрицы A не вырождены. Тогда последовательность матриц A_k при $k \to \infty,$ сходится по форме к клеточному правому треугольному виду, соответствующему клеткам с одинаковыми по модулю собственными значениями.^[1]

Для того, чтобы получить собственные векторы матрицы, нужно перемножить все матрицы Q_k.

Алгоритм считается вычислительно устойчивым, т. к. производится ортогональными преобразованиями подобия.

Доказательство для симметричной положительно определённой матрицы[править]

Будем считать, что собственные числа положительно-определённой матрицы A упорядочены по убыванию:

λ_{1} > λ_{2} > . . . > λ_{n} > 0 .

Пусть

Λ = d i a g (λ_{1}, . . ., λ_{n}),

а S — матрица, составленная из собственных векторов матрицы A. Тогда, матрица A может быть записана в виде спектрального разложения

A = S Λ S^{T} .

Найдём выражение для степеней исходной матрицы через матрицы Q_k и R_k. С одной стороны, по определению QR-алгоритма:

A^{k} = A_{1}^{k} = {(Q_{1} R_{1})}^{k} = Q_{1} {(R_{1} Q_{1})}^{k - 1} R_{1} = Q_{1} A_{2}^{k - 1} R_{1} .

Применяя это соотношение рекурсивно, получим:

A^{k} = Q_{1} \cdot . . . \cdot Q_{k} \cdot R_{k} \cdot . . . \cdot R_{1}

Введя следующие обозначения:

S_{k} = Q_{1} \cdot . . . \cdot Q_{k},

T_{k} = R_{k} \cdot . . . \cdot R_{1},

получим

A^{k} = S_{k} T_{k} .

С другой стороны:

A^{k} = S Λ^{k} S^{T} .

Приравнивая правые части последних двух формул, получим:

S Λ^{k} S^{T} = S_{k} T_{k} .

Предположим, что существует LU-разложение матрицы S^T:

S^{T} = L U,

тогда

S Λ^{k} L U = S_{k} T_{k} .

Умножим справа на обратную к U матрицу, а затем — на обратную к Λ^k:

S Λ^{k} L = S_{k} T_{k} U^{- 1},

S Λ^{k} L Λ^{- k} = S_{k} T_{k} U^{- 1} Λ^{- k} .

Можно показать, что

Λ^{k} L Λ^{- k} \to d i a g (l_{11}, . . ., l_{n n}) = L^{'} .

При $k \to \infty$ без ограничения общности можно считать, что на диагонали матрицы L стоят единицы, поэтому

S_{k} T_{k} U^{- 1} Λ^{- k} \to S .

Обозначим

P_{k} = T_{k} U^{- 1} Λ^{- k},

причём матрица P_k является верхнетреугольной, как произведение верхнетреугольных и диагональных матриц.

Таким образом, мы доказали, что

S_{k} P_{k} \to S

.

Из единственности QR-разложения следует, что если произведение ортогональной и треугольной матрицы сходится к ортогональной матрице, то треугольная матрица сходится к единичной матрице. Из сказанного следует, что

S_{k} = Q_{1} \cdot . . . \cdot Q_{k} \to S .

То есть матрицы S_k сходятся к матрице собственных векторов матрицы A.

Так как

A_{k + 1} = Q_{k}^{T} A_{k} Q_{k} = . . . = (Q_{k}^{T} \cdot . . . \cdot Q_{1}^{T}) A_{1} (Q_{1} \cdot . . . \cdot Q_{k}) = {(Q_{1} \cdot . . . \cdot Q_{k})}^{T} A (Q_{1} \cdot . . . \cdot Q_{k}),

то

A_{k + 1} = S_{k}^{T} A S_{k} .

Переходя к пределу, получим:

\lim_{k \to \infty} A_{k} = \lim_{k \to \infty} A_{k + 1} = S^{T} A S = S^{T} S Λ S^{T} S = Λ .

Итак, мы доказали, что QR-алгоритм позволяет решить полную проблему собственных значений для симметричной положительно-определённой матрицы.

Реализация QR-алгоритма[править]

При определенных условиях последовательность матриц $A_{k}$ сходится к треугольной матрице, разложению Шура матрицы $A$ . В этом случае собственные числа треугольной матрицы располагаются на ее диагонали, и задача нахождения собственных чисел считается решенной. В тестах на сходимость не практично требовать точных нулей в нулевой части матрицы, но можно воспользоваться теоремой о кругах Гершгорина, задающей пределы ошибок.

В исходном состоянии матрицы (без дополнительных преобразований) стоимость итераций относительно высока. Стоимость алгоритма можно уменьшить, сначала приведя матрицу $A$ к форме верхней матрицы Хессенберга (стоимость получения которой методом, основанном на преобразовании Хаусхолдера, оценивается как $\frac{10}{3} n^{3} + O (n^{2})$ арифметических операций), и затем используя конечную последовательность ортогональных преобразований подобия. Данный алгоритм чем-то похож на двухсторонюю QR-декомпозицию. (В обычной QR-декомпозиции матрица отражений Хаусхолдера умножается на исходную только слева, тогда как в случае использования формы Хессенберга матрица отражений умножается на исходную и слева, и справа.) Нахождение QR-декомпозиции верхней матрицы Хессенберга оценивается как $6 n^{2} + O (n)$ арифметических операций. Из-за того, что форма матрицы Хессенберга почти верхнетреугольная (у неё только один поддиагональный элемент не равен нулю), удается сразу снизить число итераций требуемых для схождения QR-алгоритма.

В случае, если исходная матрица симметричная, верхняя матрица Хессенберга также симметричная и поэтому является трехдиагональной. Этим же свойством обладает вся последовательность матриц $A_{k}$ . В этом случае стоимость процедуры оценивается как $\frac{4}{3} n^{3} + O (n^{2})$ арифметических операций с использованием метода, основанного на преобразовании Хаусхолдера. Нахождение QR-разложения симметричной трехдиагональной матрицы оценивается как $O (n)$ операций.

Скорость сходимости зависит от степени разделения собственных чисел, и в практической реализации в явном или неявном виде используюся «сдвиги» для усиления разделения собственных чисел и для ускорения сходимости. В типичном виде для симметричных матриц QR-алгоритм точно находит одно собственное число (уменьшая размерность матрицы) за одну или две итерации, что делает этот подход как эффективным, так и надежным.

Реализация QR-алгоритма в неявном виде[править]

В современной вычислительной практике QR-алгоритм реализуется с использованием его неявной версии, что упрощает добавление множественных «сдвигов». Исходно матрица приводится к форме верхней матрицы Хессенберга $A_{0} = Q A Q^{T}$ , так же как и в явной версии. Затем, на каждом шаге первая колонка $A_{k}$ преобразуется через малоразмерное преобразование подобия Хаусхолдера к первой колонке $p (A_{k})$ (или $p (A_{k}) e_{1}$ ), где $p (A_{k})$ — это полином степени $r$ , который определяет стратегию «сдвигов» (обычно $p (x) = (x - λ) (x - \bar{λ})$ , где $λ$ и $\bar{λ}$ — это два собственных числа остаточной подматрицы $A_{k}$ размера 2×2, это так называемый неявный двойной сдвиг). Затем последовательные преобразования Хаусхолдера размерности $r + 1$ производятся с целью вернуть рабочую матрицу $A_{k}$ к форме верхней матрицы Хессенберга.

Примечания[править]

↑ Численные методы / Н. С. Бахвалов, Н. П. Жидков, Г. М. Кобельков. — 3-е изд. — М.: БИНОМ, Лаборатория знаний, 2004. — С. 321. — 636 с. — ISBN 5-94774-175-X.

Ссылки[править]

Заметки Питера Олвера о ортогональных базисах и QR-алгоритме (англ.)

[1] Численные методы / Н. С. Бахвалов, Н. П. Жидков, Г. М. Кобельков. — 3-е изд. — М.: БИНОМ, Лаборатория знаний, 2004. — С. 321. — 636 с. — ISBN 5-94774-175-X.

[1]