Участник:Yastrebovks/Алгоритм Ланцоша с полной переортогонализацией: различия между версиями
| Строка 9: | Строка 9: | ||
Авторы: Алексейчук Н.Н 616, Ястребов К.С. 609. | Авторы: Алексейчук Н.Н 616, Ястребов К.С. 609. | ||
| − | == Общее описание алгоритма. == | + | ==Свойства и структура алгоритма == |
| + | === Общее описание алгоритма. === | ||
Алгоритм Ланцоша ищет собвственные значения и собственные векторы для симетричной матрицы A вещественных чисел. Является итерацонным алгоритмом. Алгоритм Ланцоша использует степенной метод (<math> b_{k+1} = \frac{Ab_k}{||Ab_k||} </math>) для поиска наибольших собственных значений и векторов матриц. | Алгоритм Ланцоша ищет собвственные значения и собственные векторы для симетричной матрицы A вещественных чисел. Является итерацонным алгоритмом. Алгоритм Ланцоша использует степенной метод (<math> b_{k+1} = \frac{Ab_k}{||Ab_k||} </math>) для поиска наибольших собственных значений и векторов матриц. | ||
| Строка 18: | Строка 19: | ||
| − | == Математическое описание алгоритма == | + | === Математическое описание алгоритма === |
Памятка: | Памятка: | ||
| Строка 34: | Строка 35: | ||
Из-за ошибок округления вектора <math> A^{k-1}x_1 </math>, формирующие подпространство Крылова, становятся неортогональными. Чтобы решить данную прорблему, проводят переортогонализацию методом Грамма-Шмидта. | Из-за ошибок округления вектора <math> A^{k-1}x_1 </math>, формирующие подпространство Крылова, становятся неортогональными. Чтобы решить данную прорблему, проводят переортогонализацию методом Грамма-Шмидта. | ||
| − | == Вычислительное ядро алгоритма == | + | === Вычислительное ядро алгоритма === |
Вычислительным ядром данного алгоритма являются следующие шаги: | Вычислительным ядром данного алгоритма являются следующие шаги: | ||
| Строка 43: | Строка 44: | ||
В данном случае первая операция выполняет умножение симметричной матрицы <math>A</math> размерности <math>n</math>X<math>n</math> на <math>n</math>-мерный вектор <math>q</math>, вследствие чего вычислительная сложность выполнения заключается в <math>n^2</math> умножений и <math>n^2-n</math> сложений. Второе действие является процессом ортогонализации Грама-Шмидта. В последнем действии вычисляются <math>k^2n+k(n+2)</math> умножений и <math>k^2n + k(n + 1) + 2</math> операций сложения. | В данном случае первая операция выполняет умножение симметричной матрицы <math>A</math> размерности <math>n</math>X<math>n</math> на <math>n</math>-мерный вектор <math>q</math>, вследствие чего вычислительная сложность выполнения заключается в <math>n^2</math> умножений и <math>n^2-n</math> сложений. Второе действие является процессом ортогонализации Грама-Шмидта. В последнем действии вычисляются <math>k^2n+k(n+2)</math> умножений и <math>k^2n + k(n + 1) + 2</math> операций сложения. | ||
| − | == Макроструктура алгоритма == | + | === Макроструктура алгоритма === |
*1. Скалярное произведение: <math> (x,y) </math>. | *1. Скалярное произведение: <math> (x,y) </math>. | ||
| Строка 49: | Строка 50: | ||
| − | == Схема реализации последовательного алгоритма == | + | === Схема реализации последовательного алгоритма === |
Input: A, b (random vector with unit norm) | Input: A, b (random vector with unit norm) | ||
: <math> | : <math> | ||
| Строка 68: | Строка 69: | ||
</math> | </math> | ||
| − | == Последовательная сложность алгоритма== | + | === Последовательная сложность алгоритма === |
Количество операций складывается из количества операций для классического метода Ланцоша[https://books.google.ru/books?id=T8EPBwAAQBAJ&pg=PA242&lpg=PA242&dq=lanczos+algorithm+complexity&source=bl&ots=kWVnVP_uJp&sig=a-4GTEt6vRE8YhR1R6dOPwvOSIk&hl=en&sa=X&ved=0ahUKEwiQ2u3U4-LPAhWvbZoKHe0QBaw4ChDoAQg7MAY#v=onepage&q=lanczos%20algorithm%20complexity&f=false] и количества операций для переортогонализации методом грамма-Шмидта[http://stackoverflow.com/questions/27986225/computational-complexity-of-gram-schmidt-orthogonalization-algorithm]. | Количество операций складывается из количества операций для классического метода Ланцоша[https://books.google.ru/books?id=T8EPBwAAQBAJ&pg=PA242&lpg=PA242&dq=lanczos+algorithm+complexity&source=bl&ots=kWVnVP_uJp&sig=a-4GTEt6vRE8YhR1R6dOPwvOSIk&hl=en&sa=X&ved=0ahUKEwiQ2u3U4-LPAhWvbZoKHe0QBaw4ChDoAQg7MAY#v=onepage&q=lanczos%20algorithm%20complexity&f=false] и количества операций для переортогонализации методом грамма-Шмидта[http://stackoverflow.com/questions/27986225/computational-complexity-of-gram-schmidt-orthogonalization-algorithm]. | ||
Итоговая сложность составляет <math>O(n^2)+O(nk^2) </math>. | Итоговая сложность составляет <math>O(n^2)+O(nk^2) </math>. | ||
| − | == Информационный граф == | + | === Информационный граф === |
[[Файл:Graph.jpg]] | [[Файл:Graph.jpg]] | ||
| − | == Входные и выходные данные алгоритма == | + | Вданной блок-схеме в качестве многоточий |
| + | === Входные и выходные данные алгоритма === | ||
На вход принимается сама матрица <math> A \in R^{n \times n}</math>, случайный вектор <math>b</math> (впрочем возможен вариант, при котором этот вектор генерируется случайным образом). | На вход принимается сама матрица <math> A \in R^{n \times n}</math>, случайный вектор <math>b</math> (впрочем возможен вариант, при котором этот вектор генерируется случайным образом). | ||
Ввиду симметричности входной матрицы достаточно передать лишь ее верхнюю треугольную матрицу, что дает объем входных данных <math>\frac{n(n + 1)}{2}</math>. | Ввиду симметричности входной матрицы достаточно передать лишь ее верхнюю треугольную матрицу, что дает объем входных данных <math>\frac{n(n + 1)}{2}</math>. | ||
| Строка 84: | Строка 86: | ||
Выходными данными для <math>j</math>ой итерации являются вектор <math> Aq_j </math> и собственное число <math> \lambda </math>. Общий объем выходных данных <math>k(n+1)</math>. | Выходными данными для <math>j</math>ой итерации являются вектор <math> Aq_j </math> и собственное число <math> \lambda </math>. Общий объем выходных данных <math>k(n+1)</math>. | ||
| − | ==Свойства алгоритма== | + | === Свойства алгоритма === |
Если <math>A</math> эрмитова матрица, то алгоритм Ланцоша и Bi-Lanczos сходятся к одинаковым трехдиагональным матрицам ритца[http://ta.twi.tudelft.nl/nw/users/vuik/papers/DUT-TWI-96-44.pdf]. | Если <math>A</math> эрмитова матрица, то алгоритм Ланцоша и Bi-Lanczos сходятся к одинаковым трехдиагональным матрицам ритца[http://ta.twi.tudelft.nl/nw/users/vuik/papers/DUT-TWI-96-44.pdf]. | ||
| Строка 91: | Строка 93: | ||
Алгоритм может завершить свою работу досрочно, когда найденные собственные значения будут достаточно близки к целевым. | Алгоритм может завершить свою работу досрочно, когда найденные собственные значения будут достаточно близки к целевым. | ||
| − | == | + | === Ресурс параллелизма алгоритма === |
Хотя алгоритм является итерационным, возможно распараллелить внутри каждой итерации умножение матрицы на вектор и процесс переортогонализации Грамма-Шмидта. | Хотя алгоритм является итерационным, возможно распараллелить внутри каждой итерации умножение матрицы на вектор и процесс переортогонализации Грамма-Шмидта. | ||
| Строка 99: | Строка 101: | ||
| − | == Библиотеки реализующие алгоритм== | + | === Библиотеки реализующие алгоритм === |
The IETL Project http://www.comp-phys.org/software/ietl/ C++ | The IETL Project http://www.comp-phys.org/software/ietl/ C++ | ||
Версия 23:36, 30 октября 2016
| Алгоритм Ланцоша с полной переортогонализацией | |
| Последовательный алгоритм | |
| Последовательная сложность | [math]O(n^2)+O(nk^2) [/math] |
| Объём входных данных | [math]\frac{n(n + 1)}{2}[/math] |
| Объём выходных данных | [math]k(n + 1)[/math] |
Авторы: Алексейчук Н.Н 616, Ястребов К.С. 609.
Содержание
- 1 Свойства и структура алгоритма
- 1.1 Общее описание алгоритма.
- 1.2 Математическое описание алгоритма
- 1.3 Вычислительное ядро алгоритма
- 1.4 Макроструктура алгоритма
- 1.5 Схема реализации последовательного алгоритма
- 1.6 Последовательная сложность алгоритма
- 1.7 Информационный граф
- 1.8 Входные и выходные данные алгоритма
- 1.9 Свойства алгоритма
- 1.10 Ресурс параллелизма алгоритма
- 1.11 Библиотеки реализующие алгоритм
- 2 Литература
1 Свойства и структура алгоритма
1.1 Общее описание алгоритма.
Алгоритм Ланцоша ищет собвственные значения и собственные векторы для симетричной матрицы A вещественных чисел. Является итерацонным алгоритмом. Алгоритм Ланцоша использует степенной метод ([math] b_{k+1} = \frac{Ab_k}{||Ab_k||} [/math]) для поиска наибольших собственных значений и векторов матриц.
В отличие от прямых алгоритмов требует мешьше памяти и мощности, что является несомненным плюсом для больших матриц.
Несмотря на свою скорость работы и экономию памяти, сначала не был популярным алгоритмом из – за недостаточной вычислительной устойчивости. В 1970 году Ojalvo и Newman [1] показали способ сделать алгоритм достаточно устойчивым. В этой же статье алгоритм был применен к расчету инженерной конструкции с большим количеством узлов, которые подвергались динамической нагрузке.
1.2 Математическое описание алгоритма
Памятка: Степенной метод нахождения наибольшего собственного числа матрицы можно сформулировать в предельном виде: если [math] b_0 [/math] – случайный вектор, и [math] b_n+1 = Ab_n [/math], тогда для больших чисел n предел [math]x_n/||x_n|| [/math] стремится к нормированному наибольшему собственному вектору.
Алгоритм Ланцоша комбинирует метод Ланцоша для нахождения крыловского подпространства и метод Релэя – Ритца.
Подпространство Крылова для степенного метода: [math] K_m(v,A) = span[x_1, Ax_1, A^2x_1, ..., A^{k-1}x_1] [/math]
В качестве входных данных для алгоритма Ланцоша подаются квадратная матрица размерности [math]n[/math]X[math]n[/math]: [math]A=A^T[/math]; а так же вектор начального приближения [math]b[/math].
Метод осуществляет поиск трехдиагональной симметричной матрицы [math]T_k=Q_k^TAQ_k[/math]. Причем собственные значения [math]T_k[/math] таковы, что приближают собственные значения исходной матрицы [math]A[/math]. То есть на каждом [math]k[/math]-м шаге из ортонормированных векторов Ланцоша строится матрица [math]Q_k = [q_1,q_2,...,q_k][/math] и в качестве приближенных собственных значений матрицы [math]A[/math] принимаются числа Ритца.
Из-за ошибок округления вектора [math] A^{k-1}x_1 [/math], формирующие подпространство Крылова, становятся неортогональными. Чтобы решить данную прорблему, проводят переортогонализацию методом Грамма-Шмидта.
1.3 Вычислительное ядро алгоритма
Вычислительным ядром данного алгоритма являются следующие шаги:
- 1. [math]Aq=( \sum\nolimits_{i=^n}a_{1i}q_i, \sum\nolimits_{i=2}^na_{2i}q_i, ..., \sum\nolimits_{i=1}^na_{ni}q_i)[/math].
- 2. [math]z=z-\sum\nolimits_{i=1}^{k}(z^Tq_i)q_i.[/math]
В данном случае первая операция выполняет умножение симметричной матрицы [math]A[/math] размерности [math]n[/math]X[math]n[/math] на [math]n[/math]-мерный вектор [math]q[/math], вследствие чего вычислительная сложность выполнения заключается в [math]n^2[/math] умножений и [math]n^2-n[/math] сложений. Второе действие является процессом ортогонализации Грама-Шмидта. В последнем действии вычисляются [math]k^2n+k(n+2)[/math] умножений и [math]k^2n + k(n + 1) + 2[/math] операций сложения.
1.4 Макроструктура алгоритма
- 1. Скалярное произведение: [math] (x,y) [/math].
- 2. Суммирование: [math] x+\alpha y [/math].
1.5 Схема реализации последовательного алгоритма
Input: A, b (random vector with unit norm)
- [math] \begin{align} q_1 = b/||b||_2, \beta_0 = 0, q_0 = 0 \\ j = 1 ,...,k \\ q_1&=b/||b||,\beta_0=0,q_o=0. \\ z&=Aq_j, \\ \alpha_j&=q_j^Tz, \\ z&=z-\alpha_jq_j-\beta_{j-1}q_{j-1}, \\ \beta&=||z||,\\ q_{j+1}&=z/\beta_j, \quad j \in [1, k]. \end{align} [/math]
1.6 Последовательная сложность алгоритма
Количество операций складывается из количества операций для классического метода Ланцоша[2] и количества операций для переортогонализации методом грамма-Шмидта[3].
Итоговая сложность составляет [math]O(n^2)+O(nk^2) [/math].
1.7 Информационный граф
Вданной блок-схеме в качестве многоточий
1.8 Входные и выходные данные алгоритма
На вход принимается сама матрица [math] A \in R^{n \times n}[/math], случайный вектор [math]b[/math] (впрочем возможен вариант, при котором этот вектор генерируется случайным образом). Ввиду симметричности входной матрицы достаточно передать лишь ее верхнюю треугольную матрицу, что дает объем входных данных [math]\frac{n(n + 1)}{2}[/math].
Выходными данными для [math]j[/math]ой итерации являются вектор [math] Aq_j [/math] и собственное число [math] \lambda [/math]. Общий объем выходных данных [math]k(n+1)[/math].
1.9 Свойства алгоритма
Если [math]A[/math] эрмитова матрица, то алгоритм Ланцоша и Bi-Lanczos сходятся к одинаковым трехдиагональным матрицам ритца[4].
При реализации классического алгоритма Ланцоша возникает большая погрешность при округлении. Вариант с полной переортогонализацией позволяет избегать больших погрешностей, однако является более ресурсоемким. Существует промежуточный вариант с частичной переортогонализацией.
Алгоритм может завершить свою работу досрочно, когда найденные собственные значения будут достаточно близки к целевым.
1.10 Ресурс параллелизма алгоритма
Хотя алгоритм является итерационным, возможно распараллелить внутри каждой итерации умножение матрицы на вектор и процесс переортогонализации Грамма-Шмидта.
Процесс умножения матриц matvec можно распараллелить несколькими способами[5].
1.11 Библиотеки реализующие алгоритм
The IETL Project http://www.comp-phys.org/software/ietl/ C++
NAG Library http://www.nag.com/content/nag-library C, C++, Fortran, C#, MATLAB, R
ARPACK https://people.sc.fsu.edu/~jburkardt/m_src/arpack/arpack.html MATLAB
GrapLab https://turi.com/products/create/open_source.html C++
С частичной переортаганализацией
LANSO/PLANSO http://web.cs.ucdavis.edu/~bai/ET/lanczos_methods/overview_PLANSO.html Fortran (уже распараллелена)
2 Литература
Ojalvo, I.U. and Newman, M., "Vibration modes of large structures by an automatic matrix-reduction method", AIAA J., 8 (7), 1234–1239 (1970).
