User:GrishinaAnna/Метод Якоби вычисления сингулярных чисел и векторов
1 Свойства и структура алгоритма
1.1 Общее описание алгоритма
Метод Якоби позволяет получить разложение произвольной матрицы [math]G[/math] размера [math]n[/math]×[math]n[/math] в виде [math]G=USV^*[/math], где [math]U[/math] и [math]V[/math] — унитарные матрицы, [math]S[/math] — диагональная матрица с вещественными положительными числами на диагонали. Диагональные элементы матрицы называются сингулярными числами матрицы [math]G[/math], а столбцы матриц и левыми и правыми сингулярными векторами соответственно.
В основе метода лежит утверждение о том, что собственными значениями симметричной матрицы [math]G^TG[/math] являются квадраты сингулярных чисел, а соответствующие ортонормированные собственные векторы – правыми сингулярными векторами. Для вычисления собственных значений и собственных векторов симметричной матрицы [math]G^TG[/math] используется метод односторонних вращений (неявный метод Якоби), на каждом шаге которого вычисляется вращение Якоби [math]J[/math]. С его помощью матрица [math]G^TG[/math] неявно пересчитывается в [math]J^TG^TGJ[/math]. Вращение выбрано так, чтобы пара внедиагональных элементов из [math]G^TG[/math] обратилась в нуль в матрице [math]J^TG^TGJ[/math]. Однако, ни [math]G^TG[/math], ни [math]J^TG^TGJ[/math] не вычисляются в явном виде; вместо них вычисляется матрица [math]GJ[/math], которая подается на вход на следующем шаге алгоритма. Алгоритм выполняется, пока [math]G^TG[/math] не станет достаточно близка к диагональной матрице, после чего вычисляются собственные числа матрицы [math]G[/math], полученной на последнем шаге, которые принимаются за сингулярные числа. Правые сингулярные векторы получаются в результате накопления произведений вращений Якоби.
Метод Якоби для вычисления сингулярных значений и векторов является самым медленным из имеющихся, но тем не менее, интерес к нему сохраняется, потому что для некоторых типов матриц [math]G[/math] он способен вычислять числа и сингулярные векторы намного точнее, чем другие методы. Так, например, метод Якоби вычисляет сингулярные числа матрицы [math]G[/math] с высокой точностью, если [math]G[/math] может быть представлена в виде [math]G=DX[/math], где [math]D[/math] – диагональная матрица, а [math]X[/math] – хорошо обусловлена. Выигрыш метода Якоби в этом случае перед другими алгоритмами объясняется следующим образом: в методе Якоби заданная матрица обрабатывается без предварительного приведения к двухдиагональному виду, в то время как другие алгоритмы включают в себя приведение матрицы к двухдиагональной форме, из-за чего и теряют все верные разряды во всех сингулярных числах, кроме старшего.
Следует иметь в виду, что метод Якоби эффективен при не очень больших размерах матрицы – максимум до [math]n[/math] ~ [math]100[/math]. Дело не только в том, что метод не будет сходиться за конечное время, дело еще и в излишне большом количестве операций, которое приведет к ошибкам округления и потере точности.
1.2 Математическое описание алгоритма
Исходные данные: матрица [math]G[/math] размера [math]n[/math]×[math]n[/math].
Вычисляемые данные: Матрица [math]S=[/math]diag([math]\sigma_{i}[/math]), где [math]\sigma_{i}[/math] - сингулярные числа, матрица [math]U[/math] левых сингулярных векторов и матрица [math]V[/math] правых сингулярных векторов.
Описание используемой в алгоритме процедуры одностороннего вращения входной матрицы [math]G[/math] One-Sided-Jacobi-Rotation [math](G, j, k)[/math]:
proc One-Sided-Jacobi-Rotation [math](G,j,k)[/math] вычислить [math]a_{jj} = (G^TG)_{jj}, a_{jk} = (G^TG)_{jk}[/math] и [math]a_{kk} = (G^TG)_{kk}[/math] if [math]|a_{jk}|[/math] не слишком мал [math] \begin{align} \tau &= {a_{jj} - a_{kk}}/{2a_{jk}} \\ t &= {sign(\tau)}/(|\tau| + \sqrt{1 + \tau^2}) \\ c &= 1/(\sqrt{1 + t^2}) \\ s &= ct \\ G &= GR(j,k,\theta) \ \dots c=cos\theta, s=sin\theta \end{align}[/math] if нужны правые сингулярные векторы [math]J = JR(j,k,\theta)[/math] end if end if
Где [math]R(j,k,\theta)[/math] – матрица вращения Якоби, которая выбирается так, чтобы обнулить пару внедиагональных элементов.
[math]R(j,k,\theta) =\begin{pmatrix} 1 & 0 & \cdots & 0\\ 0 & 1 & \cdots & 0\\ \vdots & \vdots & \ddots & \vdots\\ \0 & \cos (\theta) & \sin(\theta) & 0 \end{pmatrix} [/math]