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

Материал из Алговики
Перейти к навигации Перейти к поиску
Строка 35: Строка 35:
 
== Математическое описание алгоритма ==
 
== Математическое описание алгоритма ==
  
−
  <math> q_{1} = b_{j}/\|b\|_2, \beta_0=0, q_0=0</math>
+
== Вычислительное ядро алгоритма ==
−
  <math> for\, j=1\,\, to\, \, k</math>
+
 
−
    <math>z=Aq_j,</math>
+
* Умножение матрицы на вектор, <math>z=Aq_j,  </math>.
−
    <math>\alpha_j=q_j^Tz,</math>
+
 
−
    <math>z=z-\alpha_jq_j-\beta_{j-1}q_{j-1},</math>
+
* Ортогонализация по отношению к сошедшимся векторам <math>z = z-(y^T_{i,k},z)y_{i,k} </math> для <math>i  = 1 \dots t</math>
−
    /* Провести выборочную ортогонализацию по отношению
+
 
−
    к сошедшимся векторам Ритца */
+
== Макроструктура алгоритма ==
−
    <math>for\, i \leqslant k, \,</math>таких, что<math \, \beta_{t}|v_i(t)| \leqslant \sqrt{\varepsilon}\|T_{k}\| \, </math>
+
 
−
    <math>z = z-(y^T_{i,k},z)y_{i,k}</math>
+
Для алгоритма Ланцоша можно выделить следующие макрооперации:
−
    <math>\text{end for}</math>
+
 
−
    <math>\beta_{j}=\|z\|_2</math>
+
* Умножение матрицы на вектор. Состоит из операций умножения вектора на число и сложения векторов.
−
    <math>q_{j+1}=z/\beta_{j}, </math>
+
 
−
    Вычислить собственные значения и собственные векторы
+
<math>z=Aq_j,</math>.
−
    матрицы<math>\, \, T_{j} \, \,</math>и оценки погрешности в них
+
 
−
<math>\text{end for}</math>
+
* Вычисление вектора <math>q_{j+1}</math> как скалярное произведение двух других векторов.
 +
 
 +
<math>\alpha_j=q_j^Tz,</math>
 +
 
 +
<math>z=z-\alpha_jq_j-\beta_{j-1}q_{j-1},</math>
 +
 
 +
<math>q_{j+1}=z/\|z\|_2</math>
 +
 
 +
* Ортогонализация вектора Ланцоша с помощью скалярного произведения:
 +
 
 +
<math>z = z-(y^T_{i,t},z)y_{i,t}</math>
 +
 
 +
* Алгоритм "разделяй и властвуй" для вычисления собственных значений матрицы.
 +
 
 +
== Схема реализации последовательного алгоритма ==
 +
 
 +
Алгоритм выполняется в следующей последовательности:
 +
 +
<math>1.\, \, \beta_0 = 0,\; q_0 = 0</math> #Инициализируем векторы
 +
 +
  <math>2.\, \, \|b\|_2 = \sqrt{\sum\limits_{j=1}^{n} b_j^2}</math> #Вычисляем норму вектора начального приближения.
 +
 +
<math>3.\, \, q_{1_{j}} = \frac{b_{j}}{\|b\|_2}, \; j = 1,\, \dots\,, n</math> #Нормализуем вектор начального приближения.
 +
 +
<math>4.\, \,</math> Начинается цикл, в котором не более чем <math>k</math>. На <math>i</math>-й итерации производятся следующие вычисления:
 +
 +
<math>4.1.\, \, z_j = \sum\limits_{m=1}^{n} a_{jm} q_{i_m}, \; j = 1,\,\dots\,, n</math> #Считаем результат применения линейного оператора <math>A</math> к вектору <math>q_i</math>.
 +
 +
<math>4.2.\, \, \alpha_i = \sum\limits_{j=1}^{n}q_{i_j} z_j</math> #Получаем результат скалярного произведения векторов <math>q_i</math> и <math>z</math>.
 +
 +
<math>4.3.\, \, z_j = z_j - \alpha_i q_{i_j} - \beta_{i-1}q_{i-1_j}, \, j = 1,\,\dots\,, n</math> #Вычисляем линейную комбинацию векторов.
 +
 +
<math>4.4.\, \,</math> Начинается цикл, в котором производится выборочная ортогонализация векторов Ланцоша
 +
 +
<math>4.4.1.\, \,</math> Вычисляется вектор Ритца <math>y_{i,k_d} = \sum\limits_{x=1}^{k} q_{x_d} v_{i_x}, \; d = 1,\,\dots\,, n,</math>, где  <math>y_{i,k_d}</math> -  <math>d</math>-я компонента вектора Ритца <math>y,</math>
 +
<math>4.4.2.\, \,</math> Производится скалярное умножение <math>y_{i,k}^T</math>  и <math> z: \; \; \gamma = \sum\limits_{d=1}^{n}y_{i,k_d} z_d,</math>
 +
<math>4.4.3.\, \,</math> Ищется новое собственное значение <math> \theta_j </math> и собственный вектор <math> v_j </math> для полученной матрицы <math> T_j,</math>
 +
 +
<math>4.5.\, \, \beta_i = \|z\|_2 = \sqrt{\sum\limits_{j=1}^{n} z_j^2}</math> #Считаем норму вектора <math>z</math>.
 +
 +
<math>4.6.</math> Проверка равенства <math>\beta_i == 0</math> # Если норма оказалась равной нулю, то завершаем итерации и переходим к вычислению собственных векторов и собственных значений полученной матрицы. В обратном случае, продолжаем выполнения итераций.
 +
 +
<math>4.7.\, \, q_{i+1_j} = \frac{z_j}{\beta_i}, \; j = 1,\, \dots \,, n</math> #Нормируем вектор <math>z</math>.
 +
 +
<math>4.8.\,</math> Если выполнили <math>k</math> итераций, то завершаем выполнение итераций, переходим к следующему шагу. Иначе начинаем последующую итерацию цикла.
 +
 +
<math>5.</math> Вычисляем собственные значения и собственные вектора полученной матрицы <math>T_k</math>.
 +
 +
== Последовательная сложность алгоритма ==
 +
 
 +
* Вычисление вектора <math>z</math> - <math>n^2</math> операций умножения и <math>n(n-1)</math> операций сложения
 +
* Вычисление скалярного произведения векторов - <math>n</math> операций умножения и <math>(n-1)</math> сложений
 +
* Вычисление вектора <math>z</math> - <math>2n</math> умножений и <math>2n</math> сложений
 +
* Оценка погрешностей - <math>kn</math> операций умножения
 +
* Выборочная переортогонализация - <math>O(2nj)</math> умножений и <math>O(2j(n-1) + n)</math> сложений
 +
* Вычисление нормы вектора - <math>n</math> операций умножения и <math>(n-1)</math> сложений
 +
* Вычисление нового столбца матрицы <math>Q</math> - <math>n</math> операций умножения
 +
 
 +
Итого на каждой итерации:
 +
 
 +
* Умножений: <math>n^2 + n + 2n + kn + O(2nj) + n + n = n^2 + 5n + kn + O(2nj)</math>
 +
* Сложений: <math>n(n-1) + n-1 + 2n + O(2j(n-1) + n) + n-1 = n^2 + 3n + O(2j(n-1) + n) -2</math>
 +
 
 +
Итого за все время выполнения алгоритма:
 +
 
 +
* Умножений: <math>kn^2 + 5kn + k^2n + O(k^2n + kn)</math>
 +
* Сложений: <math>kn^2 + 3kn + O((k^2+k)(n-1)+kn) = kn^2 + 3kn + O(k^2(n-1)+2kn-k)</math>
 +
 
 +
== Информационный граф ==
 +
 
 +
[[File:anczos_graph.png|1024px|center|рис. 1: Информационный граф]]
 +
Здесь:
 +
<math>Ax</math>  — перемножение матрицы и вектора,
 +
<math>\|x\|</math>  — макрооперация взятия нормы,
 +
<math>\alpha x</math>  — умножение скаляра на вектор,
 +
<math>\cdot </math>  — операция скалярного умножения векторов,
 +
<math>\mathsf{F}</math>  — операция <math>z-\sum\nolimits_{i=1}^{j-1}(z^Tq_i)q_i</math>,
 +
<math>\mathsf{f}</math>  — вычисление линейной комбинации векторов
 +
 
 +
== Ресурс параллелизма алгоритма ==
 +
 
 +
На каждой из итераций алгоритма можно получить выигрыш за счет распараллеливания следующих шагов:
 +
 
 +
* Вычисление вектора <math>z</math> - умножение матрицы размера <math>n \times n</math> на вектор длины <math>n</math> - <math>n</math> ярусов сложений, <math>n</math> операций умножений в каждом
 +
* Вычисление коэффициента <math>\alpha_j</math> - скалярное произведение векторов - <math>n</math> ярусов сложений, 1 операция умножения в каждом
 +
* Вычисление нового значения вектора <math>z</math>: 2 яруса умножений длины <math>n</math> (умножение вектора на число), 2 яруса вычитаний длины <math>n</math>
 +
* Выборочная переортогонализация: последовательно для всех <math>i \le k</math> выполняется <math>j</math> ярусов сложений, <math>n</math> операций умножений в каждом, <math>n</math> ярусов сложений, <math>j</math> операций умножений в каждом, один ярус вычитаний размера <math>n</math>
 +
* Вычисление <math>\beta_j</math> - скалярное произведение векторов - <math>n</math> ярусов сложений, 1 операция умножения в каждом
 +
* Вычисление <math>q_j</math> - деление вектора на число - один ярус делений размера <math>n</math>
 +
 
 +
Сложность алгоритма по ширине ЯПФ - <math>O(n)</math>.
 +
 
 +
== Входные и выходные данные алгоритма ==
 +
 
 +
Входные данные: вещественная матрица <math>A</math> размера <math>n \times n</math>, вектор <math>b</math> длины <math>n</math>, число <math>k</math>
 +
 
 +
Объем входных данных: <math>n^2 + n + 1</math>
 +
 
 +
Выходные данные: матрица <math>Q</math> размера <math>n \times k</math>
 +
 
 +
Объем выходных данных: <math>nk</math>
 +
 
 +
== Свойства алгоритма ==
 +
 
 +
= Программная реализация =
 +
 
 +
== Особенности реализации последовательного алгоритма ==
 +
 
 +
== Локальность данных и вычислений ==
 +
 
 +
== Возможные способы и особенности параллельной реализации алгоритма ==
 +
 
 +
== Масштабируемость алгоритма и его реализации ==
 +
 
 +
Исследование алгоритма на масштабируемость производилось на суперкомпьютере [http://parallel.ru/cluster "Ломоносов"] с использованием технологии MPI. Для тестирования была выбрана [http://pastebin.com/SbNrh7nz вот эта] реализация алгоритма.
 +
 
 +
Набор параметров, которые изменялись при тестировании:
 +
*Число процессоров от 16 до 128.
 +
*Размерность матрицы от 5000 до 50000 с шагом 5000.
 +
 
 +
{| class="wikitable" style="border: none; background: none;"
 +
|+ align="bottom" style="caption-side: bottom" | Табл. 1. Результаты тестирования. В ячейках таблицы -- время работы алгоритма в секундах.
 +
! colspan="2" rowspan="2" style="border: none; background: none;"|
 +
! colspan="4"| Число процессоров
 +
|-
 +
! 16 !! 32 !! 64 !! 128
 +
|-
 +
! rowspan="10"| Размер матрицы
 +
! 5000
 +
| 1,07 || 0,71 || 0,62 || 0,54
 +
|-
 +
! 10000
 +
| 2,39 || 1,94 || 1,17 || 0,79
 +
|-
 +
! 15000
 +
| 4,35 || 3,38 || 1,89 || 1,59
 +
|-
 +
! 20000
 +
| 7,63 || 5,19 || 3,78 || 2,35
 +
|-
 +
! 15000
 +
| 12,41 || 7,33 || 4,32 || 3,44
 +
|-
 +
! 30000
 +
| 16,68 || 9,36 || 5,02 || 4,68
 +
|-
 +
! 35000
 +
| 22,58 || 11,92 || 8,04 || 6,20
 +
|-
 +
! 40000
 +
| 29,31 || 15,90 || 8,47 || 7,43
 +
|-
 +
! 45000
 +
| 37,04 || 19,42 || 10,55 || 9,04
 +
|-
 +
! 50000
 +
| 45,64 || 23,39 || 16,09 || 10,25
 +
|}
 +
 
 +
[[File:grph.png|1024px|center|рис. 1: Информационный граф]]
 +
 
 +
== Динамические характеристики и эффективность реализации алгоритма ==
 +
 
 +
== Выводы для классов архитектур ==
 +
 
 +
== Существующие реализации ==
 +
 
 +
Реализация в проекте Vienna VL<ref>https://github.com/viennacl/viennacl-dev</ref>
 +
 
 +
Реализация в проекте IETL<ref>http://www.comp-phys.org/software/ietl/lanczos.html</ref>
 +
 
 +
Реализация в проекте ARPACK <ref>http://www.caam.rice.edu/software/ARPACK/</ref>
 +
 
 +
Lanczos Plus Plus <ref>https://github.com/g1257/LanczosPlusPlus</ref>
 +
 
 +
ED Lanczos <ref>https://github.com/henhans/ED_lanczos</ref>
 +
 
 +
Paralleles Rechnen 2 <ref>https://github.com/soneyworld/ParallelesRechnen2</ref>
 +
 
 +
Cuda-Arnoldi <ref>https://github.com/trantalaiho/Cuda-Arnoldi/</ref>
 +
 
 +
= Литература =

Версия 12:17, 19 января 2017

Шаблон:Assignmenta Авторы: Хакимов А. С.

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

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

Алгоритм Ланцоша служит для нахождения собственных значений и собственных векторов для больших разреженных матриц, к которым нельзя применить прямые методы из-за больших требований к памяти и времени. Он был опубликован Корнелием Ланцошем в 1950 году. Его эффективность обусловлена экономией памяти для хранения матриц и экономией вычислительных ресурсов. Алгоритм итерационный и использует степенной метод для поиска наибольших собственных значений и векторов матриц. Основной недостаток алгоритма заключается в накоплении ошибок округления, для решения которых появились методы поддержания ортогонализации т.н. векторов Ланцоша. Здесь мы рассмотрим выборочный метод поддержания ортогонализации, который существенно экономит процессорное время.


На вход алгоритма подаётся [math]A = A^T[/math],


[math] A = \begin{pmatrix} a_{11} & a_{12} & a_{13} & \cdots & a_{1\ n-1} & a_{1\ n} \\ a_{12} & a_{22} & a_{23} & \cdots & a_{2\ n-1} & a_{2\ n} \\ a_{13} & a_{23} & a_{33} & \cdots & a_{3\ n-1} & a_{3\ n} \\ \vdots & \vdots & \ddots & \ddots & \ddots & \vdots \\ a_{1\ n-1} & \cdots & \cdots & a_{n-2\ n-1} & a_{n-1\ n-1} & a_{n-1\ n} \\ a_{1\ n} & \cdots & \cdots & a_{n-2\ n} & a_{n-1\ n} & a_{n\ n} \\ \end{pmatrix} [/math] [math] ,\, \;[/math]


случайный вектор [math]b[/math], как первое приближение собственного вектора матрицы и [math]k [/math] - количество собственных значений и собственных векторов, которые требуется найти.

Матрица [math]Q_j = [q_1, q_2, \dots, q_j][/math] размерности [math]n \times j[/math] строится на каждой итерации и состоит из ортонормированных векторов Ланцоша. А в качестве приближенных собственных значений берутся числа Ритца [math]\theta_i [/math], - собственные значения симметричной трехдиагональной матрицы [math]T_j = Q^T_j A Q_j[/math] размерности [math]j \times j[/math].

[math] T_j = \begin{pmatrix} \alpha_1 & \beta_1 \\ \beta_1 & \alpha_2 & \beta_2 \\ & \beta_2 & \ddots & \ddots \\ & & \ddots & \ddots & \beta_{j-1} \\ & & & \beta_{j-1} & \alpha_j \end{pmatrix}\; (2). [/math]


Однако, векторы [math]q_j [/math] теряют ортогональность вследствие приобретения больших компонент в направлениях векторов Ритца [math]y_{i,j} = Q_j v_i [/math], отвечающих сошедшимся числам Ритца [math] \theta_i [/math]. Поэтому чтобы построить [math]q_j [/math], предлагается на каждом шаге следить за оценками погрешностей [math]\beta_{t}|v_i(t)|, i = 1 \dots t, t = j - 1 [/math], где [math]v_i(t) [/math] - [math]t[/math]-я компонента собственного вектора [math]v_i [/math]. И когда какая-то оценка становится слишком малой, проводить ортогонализацию вектора Ланцоша [math]z [/math]. Величина [math]\beta_{t}|v_i(t)| [/math] считается малой, если она меньше, чем [math]\sqrt{\varepsilon}||T_{t}|| [/math], где [math]\varepsilon[/math] - доступная машинная точность чисел.

После следует вычисление собственных значений [math] \theta_j [/math] и собственных векторов [math]v_j [/math] полученной трехдиагональной матрицы [math]T_j[/math], например, с помощью метода "разделяй и властвуй"[1]

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

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

  • Умножение матрицы на вектор, [math]z=Aq_j, [/math].
  • Ортогонализация по отношению к сошедшимся векторам [math]z = z-(y^T_{i,k},z)y_{i,k} [/math] для [math]i = 1 \dots t[/math]

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

Для алгоритма Ланцоша можно выделить следующие макрооперации:

  • Умножение матрицы на вектор. Состоит из операций умножения вектора на число и сложения векторов.

[math]z=Aq_j,[/math].

  • Вычисление вектора [math]q_{j+1}[/math] как скалярное произведение двух других векторов.

[math]\alpha_j=q_j^Tz,[/math]

[math]z=z-\alpha_jq_j-\beta_{j-1}q_{j-1},[/math]

[math]q_{j+1}=z/\|z\|_2[/math]

  • Ортогонализация вектора Ланцоша с помощью скалярного произведения:

[math]z = z-(y^T_{i,t},z)y_{i,t}[/math]

  • Алгоритм "разделяй и властвуй" для вычисления собственных значений матрицы.

1.5 Схема реализации последовательного алгоритма

Алгоритм выполняется в следующей последовательности: 

[math]1.\, \, \beta_0 = 0,\; q_0 = 0[/math] #Инициализируем векторы

[math]2.\, \, \|b\|_2 = \sqrt{\sum\limits_{j=1}^{n} b_j^2}[/math] #Вычисляем норму вектора начального приближения.

[math]3.\, \, q_{1_{j}} = \frac{b_{j}}{\|b\|_2}, \; j = 1,\, \dots\,, n[/math] #Нормализуем вектор начального приближения.

[math]4.\, \,[/math] Начинается цикл, в котором не более чем [math]k[/math]. На [math]i[/math]-й итерации производятся следующие вычисления:

[math]4.1.\, \, z_j = \sum\limits_{m=1}^{n} a_{jm} q_{i_m}, \; j = 1,\,\dots\,, n[/math] #Считаем результат применения линейного оператора [math]A[/math] к вектору [math]q_i[/math].

[math]4.2.\, \, \alpha_i = \sum\limits_{j=1}^{n}q_{i_j} z_j[/math] #Получаем результат скалярного произведения векторов [math]q_i[/math] и [math]z[/math].

[math]4.3.\, \, z_j = z_j - \alpha_i q_{i_j} - \beta_{i-1}q_{i-1_j}, \, j = 1,\,\dots\,, n[/math] #Вычисляем линейную комбинацию векторов.

[math]4.4.\, \,[/math] Начинается цикл, в котором производится выборочная ортогонализация векторов Ланцоша

[math]4.4.1.\, \,[/math] Вычисляется вектор Ритца [math]y_{i,k_d} = \sum\limits_{x=1}^{k} q_{x_d} v_{i_x}, \; d = 1,\,\dots\,, n,[/math], где  [math]y_{i,k_d}[/math] -  [math]d[/math]-я компонента вектора Ритца [math]y,[/math]
[math]4.4.2.\, \,[/math] Производится скалярное умножение [math]y_{i,k}^T[/math]  и [math] z: \; \; \gamma = \sum\limits_{d=1}^{n}y_{i,k_d} z_d,[/math]
[math]4.4.3.\, \,[/math] Ищется новое собственное значение [math] \theta_j [/math] и собственный вектор [math] v_j [/math] для полученной матрицы [math] T_j,[/math]

[math]4.5.\, \, \beta_i = \|z\|_2 = \sqrt{\sum\limits_{j=1}^{n} z_j^2}[/math] #Считаем норму вектора [math]z[/math].

[math]4.6.[/math] Проверка равенства [math]\beta_i == 0[/math] # Если норма оказалась равной нулю, то завершаем итерации и переходим к вычислению собственных векторов и собственных значений полученной матрицы. В обратном случае, продолжаем выполнения итераций.

[math]4.7.\, \, q_{i+1_j} = \frac{z_j}{\beta_i}, \; j = 1,\, \dots \,, n[/math] #Нормируем вектор [math]z[/math].

[math]4.8.\,[/math] Если выполнили [math]k[/math] итераций, то завершаем выполнение итераций, переходим к следующему шагу. Иначе начинаем последующую итерацию цикла.

[math]5.[/math] Вычисляем собственные значения и собственные вектора полученной матрицы [math]T_k[/math].

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

  • Вычисление вектора [math]z[/math] - [math]n^2[/math] операций умножения и [math]n(n-1)[/math] операций сложения
  • Вычисление скалярного произведения векторов - [math]n[/math] операций умножения и [math](n-1)[/math] сложений
  • Вычисление вектора [math]z[/math] - [math]2n[/math] умножений и [math]2n[/math] сложений
  • Оценка погрешностей - [math]kn[/math] операций умножения
  • Выборочная переортогонализация - [math]O(2nj)[/math] умножений и [math]O(2j(n-1) + n)[/math] сложений
  • Вычисление нормы вектора - [math]n[/math] операций умножения и [math](n-1)[/math] сложений
  • Вычисление нового столбца матрицы [math]Q[/math] - [math]n[/math] операций умножения

Итого на каждой итерации:

  • Умножений: [math]n^2 + n + 2n + kn + O(2nj) + n + n = n^2 + 5n + kn + O(2nj)[/math]
  • Сложений: [math]n(n-1) + n-1 + 2n + O(2j(n-1) + n) + n-1 = n^2 + 3n + O(2j(n-1) + n) -2[/math]

Итого за все время выполнения алгоритма:

  • Умножений: [math]kn^2 + 5kn + k^2n + O(k^2n + kn)[/math]
  • Сложений: [math]kn^2 + 3kn + O((k^2+k)(n-1)+kn) = kn^2 + 3kn + O(k^2(n-1)+2kn-k)[/math]

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

рис. 1: Информационный граф
Здесь:
[math]Ax[/math]  — перемножение матрицы и вектора,
[math]\|x\|[/math]  — макрооперация взятия нормы,
[math]\alpha x[/math]  — умножение скаляра на вектор,
[math]\cdot [/math]  — операция скалярного умножения векторов,
[math]\mathsf{F}[/math]  — операция [math]z-\sum\nolimits_{i=1}^{j-1}(z^Tq_i)q_i[/math],
[math]\mathsf{f}[/math]  — вычисление линейной комбинации векторов

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

На каждой из итераций алгоритма можно получить выигрыш за счет распараллеливания следующих шагов:

  • Вычисление вектора [math]z[/math] - умножение матрицы размера [math]n \times n[/math] на вектор длины [math]n[/math] - [math]n[/math] ярусов сложений, [math]n[/math] операций умножений в каждом
  • Вычисление коэффициента [math]\alpha_j[/math] - скалярное произведение векторов - [math]n[/math] ярусов сложений, 1 операция умножения в каждом
  • Вычисление нового значения вектора [math]z[/math]: 2 яруса умножений длины [math]n[/math] (умножение вектора на число), 2 яруса вычитаний длины [math]n[/math]
  • Выборочная переортогонализация: последовательно для всех [math]i \le k[/math] выполняется [math]j[/math] ярусов сложений, [math]n[/math] операций умножений в каждом, [math]n[/math] ярусов сложений, [math]j[/math] операций умножений в каждом, один ярус вычитаний размера [math]n[/math]
  • Вычисление [math]\beta_j[/math] - скалярное произведение векторов - [math]n[/math] ярусов сложений, 1 операция умножения в каждом
  • Вычисление [math]q_j[/math] - деление вектора на число - один ярус делений размера [math]n[/math]

Сложность алгоритма по ширине ЯПФ - [math]O(n)[/math].

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

Входные данные: вещественная матрица [math]A[/math] размера [math]n \times n[/math], вектор [math]b[/math] длины [math]n[/math], число [math]k[/math]

Объем входных данных: [math]n^2 + n + 1[/math]

Выходные данные: матрица [math]Q[/math] размера [math]n \times k[/math]

Объем выходных данных: [math]nk[/math]

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

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

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

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

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

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

Исследование алгоритма на масштабируемость производилось на суперкомпьютере "Ломоносов" с использованием технологии MPI. Для тестирования была выбрана вот эта реализация алгоритма.

Набор параметров, которые изменялись при тестировании:

  • Число процессоров от 16 до 128.
  • Размерность матрицы от 5000 до 50000 с шагом 5000.
Табл. 1. Результаты тестирования. В ячейках таблицы -- время работы алгоритма в секундах.
Число процессоров
16 32 64 128
Размер матрицы 5000 1,07 0,71 0,62 0,54
10000 2,39 1,94 1,17 0,79
15000 4,35 3,38 1,89 1,59
20000 7,63 5,19 3,78 2,35
15000 12,41 7,33 4,32 3,44
30000 16,68 9,36 5,02 4,68
35000 22,58 11,92 8,04 6,20
40000 29,31 15,90 8,47 7,43
45000 37,04 19,42 10,55 9,04
50000 45,64 23,39 16,09 10,25
рис. 1: Информационный граф

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

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

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

Реализация в проекте Vienna VL[2]

Реализация в проекте IETL[3]

Реализация в проекте ARPACK [4]

Lanczos Plus Plus [5]

ED Lanczos [6]

Paralleles Rechnen 2 [7]

Cuda-Arnoldi [8]

3 Литература

  1. ↑ Дж. Деммель «Вычислительная линейная алгебра», c. 232, алгоритм 5.2
  2. ↑ https://github.com/viennacl/viennacl-dev
  3. ↑ http://www.comp-phys.org/software/ietl/lanczos.html
  4. ↑ http://www.caam.rice.edu/software/ARPACK/
  5. ↑ https://github.com/g1257/LanczosPlusPlus
  6. ↑ https://github.com/henhans/ED_lanczos
  7. ↑ https://github.com/soneyworld/ParallelesRechnen2
  8. ↑ https://github.com/trantalaiho/Cuda-Arnoldi/
Источник — https://algowiki-project.org/w/ru/index.php?title=Участник:ArtyomKhakimov/Алгоритм_Ланцоша_для_арифметики_с_плавающей_точкой_с_выборочной_ортогонализацией&oldid=21573