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

Авторы: Хакимов А. С.

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

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

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

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

[math]\displaystyle{ 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]\displaystyle{ ,\, \; }[/math]

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

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

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

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

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

[math]\displaystyle{  q_{1} = b_{j}/\|b\|_2, \beta_0=0, q_0=0 }[/math]
[math]\displaystyle{  for\, j=1\,\, to\, \, k }[/math]
    [math]\displaystyle{ z=Aq_j, }[/math]
    [math]\displaystyle{ \alpha_j=q_j^Tz, }[/math]
    [math]\displaystyle{ z=z-\alpha_jq_j-\beta_{j-1}q_{j-1}, }[/math]
    [math]\displaystyle{ \text{/* Провести выборочную ортогонализацию по отношению},   }[/math]
    [math]\displaystyle{ \text{   к сошедшимся векторам Ритца */} }[/math]
    [math]\displaystyle{ for\, i \leqslant k, \text{таких, что} \, \beta_{t}|v_i(t)| \leqslant \sqrt{\varepsilon}\|T_{k}\| \,  }[/math]
    [math]\displaystyle{ z = z-(y^T_{i,k},z)y_{i,k} }[/math]
    [math]\displaystyle{ end \, for }[/math]
    [math]\displaystyle{ \beta_{j}=\|z\|_2 }[/math]
    [math]\displaystyle{ q_{j+1}=z/\beta_{j},  }[/math]
    [math]\displaystyle{ \text{Вычислить собственные значения и собственные векторы} }[/math]
    [math]\displaystyle{ \text{матрицы} \, \, T_{j} \, \text{и оценки погрешности в них} }[/math]
[math]\displaystyle{ end \, for }[/math]