Участник:VolkovNikita94/Алгоритм Ланцоша для точной арифметики (без переортогонализации): различия между версиями

Материал из Алговики
Перейти к навигации Перейти к поиску
Строка 18: Строка 18:
  
 
Вторая часть алгоритма представляет собой поиск собственных значений и собственных векторов симметричной трёхдиагональной матрицы <math>T_{j}</math>, что реализуется отдельным алгоритмом ( QR-алгоритм, метод «Разделяй-и-властвуй», итерации Арнольди и.т.д. ), который мы будем рассматривать как отдельную макрооперацию.
 
Вторая часть алгоритма представляет собой поиск собственных значений и собственных векторов симметричной трёхдиагональной матрицы <math>T_{j}</math>, что реализуется отдельным алгоритмом ( QR-алгоритм, метод «Разделяй-и-властвуй», итерации Арнольди и.т.д. ), который мы будем рассматривать как отдельную макрооперацию.
 
===Схема реализации последовательного алгоритма===
 
Каждая итерация алгоритма Ланцоща выполняется строго последовательно. Однако внутри итерации возможно распараллеливание алгоритма.
 
 
*Умножение квадратной матрицы порядка <math>n</math> на вектор длины <math>n</math> требует последовательного выполнения <math>n</math> ярусов умножений и сложений.
 
*Все остальные операции в рамках одной итерации, ведущие к вычислению матрицы <math>Tj</math> выполняются строго последовательно, но каждая из них требует последовательного выполнения только одного яруса операций ( умножений, сложений ).
 
*Ресурс параллелизма алгоритма вычисления собственных значений зависит от используемого алгоритма и рассматривается в соответствующей статье.
 
 
Таким образом, при классисифкации по ширине ЯПФ алгоритм обладает <i>линейной</i> сложностью. В случае классификации по высоте ЯПФ алгоритм также имеет <i>линейную</i> сложность.
 
  
 
===Последовательная сложность алгоритма===
 
===Последовательная сложность алгоритма===

Версия 15:41, 15 октября 2016

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

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

Метод Ланцоша — процедура, применяемая для нахождения [math]k[/math] собственных значений симметричной матрицы произвольного размера путём вычисления их приближений. Суть метода заключается в итеративном построении специальной матрицы [math]T_{j}[/math], а затем вычислении её собственных значений и собственных векторов.

Здесь мы рассмотрим упрощенный вариант алгоритма Ланцоша, который подразумевает отсутствие влияния ошибок округления на вычислительный процесс. На практике используются модификации алгоритма Ланцоша, обладающие устойчивостью и без этого допущения ( например, алгоритм Ланцоша с полной переортогонализацией ).

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

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

Вычислительное ядро алгоритма Ланцоща состоит в получении на каждой итерации очередного промежуточного вектора [math]z[/math], получаемого путём умножения исходной матрицы [math]A[/math] на вектор Ланцоша [math]q_{j}[/math], полученный на предыдущей итерации.

Также при значениях [math]k[/math], сопоставимых с [math]n[/math], процедура вычисления собственных значений и собственных векторов симметричной трёхдиагональной матрицы по некоторому выбранному алгоритму подразумевает значительный объём расчётов и может считаться вычислительным ядром. Однако на практике метод используется прежде всего при малых [math]k[/math].

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

Метод Ланцоша соединяет метод Ланцоша для построения крыловского подпространства в процедурой Рэлея-Ритца, подразумевающей вычисление собственных значений симметричной трёхдиагональной матрицы.

Первая часть алгоритма представляет собой строго последовательные итерации, на каждой из которых вначале строится очередной столбец матрицы [math]Q_{j}[/math] из первых [math]j[/math] ортонормированных векторов Ланцоша и матрица [math]T_{j} = Q^T_{j} * A * Q_{j}[/math] порядка j. Итерационный процесс завершается, когда [math]j = k[/math]

Вторая часть алгоритма представляет собой поиск собственных значений и собственных векторов симметричной трёхдиагональной матрицы [math]T_{j}[/math], что реализуется отдельным алгоритмом ( QR-алгоритм, метод «Разделяй-и-властвуй», итерации Арнольди и.т.д. ), который мы будем рассматривать как отдельную макрооперацию.

1.5 Последовательная сложность алгоритма

Все дальнейшие выкладки верны для наиболее быстрого последовательного варианта выполнения указываемых операций.

  • Умножение квадратной матрицы порядка [math]n[/math] на вектор длины [math]n[/math] и наоборот требует по [math]n^2[/math] умножений и сложений.
  • Перемножение векторов длины [math]n[/math] требует по [math]n[/math] умножений и сложений.
  • Поэлементное сложение векторов длины [math]n[/math] требует [math]n[/math] сложений.
  • Умножение вектора длины [math]n[/math] на число требует [math]n[/math] умножений.
  • Нахождение квадратичной нормы вектора длины [math]n[/math] требует по [math]n[/math] умножений и сложений, а также одну операцию извлечения квадратного корня.
  • Вычисление собственных значений матрицы порядка [math]k[/math] QR-алгоритмом требует [math]O(k^2)[/math] операций, вычисление также и собственных векторов требует примерно [math]6 * k ^ 3[/math], то есть [math]O(k^3)[/math] операций; при использовании метода «Разделяй-и-властвуй» для вычисления собственных значений и векторов аналогичной матрицы в среднем затрачивается [math]O(k^{2.3})[/math] операций.

Таким образом, для выполнения одной итерации метода Ланцоща требуется 1 операция вычисления квадратного корня, [math]n[/math] умножений, а также по [math]n ^ 2 + n + 2 * n + n[/math] сложений и умножений, а также суммарно [math]O(k ^ 3)[/math] операций, требуемых для поиска собственных значений и собственных векторов матрицы [math]T_{j}[/math]. То есть последовательная сложность алгоритма Ланцоша составляет [math]O(n ^ 2) + O(k ^ 3)[/math]. Это, очевидно, меньше, чем [math]O(n ^ 3)[/math] операций, требуемых в алгоритмах вычисления всех собственных значений произвольных симметричных матриц, в чём и заключается выгода метода Ланцоша.

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

1.7 Ресурс параллелизма алгоритма

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

Входные данные: квадратная матрица [math]A[/math] ( элементы [math]a_{ij}[/math] ), вектор [math]b[/math], число [math]k[/math]
Объём входных данных: [math]n^2 + n + 1[/math]
Выходные данные: вектор [math]L[/math], матрица [math]E[/math] ( элементы [math]e_{pq}[/math] )
Объём выходных данных: [math]k + k * n[/math]

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

В случае неограниченности ресурсов соотношение последовательной и параллельной сложности алгоритма Ланцоща по сумме всех итераций представляет собой [math]k ^ 3 + k * n ^ 2[/math] к [math]X + k * n[/math], где [math]X[/math] - параллельная сложность алгоритма, используемого для вычисления собственных значений. В случае использования метода вычисления собственных значений с линейной сложностью при классисифкации по высоте ЯПФ это соотношение будет равно [math](k ^ 2 + n ^ 2) / (k + n)[/math].

Алгоритм Ланцоша не является полностью детерминированным в том смысле, что возможно выполнение числа итераций алгоритма, меньшего, чем заданное ( из-за того, что вычислены все ненулевые собственные значения матрицы [math]A[/math] ).

Важное свойство метода Ланцоша состоит в том, что первыми в матрице [math]T_{j}[/math] появляются собственные значения с максимальной величиной по модулю. Таким образом, метод особенно хорошо подходит для вычисления собственных значений матрицы [math]A[/math], находящихся на краях её спектра.

2 Программная реализация алгоритма

2.1 Особенности реализации последователього алгоритма

2.2 Локальность данных и вычислений

2.3 Возможные способы и особенности параллельной реализации алгоритма

2.4 Масштабируемость алгоритма и его реализации

Для проверки масштабируемости алгоритма из реализации была исключена непосредственно задача поиска собственных значений матрицы [math]T_{j}[/math] на каждой итерации, так как это является предметом отдельных статей ( QR-алгоритм, метод «Разделяй-и-властвуй», итерации Арнольди и.т.д. ). Таким образом, проверяемая на масштабируемость задача сводится к итерационному вычислению коэффициентов матрицы [math]T_{j}[/math]. Матрица [math]A[/math] в рассматриваемой реализации хранится построчно и построчно же разделяется между процессорами.

2.5 Динамические характеристики и эффективность реализации алгоритма

2.6 Выводы для классов архитектур

2.7 Существующие реализации алгоритма

Существует несколько как последовательных, так и параллельных реализаций алгоритма Ланцоша, включенных в программные библиотеки для поиска собственных значений матриц. Свойства этих реализаций можно увидеть в таблицах ниже. Первая таблица — относительно недавние реализации, вторая - «музейные экспонаты». Следует уделить особенное внимание столбцам «тип метода» и «параллелизация». В первом из них значение только N соответствует описанному варианту без реортогонализации, F – полной реортогонализации, P – частичной реортогонализации, S – выборочной реортогонализации. Во втором значние «none» соответствует отсутствию параллельной реализации, M – параллельной реализации посредством MPI, O – параллельной реализации посредством OpenMP. Их приведение целесообразно, так как на практике алгоритм Ланцоша без переортогонализации неустойчив. Также указываются библиотеки, в которых реализованы более глубокие модификации метода Ланцоша, с указанием изменений в графе «тип метода».

Таблица 1 - «современные» реализации.
Название библиотеки Язык программирования Дата появления, версия библиотки Тип метода Параллелизация
BLKLAN C/Matlab 2003 P none
BLZPACK F77 2000, 04/00 P + S M
IETL C++ 2006, 2.2 N none
SLEPc C/F77 2009, 3.0.0 All M
TRLAN F90 2006 Dynamic thick-restart M
PROPACK F77/Matlab 2005, 2.1 / 1.1 P, finds SVD O
IRBLEIGS Matlab 2000, 1.0 Indefinitie symmetric none
Таблица 2 - «исторические» реализации.
Название библиотеки Язык программирования Дата появления, версия библиотки Тип метода Параллелизация
ARPACK F77 1995, 2.1 Arnoldi iterations, impliicit restart M
ARPACK++ C++ 1998, 1.1 Arnoldi iterations, impliicit restart none
LANCZOS F77 1992 N none
LANZ F77 1991, 1.0 P none
LASO F77 1983, 2 S none
NAPACK F77 1987 N none
QMRPACK F77 1996 Designed for nonsymmetric matrices ( lookahead ) none
SVDPACK C/F77 1992 P, finds SVD none
Underwood F77 1975 F, block version none