Участник:ArtyomKhakimov/Алгоритм Ланцоша для арифметики с плавающей точкой с выборочной ортогонализацией
| Файл:Symbol wait.svg | Эта работа ждет рассмотрения преподавателем Дата последней правки страницы: 19.01.2017 Авторы этой статьи считают, что задание выполнено. |
Авторы: Хакимов А. С.
Свойства и структура алгоритма
Общее описание алгоритма
Алгоритм Ланцоша служит для нахождения собственных значений и собственных векторов для больших разреженных матриц, к которым нельзя применить прямые методы из-за больших требований к памяти и времени. Он был опубликован Корнелием Ланцошем в 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], например, с помощью метода "разделяй и властвуй"[1]
Математическое описание алгоритма
Вычислительное ядро алгоритма
- Умножение матрицы на вектор, [math]\displaystyle{ z=Aq_j, }[/math].
- Ортогонализация по отношению к сошедшимся векторам [math]\displaystyle{ z = z-(y^T_{i,k},z)y_{i,k} }[/math] для [math]\displaystyle{ i = 1 \dots t }[/math]
Макроструктура алгоритма
- ↑ Дж. Деммель «Вычислительная линейная алгебра», c. 232, алгоритм 5.2