Участник:Макарова Алёна/Итерационный метод решения системы линейных алгебраических уравнений GMRES (обобщенный метод минимальных невязок): различия между версиями

Материал из Алговики
Перейти к навигации Перейти к поиску
(Описание)
(сходимость)
Строка 72: Строка 72:
  
 
</math>.
 
</math>.
 +
 +
==== Сходимость метода ====
 +
 +
Метод минимизации невязок на подпространствах Крылова за конечное число итераций приводит к точному решению и по этой причине является прямым методом. Но чаще его рассматривают как итерационный метод - с типичными для итерационных методов оценками сходимости. Идея заключается в том, что спустя некоторое количество итераций, полученное приближение уже является хорошей оценкой ответа.
 +
 +
Но об этом не верно говорить в общем. Согласно теореме (Greenbaum, Pták and Strakoš), для каждой монотонной последовательности <math> a_1,..., a_{m-1}, a_m = 0 </math>, найдётся матрица <math>A</math> такая что <math>||r_n|| = a_n</math> для всех <math>n</math>, где <math>r_n</math> малое большего порядка (is the residual defined above. ) ????
 +
 +
На практике для матрицы порядка  <math>n</math> число шагов может оказаться равным <math>n</math>, а норма невязки может оставаться равной норме начальной невязки на всех щагах, кроме последнего.
 +
 +
'''Пример'''
 +
 +
:<math>
 +
 +
\begin{bmatrix} 0 1 & ... \\ 0 0 1 & ... \\...\\... & 01 \\1 & ... \end{bmatrix} \begin{bmatrix} 1 \\ 0 \\...\\ 0 \\ 0 \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\...\\ 0 \\ 1 \end{bmatrix}
 +
 +
</math>
 +
 +
Если <math> x_0 = 0 </math>, то <math> r_0 = b</math>. Подпространства Крылова получаются такие:
 +
 +
:<math>
 +
 +
\mathcal{K}_i = span \{ e_n, e_{n-1}, ..., e_{n-i+1}\}, i \le i \le n
 +
 +
</math>
 +
 +
Точное решение системы имеет вид <math>x = e_1</math>, а на итерациях получаются следующие векторы:
 +
 +
:<math>
 +
 +
x_0 = x_1 = ... = x_{n-1} = 0, x_n = e_1
 +
 +
</math>.
 +
 +
Для получения оценок требуются дополнительные предположения относительно матрицы коэффициентов.
 +
 +
===== Сходимость для положительно определенной матрицы =====
 +
 +
Если матрица <math>A</math> положительно определена, то
 +
:<math>
 +
\|r_n\| \leq \left( 1-\frac{\lambda_{\mathrm{min}}^2(1/2(A^T + A))}{ \lambda_{\mathrm{max}}(A^T A)} \right)^{n/2} \|r_0\|,
 +
</math>,
 +
где <math>\lambda_{\mathrm{min}}(M)</math> и <math>\lambda_{\mathrm{max}}(M)</math> обозначают наибольшее и наименьшее собственное значение матрицы <math>M</math>, соответственно. <ref>Eisenstat, Elman & Schultz, Thm 3.3. NB all results for GCR also hold for GMRES, cf. Saad & Schultz</ref>
 +
 +
===== Сходимость для симметричной матрицы =====
 +
Если матрица <math>A</math> является симметрично и положительно определена, тогда
 +
 +
:<math>
 +
 +
\|r_n\| \leq \left( \frac{\kappa_2(A)^2-1}{\kappa_2(A)^2} \right)^{n/2} \|r_0\|
 +
 +
</math>,
 +
где <math>\kappa_2(A)</math> число обусловленности матрицы <math>A</math> в Евклидовой норме.
 +
 +
Если не требовать положительной определённости
 +
:<math>
 +
\|r_n\| \le \inf_{p \in P_n} \|p(A)\| \le \kappa_2(V) \inf_{p \in P_n} \max_{\lambda \in \sigma(A)} |p(\lambda)| \|r_0\|, \,
 +
</math>
 +
где <math>P_n</math> обозначает множество полиномов степени не выше <math>n</math> и <math>p(0) = 1</math>, матрица <math>V</math> появляется из спектрального разложения матрицы <math>A</math> и <math>\sigma(A)</math> - спектр матрицы <math>A</math>. Грубо говоря, это означает, что быстрая сходимость имеет место, когда собственные значения матрицы <math>A</math> группируются от происхождения и <math>A</math> не сильно отклоняется от нормальной матрицы (are clustered away from the origin and ''A'' is not too far from).<ref>Trefethen & Bau, Thm 35.2</ref>
 +
 +
Все эти неравенства связаны с невязкой вместо фактической ошибки, т.е. расстоянием между текущей итерацией и точным решением.
  
 
=== Вычислительное ядро алгоритма ===
 
=== Вычислительное ядро алгоритма ===

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

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

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

Метод минимальных невязок представляет собой итерационный метод для численного решения несимметричной системы линейных уравнений. Метод приближает решение вектором в подпространстве Крылова с минимальной невязкой. Чтобы найти этот вектор используется итерация Арнольди.

Впервые метод минимальных невязок был описан в книге Марчука и Кузнецова "Итерационные методы и квадратичные функционалы" [1] и представлял собой геометрическую реализацию метода минимальных невязок.

Алгебраическая реализация метода минимальных невязок была предложена Саадом и Шульцем в 1986 году [2]. С их подачи за методом закрепилась аббревиатура GMRES (обобщенный метод минимальных невязок).

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

Для решения системы [math]Ax = b[/math] с невырожденной матрицей [math]A[/math], выбирается начальное приближение [math]x_0[/math], затем решается редуцированная система [math]Au = r_0[/math], где [math]r_0 = b - Ax_0, x = x_0 + u[/math]. Подпространства Крылова строятся для редуцированной системы (в предположении, что [math]r_0 \neq 0[/math]).

Пусть [math]x_i = x_0 + y[/math], где [math] y \in \mathcal{K}_i [/math]. Невязка имеет вид [math]r_i = r_0 - Ay[/math], а её длина минимальна, в силу теоремы Пифагора, в том и только том случае, когда

[math] r_i \perp A\mathcal{K}_i [/math].

Таким образом, для реализации [math]i[/math]-го шага требуется опустить перпендикуляр из вектора [math]r_0[/math] на подпространство [math]A\mathcal{K}_i[/math]. Проще всего это сделать, если в данном подпространстве уже найден ортогональный базис.

1.2.1 Геометрическая реализация

Геометрическая реализация метода минимальных невязок заключается в построении последовательности векторов [math]q_1, q_2,...[/math] таким образом, что [math]q_1, q_2,.., q_i[/math] дают базис в подпространстве Крылова [math]\mathcal{K}_i[/math] и при этом вектора [math]p_1 = Aq_1, ..., p_i = Aq_i[/math] образуют ортогональный базис в подпространстве [math]A\mathcal{K}_i[/math]. Вектор [math]q_{i+1} \notin \mathcal{K}_i[/math] должен обладать следующими свойствами:

[math] q_{i+1} \in \mathcal{K}_{i+1}, p_{i+1} = Aq_{i+1} \perp A\mathcal{K}_i [/math].

Его можно получить с помощью процесса ортогонализации, применённого к вектору [math]p = Aq[/math], где [math]q = Aq_i[/math]. В геометрической реализации необходимо хранить две системы векторов: [math]q_1, ..., q_i[/math] и [math]p_1, ..., p_i[/math].

1.2.2 Алгебраическая реализация

Алгебраическая реализация метода минимальных невязок использует лишь одну систему векторов, образующих ортогональные базисы в подпространствах [math]\mathcal{K}_i[/math].

[math]q_1 = r_0/||r_0||_2[/math]. Чтобы получить ортонормированный базис [math]q_1, ..., q_{i+1}[/math] в [math]\mathcal{K}_{i+1} = \mathcal{K}_{i+1}(r_0, A)[/math], проводим ортогонализацию вектора [math]Aq_i[/math] к векторам [math]q_1, ..., q_i[/math]. Если [math]Q_i \equiv [q_1, ..., q_i] \in \mathbb{C}^{n*i}[/math], то получим

[math] AQ_i = Q_{i+1}\hat{H_i}, \hat{H_i}= \begin{bmatrix} H_i \\ 0 ... 0 & h_{i+1 i} \end{bmatrix}, [/math]

где [math]H_i[/math] - верхняя хессенбергова матрица порядка [math]i[/math].

Далее рассмотрим [math]QR[/math]-разложение прямоугольной матрицы [math]\hat{H_i}[/math]:

[math] \hat{H_i} = U_iR_i, R_i \in \mathbb{C}^{i*i} [/math].

Тогда минимум невязки [math]||r_0 - AQ_iy||_2[/math] по всем [math]y[/math] будет достигаться в том случае, если [math]y = y_i[/math] удовлетворяет уравнению

[math] R_iy_i = z_i \equiv ||r_0||_2U^*_ie_1, e_1 = \begin{bmatrix} 1 0 ... 0 \end{bmatrix}^T, [/math].

Матрица [math]R_i[/math] невырожденная, следовательно,

[math] x_i = x_0 + Q_iy_i = x_0 + Q_iR_i^{-1}z_i [/math].

1.2.3 Сходимость метода

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

Но об этом не верно говорить в общем. Согласно теореме (Greenbaum, Pták and Strakoš), для каждой монотонной последовательности [math] a_1,..., a_{m-1}, a_m = 0 [/math], найдётся матрица [math]A[/math] такая что [math]||r_n|| = a_n[/math] для всех [math]n[/math], где [math]r_n[/math] малое большего порядка (is the residual defined above. ) ????

На практике для матрицы порядка [math]n[/math] число шагов может оказаться равным [math]n[/math], а норма невязки может оставаться равной норме начальной невязки на всех щагах, кроме последнего.

Пример

[math] \begin{bmatrix} 0 1 & ... \\ 0 0 1 & ... \\...\\... & 01 \\1 & ... \end{bmatrix} \begin{bmatrix} 1 \\ 0 \\...\\ 0 \\ 0 \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\...\\ 0 \\ 1 \end{bmatrix} [/math]

Если [math] x_0 = 0 [/math], то [math] r_0 = b[/math]. Подпространства Крылова получаются такие:

[math] \mathcal{K}_i = span \{ e_n, e_{n-1}, ..., e_{n-i+1}\}, i \le i \le n [/math]

Точное решение системы имеет вид [math]x = e_1[/math], а на итерациях получаются следующие векторы:

[math] x_0 = x_1 = ... = x_{n-1} = 0, x_n = e_1 [/math].

Для получения оценок требуются дополнительные предположения относительно матрицы коэффициентов.

1.2.3.1 Сходимость для положительно определенной матрицы

Если матрица [math]A[/math] положительно определена, то

[math] \|r_n\| \leq \left( 1-\frac{\lambda_{\mathrm{min}}^2(1/2(A^T + A))}{ \lambda_{\mathrm{max}}(A^T A)} \right)^{n/2} \|r_0\|, [/math],

где [math]\lambda_{\mathrm{min}}(M)[/math] и [math]\lambda_{\mathrm{max}}(M)[/math] обозначают наибольшее и наименьшее собственное значение матрицы [math]M[/math], соответственно. [3]

1.2.3.2 Сходимость для симметричной матрицы

Если матрица [math]A[/math] является симметрично и положительно определена, тогда

[math] \|r_n\| \leq \left( \frac{\kappa_2(A)^2-1}{\kappa_2(A)^2} \right)^{n/2} \|r_0\| [/math],

где [math]\kappa_2(A)[/math] число обусловленности матрицы [math]A[/math] в Евклидовой норме.

Если не требовать положительной определённости

[math] \|r_n\| \le \inf_{p \in P_n} \|p(A)\| \le \kappa_2(V) \inf_{p \in P_n} \max_{\lambda \in \sigma(A)} |p(\lambda)| \|r_0\|, \, [/math]

где [math]P_n[/math] обозначает множество полиномов степени не выше [math]n[/math] и [math]p(0) = 1[/math], матрица [math]V[/math] появляется из спектрального разложения матрицы [math]A[/math] и [math]\sigma(A)[/math] - спектр матрицы [math]A[/math]. Грубо говоря, это означает, что быстрая сходимость имеет место, когда собственные значения матрицы [math]A[/math] группируются от происхождения и [math]A[/math] не сильно отклоняется от нормальной матрицы (are clustered away from the origin and A is not too far from).[4]

Все эти неравенства связаны с невязкой вместо фактической ошибки, т.е. расстоянием между текущей итерацией и точным решением.

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

Главные затраты [math]i[/math]-ой итерации связаны с ортогонализацией.

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

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

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

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

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

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

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

2 Литература

  1. Г.И. Марчук, Ю.А. Кузнецов, Итерационные методы и квадратичные функционалы, Методы вычислительной математики, Новосибирск, 1975, с. 4-143.
  2. Y.Saad, M.H. Schultz, GMRES: A generalized minimal residual algorithm for solving nonsymmetric linear systems, SIAM J. Scientific and Stat. Comp. 7: 856-869 (1986).
  3. Eisenstat, Elman & Schultz, Thm 3.3. NB all results for GCR also hold for GMRES, cf. Saad & Schultz
  4. Trefethen & Bau, Thm 35.2