Участница:Ekaterina.ivkina/Метод Якоби вычисления сингулярных чисел и векторов: различия между версиями
Victoria (обсуждение | вклад) |
Victoria (обсуждение | вклад) |
||
Строка 28: | Строка 28: | ||
* <math>U</math> - матрица левых сингулярных векторов, | * <math>U</math> - матрица левых сингулярных векторов, | ||
* <math>V</math> - матрица правых сингулярных векторов. | * <math>V</math> - матрица правых сингулярных векторов. | ||
+ | |||
+ | '''Описание алгоритма:''' | ||
+ | |||
+ | '''Шаг 1''' Вычисляем матрицу <math>G'</math> применением одностороннего вращения Якоби к исходной матрице <math>G</math>. | ||
+ | |||
+ | '''Шаг 2''' Повторяем Шаг 1 пока матрица <math>G'^TG</math> не станет достаточно близка к диагональной. | ||
+ | |||
+ | '''Шаг 3''' Вычисляем по полученной матрице <math>G'</math> сингулярные числа <math>\sigma_i</math> и сингулярные векторы <math>U, V</math>. | ||
'''Основные формулы метода:''' | '''Основные формулы метода:''' | ||
Строка 37: | Строка 45: | ||
<math>V = J</math>, где <math>J</math> - накопленное произведение вращений Якоби. | <math>V = J</math>, где <math>J</math> - накопленное произведение вращений Якоби. | ||
− | + | '''Основные формулы процедуры одностороннего вращения Якоби:''' | |
− | |||
− | Основные формулы процедуры одностороннего вращения Якоби: | ||
Пусть <math>a_{ij}</math> - элемент, расположенный в <math>i</math>-ой строке и <math>j</math>-ом столбце матрицы <math>G^TG</math>. | Пусть <math>a_{ij}</math> - элемент, расположенный в <math>i</math>-ой строке и <math>j</math>-ом столбце матрицы <math>G^TG</math>. | ||
− | <math>\tau = (a_{jj} - a_{kk})/(2 \ | + | <math>\tau = (a_{jj} - a_{kk})/(2 \cdot a_{jk})</math> |
<math>t = sign(\tau)/(|\tau|+\sqrt{1+\tau^2})</math> | <math>t = sign(\tau)/(|\tau|+\sqrt{1+\tau^2})</math> | ||
− | <math>c = 1/\sqrt{1+t^2}</math> | + | <math>c = 1/\sqrt{1+t^2}</math>, где <math>c = cos\theta</math> |
+ | |||
+ | <math>s = c \cdot t</math>, где <math>s = sin\theta</math> | ||
+ | |||
+ | <math>J = J \cdot R(j,k,\theta)</math> | ||
+ | Здесь <math>R(j,k,\theta)</math> - матрица вращений Якоби, которая имеет следующий вид: | ||
− | + | | |
− | + | <math> | |
+ | \begin{matrix} | ||
+ | & & & & j & & k & & & \\ | ||
+ | \end{matrix} | ||
+ | </math> | ||
+ | <math> | ||
+ | R(j,k,\theta) = | ||
+ | \begin{matrix} | ||
+ | \\ | ||
+ | \\ | ||
+ | \\ | ||
+ | j \\ | ||
+ | \\ | ||
+ | k \\ | ||
+ | \\ | ||
+ | \\ | ||
+ | \\ | ||
+ | \end{matrix} | ||
+ | </math> | ||
+ | <math> | ||
+ | \begin{pmatrix} | ||
+ | & 1 & & & & & & & & \\ | ||
+ | & & 1 & & & & & & & \\ | ||
+ | & & & \ddots & & & & & & \\ | ||
+ | & & & & \cos(\theta) & & -\sin(\theta) & & & \\ | ||
+ | & & & & & \ddots & & & & \\ | ||
+ | & & & & \sin(\theta) & & \cos(\theta) & & & \\ | ||
+ | & & & & & & & \ddots & & \\ | ||
+ | & & & & & & & & 1 & \\ | ||
+ | & & & & & & & & & 1 \\ | ||
+ | \end{pmatrix} | ||
+ | </math> | ||
+ | |||
+ | |||
+ | Одностороннее вращение Якоби в координатной плоскости <math>i, j</math> вычисляется в том случае, если элемент <math>a_{ij}</math> удовлетворяет условию: | ||
+ | <math>|a_{ij}| \geq \epsilon \sqrt{a_{ii}a_{jj}}</math>, где <math>\epsilon</math> - | ||
+ | |||
+ | '''Критерий останова:''' | ||
+ | |||
+ | Для проверки близости матрицы <math>G'^TG</math> к диагональной используется следующее неравенство: | ||
=== Вычислительное ядро алгоритма === | === Вычислительное ядро алгоритма === | ||
+ | Вычислительное ядро метода Якоби вычисления сингулярных чисел и векторов можно составить из множественных вызовов процедуры одностороннего вращения Якоби '''''One-Sided-Jacobi-Rotation(G,j,k)'''''. | ||
+ | Основными операциями процедуры '''''One-Sided-Jacobi-Rotation(G,j,k)''''' являются: | ||
+ | * Вычисление элементов <math>a_{jj}, a_{jk}, a_{kk}</math>, где <math>a_{ij} = (G^TG)_{ij}</math>. | ||
+ | * Вычисление матрицы вращения Якоби <math>R(j,k,\theta)</math>. | ||
+ | * Применение вращения Якоби к матрицам <math>G</math> и <math>J</math>. | ||
+ | === Макроструктура алгоритма === | ||
− | + | Алгоритм вычисления сингулярных чисел и векторов Якоби состоит из множественных обращений к процедуре одностороннего вращения Якоби. На каждой итерации цикла происходит вызов процедуры '''''One-Sided-Jacobi-Rotation(G, j, k)''''' для всех <math>j \in [1, n-1]</math> и <math>k \in [j+1, n]</math>. Таким образом, на каждой итерации цикла процедура '''''One-Sided-Jacobi-Rotation(G, j, k)''''' вызывается для всех наддиагональных элементом матрицы <math>G^TG</math>, то есть всего <math>\frac{n(n-1)}{2}</math> раз. Количество итераций цикла непосредственно зависит от матрицы <math>G</math>, выход из цикла происходит по достижении критерия останова (описан в разделе ''"Математическое описание алгоритма"''.) | |
+ | === Схема реализации последовательного алгоритма === | ||
+ | Схему реализации последовательного алгоритма можно описать следующиv псевдокодом: | ||
− | === | + | ''repeat |
+ | ''for <math>j=1</math> to <math>n-1</math> | ||
+ | ''for <math>k=j+1</math> to <math>n</math> | ||
+ | '''''One-Sided-Jacobi-Rotation(G, j, k) | ||
+ | ''end for | ||
+ | ''end for | ||
+ | ''пока <math>G^TG</math> не станет достаточно близка к диагональной матрице. | ||
+ | ''Положить <math>\sigma_i = \parallel G'(:,i)\parallel_2</math> (2-я норма <math>i</math>-го столбца в <math>G'</math>). | ||
+ | ''Положить <math>U = [u_1, \dots, u_n]</math>, где <math> u_i = G'(:,i)/\sigma_i</math>. | ||
+ | ''Положить <math>V = J</math>, где <math>J</math> - накопленное произведение вращений Якоби. | ||
+ | Процедура '''''One-Sided-Jacobi-Rotation(G, j, k)''''': | ||
+ | ''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>\tau = {a_{jj} - a_{kk}}/{2a_{jk}}</math> | ||
+ | <math>t = {sign(\tau)}/(|\tau| + \sqrt{1 + \tau^2})</math> | ||
+ | <math>c = 1/(\sqrt{1 + t^2})</math> | ||
+ | <math>s = ct</math> | ||
+ | <math>G = GR(j,k,\theta) \dots </math>''где <math>c=cos\theta, s=sin\theta</math> | ||
+ | ''if нужны правые сингулярные векторы | ||
+ | <math>J = JR(j,k,\theta)</math> | ||
+ | ''end if | ||
+ | ''end if | ||
=== Последовательная сложность алгоритма === | === Последовательная сложность алгоритма === |
Версия 20:48, 20 октября 2016
Метод Якоби вычисления сингулярных чисел и векторов | |
Последовательный алгоритм | |
Последовательная сложность | [math][/math] |
Объём входных данных | [math][/math] |
Объём выходных данных | [math][/math] |
Параллельный алгоритм | |
Высота ярусно-параллельной формы | [math][/math] |
Ширина ярусно-параллельной формы | [math][/math] |
Содержание
- 1 Свойства и структура алгоритма
- 1.1 Общее описание алгоритма
- 1.2 Математическое описание алгоритма
- 1.3 Вычислительное ядро алгоритма
- 1.4 Макроструктура алгоритма
- 1.5 Схема реализации последовательного алгоритма
- 1.6 Последовательная сложность алгоритма
- 1.7 Информационный граф
- 1.8 Ресурс параллелизма алгоритма
- 1.9 Входные и выходные данные алгоритма
- 1.10 Свойства алгоритма
- 2 Программная реализация алгоритма
- 2.1 Особенности реализации последовательного алгоритма
- 2.2 Локальность данных и вычислений
- 2.3 Возможные способы и особенности параллельной реализации алгоритма
- 2.4 Масштабируемость алгоритма и его реализации
- 2.5 Динамические характеристики и эффективность реализации алгоритма
- 2.6 Выводы для классов архитектур
- 2.7 Существующие реализации алгоритма
- 3 Литература
1 Свойства и структура алгоритма
1.1 Общее описание алгоритма
Сингулярное разложение было первоначально разработано в дифференциальной геометрии при изучении свойств билинейных форм учеными Эудженио Бельтрами и Камилем Жорданом независимо в 1873 и 1874 годах соответственно. Первое доказательство сингулярного разложения для прямоугольных и комплексных матриц было осуществлено математиками Карлом Эскартом и Гэйлом Янгом в 1936 году.
Сингулярное разложение (Singular Values Decomposition, SVD) является удобным методом при работе с матрицами. Cингулярное разложение показывает геометрическую структуру матрицы и позволяет наглядно представить имеющиеся данные. Сингулярное разложение используется при решении самых разных задач — от приближения методом наименьших квадратов и решения систем уравнений до сжатия и распознавания изображений. Используются разные свойства сингулярного разложения, например, способность показывать ранг матрицы и приближать матрицы данного ранга. Так как вычисления ранга матрицы — задача, которая встречается очень часто, то сингулярное разложение является довольно популярным методом.
Сингулярным разложением матрицы [math]G (n \times n)[/math] называется разложение вида [math]G = U \Sigma V^T[/math] , где [math]U, V[/math] - унитарные матрицы [math]n \times n[/math], [math]\Sigma[/math] - диагональная матрица [math]n \times n[/math] с вещественными положительными числами на главной диагонали. Столбцы матриц [math]U, V[/math] называются соответственно левыми и правыми сингулярными векторами, а значения на диагонали матрицы [math]\Sigma[/math] - сингулярными значениями матрицы [math]G[/math].
Одним из методов нахождения сингулярного разложения матрицы является метод Якоби. Метод Якоби был предложен Карлом Густавом Якоби Якоби в 1846 году и представляет собой итерационный алгоритм вычисления собственных значений и собственных векторов симметричной матрицы. Для вычисления собственных значений необходимо неявно применить метод Якоби к симметричной матрице [math]A = G^T G[/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]. Поэтому алгоритм называется методом односторонних вращений.
1.2 Математическое описание алгоритма
Исходные данные: Матрица [math]G[/math] размера [math](n \times n)[/math] над полем вещественных или комплексных чисел.
Вычисляемые данные:
- [math]\sigma_i[/math] - сингулярные числа матрицы [math]G[/math],
- [math]U[/math] - матрица левых сингулярных векторов,
- [math]V[/math] - матрица правых сингулярных векторов.
Описание алгоритма:
Шаг 1 Вычисляем матрицу [math]G'[/math] применением одностороннего вращения Якоби к исходной матрице [math]G[/math].
Шаг 2 Повторяем Шаг 1 пока матрица [math]G'^TG[/math] не станет достаточно близка к диагональной.
Шаг 3 Вычисляем по полученной матрице [math]G'[/math] сингулярные числа [math]\sigma_i[/math] и сингулярные векторы [math]U, V[/math].
Основные формулы метода:
[math]\sigma_i = \parallel G'(:,i)\parallel_2[/math] (2-я норма [math]i[/math]-го столбца в [math]G'[/math]).
[math]U = [u_1, \dots, u_n][/math], где [math] u_i = G'(:,i)/\sigma_i[/math].
[math]V = J[/math], где [math]J[/math] - накопленное произведение вращений Якоби.
Основные формулы процедуры одностороннего вращения Якоби:
Пусть [math]a_{ij}[/math] - элемент, расположенный в [math]i[/math]-ой строке и [math]j[/math]-ом столбце матрицы [math]G^TG[/math].
[math]\tau = (a_{jj} - a_{kk})/(2 \cdot a_{jk})[/math]
[math]t = sign(\tau)/(|\tau|+\sqrt{1+\tau^2})[/math]
[math]c = 1/\sqrt{1+t^2}[/math], где [math]c = cos\theta[/math]
[math]s = c \cdot t[/math], где [math]s = sin\theta[/math]
[math]J = J \cdot R(j,k,\theta)[/math]
Здесь [math]R(j,k,\theta)[/math] - матрица вращений Якоби, которая имеет следующий вид:
[math] \begin{matrix} & & & & j & & k & & & \\ \end{matrix} [/math]
[math] R(j,k,\theta) = \begin{matrix} \\ \\ \\ j \\ \\ k \\ \\ \\ \\ \end{matrix} [/math] [math] \begin{pmatrix} & 1 & & & & & & & & \\ & & 1 & & & & & & & \\ & & & \ddots & & & & & & \\ & & & & \cos(\theta) & & -\sin(\theta) & & & \\ & & & & & \ddots & & & & \\ & & & & \sin(\theta) & & \cos(\theta) & & & \\ & & & & & & & \ddots & & \\ & & & & & & & & 1 & \\ & & & & & & & & & 1 \\ \end{pmatrix} [/math]
Одностороннее вращение Якоби в координатной плоскости [math]i, j[/math] вычисляется в том случае, если элемент [math]a_{ij}[/math] удовлетворяет условию:
[math]|a_{ij}| \geq \epsilon \sqrt{a_{ii}a_{jj}}[/math], где [math]\epsilon[/math] -
Критерий останова:
Для проверки близости матрицы [math]G'^TG[/math] к диагональной используется следующее неравенство:
1.3 Вычислительное ядро алгоритма
Вычислительное ядро метода Якоби вычисления сингулярных чисел и векторов можно составить из множественных вызовов процедуры одностороннего вращения Якоби One-Sided-Jacobi-Rotation(G,j,k).
Основными операциями процедуры One-Sided-Jacobi-Rotation(G,j,k) являются:
- Вычисление элементов [math]a_{jj}, a_{jk}, a_{kk}[/math], где [math]a_{ij} = (G^TG)_{ij}[/math].
- Вычисление матрицы вращения Якоби [math]R(j,k,\theta)[/math].
- Применение вращения Якоби к матрицам [math]G[/math] и [math]J[/math].
1.4 Макроструктура алгоритма
Алгоритм вычисления сингулярных чисел и векторов Якоби состоит из множественных обращений к процедуре одностороннего вращения Якоби. На каждой итерации цикла происходит вызов процедуры One-Sided-Jacobi-Rotation(G, j, k) для всех [math]j \in [1, n-1][/math] и [math]k \in [j+1, n][/math]. Таким образом, на каждой итерации цикла процедура One-Sided-Jacobi-Rotation(G, j, k) вызывается для всех наддиагональных элементом матрицы [math]G^TG[/math], то есть всего [math]\frac{n(n-1)}{2}[/math] раз. Количество итераций цикла непосредственно зависит от матрицы [math]G[/math], выход из цикла происходит по достижении критерия останова (описан в разделе "Математическое описание алгоритма".)
1.5 Схема реализации последовательного алгоритма
Схему реализации последовательного алгоритма можно описать следующиv псевдокодом:
repeat for [math]j=1[/math] to [math]n-1[/math] for [math]k=j+1[/math] to [math]n[/math] One-Sided-Jacobi-Rotation(G, j, k) end for end for пока [math]G^TG[/math] не станет достаточно близка к диагональной матрице. Положить [math]\sigma_i = \parallel G'(:,i)\parallel_2[/math] (2-я норма [math]i[/math]-го столбца в [math]G'[/math]). Положить [math]U = [u_1, \dots, u_n][/math], где [math] u_i = G'(:,i)/\sigma_i[/math]. Положить [math]V = J[/math], где [math]J[/math] - накопленное произведение вращений Якоби.
Процедура One-Sided-Jacobi-Rotation(G, j, k):
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]\tau = {a_{jj} - a_{kk}}/{2a_{jk}}[/math] [math]t = {sign(\tau)}/(|\tau| + \sqrt{1 + \tau^2})[/math] [math]c = 1/(\sqrt{1 + t^2})[/math] [math]s = ct[/math] [math]G = GR(j,k,\theta) \dots [/math]где [math]c=cos\theta, s=sin\theta[/math] if нужны правые сингулярные векторы [math]J = JR(j,k,\theta)[/math] end if end if
1.6 Последовательная сложность алгоритма
1.7 Информационный граф
1.8 Ресурс параллелизма алгоритма
1.9 Входные и выходные данные алгоритма
1.10 Свойства алгоритма
2 Программная реализация алгоритма
2.1 Особенности реализации последовательного алгоритма
2.2 Локальность данных и вычислений
2.3 Возможные способы и особенности параллельной реализации алгоритма
2.4 Масштабируемость алгоритма и его реализации
2.5 Динамические характеристики и эффективность реализации алгоритма
2.6 Выводы для классов архитектур
2.7 Существующие реализации алгоритма
3 Литература
<references \>