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

Основные авторы описания: Анюшева Ирина (614 группа), Исхаков Эльдар (614 группа)

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

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

Обобщенный метод минимальных невязок (GMRES) для решения системы линейных алгебраических уравнений аппроксимирует решение с помощью вектора в подпространства Крылова с минимальным остатком. Метод GMRES разработан Ю.Саадом и Мартин Х. Шульцом в 1986 году. Метод по праву считается одним из самых эффективных численных методов решения несимметричных систем. Он гарантирует неувеличение нормы невязки в ходе итерационного процесса, который строится основан на построении базиса в соответствующей системе подпространстве Крылова [math]\displaystyle{ K_j(A, r_0) }[/math]. Затем решение уточняется некоторой добавкой, представленной в виде разложения по этому базису.

Метод обобщенных минимальных невязок популярен из-за наличия ряда преимуществ: он ошибкоустойчив, допускает эффективное распараллеливание, не требует нахождения параметра релаксации, обладает суперлинейной скоростью сходимости. Однако в чистом виде для больших систем метод не используется из-за чрезвычайно больших затрат памяти. Зато широкое распространение получила перезапускаемая версия метода, подоразумевающая перезапуск метода каждые m итераций. Уменьшение размерности пространства может привести к стагнации метода – процессу, при котором с каждой следующей итерацией норма невязки уменьшается незначительно, либо не уменьшается вовсе.

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

Исходные данные: система линейных алгебраических уравнений [math]\displaystyle{ Ax=b }[/math], где:

[math]\displaystyle{ A=(a_{ij}) }[/math] [math]\displaystyle{ - }[/math] вещественная матрица размера [math]\displaystyle{ n \times n }[/math];

[math]\displaystyle{ b }[/math] и [math]\displaystyle{ x }[/math] вектора из [math]\displaystyle{ n }[/math] элементов;

[math]\displaystyle{ x^* }[/math] [math]\displaystyle{ - }[/math] точное решение системы.

Вычисляемые данные: [math]\displaystyle{ x^{(s)} - }[/math]приближенное решение системы.

Ортогонализация Арнольди

Для построения базиса [math]\displaystyle{ V }[/math] в подпространстве Крылова в методе GMRES применяется ортогонализация Арнольди.

Формулы метода:

[math]\displaystyle{ \begin{align} (h_{i+1,i}v_{i+1}) & = Av_i - \sum_{j = 1}^{i} h_{j,i}v_j, \\ v_1 &= \frac{r_0}{\parallel r_0 \parallel}, \end{align} }[/math]

В матрично-векторных обозначениях соотношение может быть записано как

[math]\displaystyle{ \begin{align} AV_i &= V_i H_i + h_{i+1,i}v_{i+1} e_i^T &=V_{i+1}\overline{H_i}, \end{align} }[/math]

[math]\displaystyle{ \begin{align} V_i &= [v_1|v_2|...|v_i], \overline{H_i} - \end{align} }[/math] верхняя матрица Хессенберга размерности [math]\displaystyle{ (i+1)\times i }[/math]

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

1. Выбрать произвольное [math]\displaystyle{ x_0 }[/math]; вычислить [math]\displaystyle{ r_0=b-Ax_0 }[/math] и [math]\displaystyle{ v_1=\frac{r_0}{\parallel r_0 \parallel} }[/math]

2. Задать матрицу [math]\displaystyle{ \overline{H}_m }[/math] размерности [math]\displaystyle{ (m+1)\times m }[/math] и всем ее элементам присвоить нулевые значения

3. Для [math]\displaystyle{ j=1,\dots,m }[/math] выполнять:

[math]\displaystyle{ \bullet \  }[/math] для [math]\displaystyle{ i=1,\dots,j }[/math] выполнять: [math]\displaystyle{ h_{i,j}=(Av_j,v_i) }[/math]
[math]\displaystyle{ \bullet \ \omega_j=Av_j-\sum_{i=1}^j h_{i,j}v_i }[/math]
[math]\displaystyle{ \bullet \ h_{j+1,j}={\parallel \omega_j \parallel} }[/math]
[math]\displaystyle{ \bullet \ v_{j+1}=\frac{\omega_j}{h_{j+1,j}} }[/math]

4. Вычислить [math]\displaystyle{ y_m=\arg\min_{y} \parallel \beta e_1 - \overline{H}_my \parallel,\ y \in \mathbb{R}^m }[/math] и [math]\displaystyle{ x_m=x_0+V_my_m }[/math]

5. Вычислить [math]\displaystyle{ r_m=f-Ax_m }[/math]. Если достигнутая точность удовлетворительна, остановиться.

Иначе принять [math]\displaystyle{ x_0=x_m }[/math] и [math]\displaystyle{ v_1=\frac{r_m}{\parallel r_m \parallel} }[/math] и вернуться к пункту 2.