Участник:Alexander34396/Обобщенный метод минимальных невязок: различия между версиями

Материал из Алговики
Перейти к навигации Перейти к поиску
Строка 86: Строка 86:
 
=== Схема реализации последовательного алгоритма ===
 
=== Схема реализации последовательного алгоритма ===
  
На практике используется перезапускаемая версия метода - GMRES(m):
+
==== Схема реализации классического GMRES ====
+
 
'''1. Подготовка перед итерационным процессом:'''
+
Для решения исходной системы методом GMRES можно воспользоваться следующим алгоритмом:
 +
 
 +
'''Подготовка перед итерационным процессом:'''
 
# Выбрать начальное приближение <math> x_0 </math>;
 
# Выбрать начальное приближение <math> x_0 </math>;
 
# Посчитать невязку <math> r_0 = b - Ax_0 </math>;
 
# Посчитать невязку <math> r_0 = b - Ax_0 </math>;
 
# Вычислить <math> v_1 = \frac{r_0}{\|r_0\|_2} </math>.
 
# Вычислить <math> v_1 = \frac{r_0}{\|r_0\|_2} </math>.
 +
 +
'''k-ая итерация:'''
 +
# <math> h_{ik} = (Av_k, v_i), \quad i=1,\ldots,k </math>;
 +
# <math> \hat{v}_{k+1} = Av_k - \sum_{i=1}^k h_{ik}v_{i} </math>;
 +
# <math> h_{k+1k} = \|\hat{v}_{k+1}\|_2 </math>;
 +
# <math> v_{k+1} = \frac{\hat{v}_{k+1}}{h_{k+1k}} </math>;
 +
# <math> x_k = x_0 + V_ky_k </math>, где <math>y_k</math> минимизирует <math>\|r_0 - AV_ky_k\|_2</math>.
 +
 +
==== Схема реализации перезапускаемого GMRES ====
 +
 +
На практике используется перезапускаемая версия метода - GMRES(m):
 +
 +
'''1 Подготовка перед итерационным процессом:'''
 +
:1.1 Выбрать начальное приближение <math> x_0 </math>;
 +
:1.2 Посчитать невязку <math> r_0 = b - Ax_0 </math>;
 +
:1.3 Вычислить <math> v_1 = \frac{r_0}{\|r_0\|_2} </math>.
  
 
'''2. Построение ортонормированного базиса <math> K_m </math>:'''
 
'''2. Построение ортонормированного базиса <math> K_m </math>:'''
  
 
: Для всех <math> j </math> от 1 до m по нарастанию выполнять:
 
: Для всех <math> j </math> от 1 до m по нарастанию выполнять:
# <math> h_{ij} := (Av_j, v_i), \quad i=1,\ldots,j </math>;
+
:2.1 <math> h_{ij} := (Av_j, v_i), \quad i=1,\ldots,j </math>;
# <math> \hat{v}_{j+1} := Av_j - \sum_{i=1}^j h_{ij}v_{i} </math>;
+
:2.2 <math> \hat{v}_{j+1} := Av_j - \sum_{i=1}^j h_{ij}v_{i} </math>;
# <math> h_{j+1j} = \|\hat{v}_{j+1}\|_2 </math>;
+
:2.3 <math> h_{j+1j} = \|\hat{v}_{j+1}\|_2 </math>;
# <math> v_{j+1} = \frac{\hat{v}_{j+1}}{h_{j+1j}} </math>.
+
:2.4 <math> v_{j+1} = \frac{\hat{v}_{j+1}}{h_{j+1j}} </math>.
  
'''3. Вычисление <math> x_m </math>:'''
+
'''3 Вычисление <math> x_m </math>:'''
# <math> x_m = x_0 + V_my_m </math>, где <math>y_m</math> минимизирует <math>\|r_0 - AV_my_m\|_2</math>;
+
:3.1 <math> x_m = x_0 + V_my_m </math>, где <math>y_m</math> минимизирует <math>\|r_0 - AV_my_m\|_2</math>;
# Вычислить <math> r_m </math>;
+
:3.2 Вычислить <math> r_m </math>;
# Если требуемая точность достигнута, остановиться.
+
:3.3 Если требуемая точность достигнута, остановиться.
  
'''4. Перезапуск:'''
+
'''4 Перезапуск:'''
# <math>x_0 = x_m</math>;
+
:4.1 <math>x_0 = x_m</math>;
# <math>v_1 = \frac{r_m}{\|r_m\|_2}</math>;
+
:4.2 <math>v_1 = \frac{r_m}{\|r_m\|_2}</math>;
# Перейти к шагу 2.
+
:4.3 Перейти к шагу 2.
  
 
=== Последовательная сложность алгоритма ===
 
=== Последовательная сложность алгоритма ===

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

Symbol wait.svgЭта работа ждет рассмотрения преподавателем
Дата последней правки страницы:
17.10.2016
Авторы этой статьи считают, что задание выполнено.

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

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

Обобщённый метод минимальных невязок (англ. Generalized minimal residual method, GMRES) - итерационный метод численного решения системы линейных алгебраических уравнений с невырожденной матрицей. Метод основан на минимизации квадратичного функционала невязки на подпространствах Крылова. Разработан Юсефом Саадом и Мартином Шульцем[1] в 1986 году как обобщение метода MINRES[2] на случай систем с несимметричными матрицами.

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

Исходные данные:

  • система линейных алгебраических уравнений вида [math] Ax = b [/math], где [math] A [/math] — невырожденная матрица размера [math]n[/math]-на-[math] n [/math];
  • [math] x_0 [/math] - начальное приближение.

Вычисляемые данные:

  • [math] x_k [/math] - приближённое решение исходной системы.

Метод GMRES приближает точное решение исходной системы [math] Ax = b [/math] вектором [math] x_k [/math], который минимизирует [math] l_2 [/math]-норму невязки [math] r_k = Ax_k - b[/math] на подпространстве Крылова:

[math] K_k = K_k(A,b) = \operatorname{span} \, \{ b, Ab, A^2b, \ldots, A^{k-1}b \} [/math].

На каждой итерации метода решение уточняется поправкой, представленной в виде разложения по [math] l_2 [/math]-ортонормированному базису пространства [math] K_k [/math]:

[math] x_k = x_0 + z_k [/math], где [math] x_0 [/math] - некоторое начальное приближение, [math] z_k \in K_k [/math] - поправка решения.

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

Для построения ортонормированного базиса [math] K_k [/math] метод использует процесс Арнольди.

Исходные данные:

  • [math] v_1 [/math], такой что [math] \|v_1\|_2 = 1 [/math];
  • матрица A размером [math] n [/math]-на-[math] n [/math].

Формула процесса:

[math] h_{j+1j}v_{j+1} = Av_j - \sum_{i=1}^j h_{ij}v_{j} [/math], где [math] h_{ij} = (Av_j, v_i)[/math], [math] \quad j=1,\ldots,k [/math].

Алгоритм процесса:

Выполнять для [math] \quad j=1,\ldots,k [/math]:

  1. [math] h_{ij} = (Av_j, v_i), \quad i=1,\ldots,j [/math];
  2. [math] \hat{v}_{j+1} = Av_j - \sum_{i=1}^j h_{ij}v_{i} [/math];
  3. [math] h_{j+1j} = \|\hat{v}_{j+1}\|_2 [/math];
  4. [math] v_{j+1} = \frac{\hat{v}_{j+1}}{h_{j+1j}} [/math].

При введении матричных обозначений можно записать:

  1. [math] z_k = V_ky_k [/math], где [math] y_k \in \mathbb{R}^k [/math] - вектор коэффициентов;
  2. [math] AV_k = V_kH_k + h_{k+1k}v_{k+1} e_k^T = V_{k+1}\overline{H}_k [/math], где [math] V_k = [v_1|v_2|...|v_k] [/math], а [math] \overline{H}_k [/math] - верхняя матрица Хессенберга размерности [math] (j+1) [/math] на [math] j [/math].

1.2.2 Минимизация функционала невязки

Для решения исходной системы GMRES вычисляет приближённое решение [math] x_k [/math], которое минимизизует функционал невязки:

[math] J(y) = \|b- A x_k\|_2 = \|b - A(x_0 + V_ky)\|_2 = \|\beta e_1 - \overline{H}_ky_k \|_2 [/math], где [math] \beta = \|r_0\|_2 [/math].

Вектор [math] y_k [/math] может найден как решение линейной задачи наименьших квадратов размера [math] k+1 [/math]-на-[math] k [/math]:

[math] y_k = \arg\min_{y}\| \beta e_1 - \overline{H}_ky\|_2 [/math].

Для решения задачи матрица [math] \overline{H}_k [/math] приводится к верхнему треугольному виду с помощью [math] k [/math] последовательных вращений Гивенса.

1.2.3 Общая схема метода

В общем виде k-aя итерация алгоритма GMRES может быть записан следующим образом:

  1. найти ортонормированный базис [math] V_k [/math] подпространства [math] K_k [/math] с помощью ортогонализации Арнольди;
  2. найти [math] y_k [/math], минимизирующий [math] \|r_k\|_2 [/math];
  3. вычислить [math] x_k = x_0 + V_ky_k [/math];
  4. вычислить [math] r_k [/math] и остановиться,если требуемая точность была достигнута, иначе повторить.

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

Вычислительное ядро последовательной версии метода GMRES состоит из двух частей:

  • Вычисление ортонормированного базиса [math] K_k [/math];
  • Формирование приближенного решения [math] x_k [/math].

Для вычисления ортонормированного базиса [math] K_k [/math] метод использует процесс Арнольди:

[math] h_{j+1j}v_{j+1} = Av_j - \sum_{i=1}^j h_{ij}v_{j} [/math], где [math] h_{ij} = (Av_j, v_i)[/math], [math] \quad j=1,\ldots,k [/math].

Для нахождения приближённого решения [math] x_k [/math] метод использует формулу:

[math] x_k = x_0 + V_ky_k [/math].

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

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

  • Умножение матрицы на вектор;
  • Вычисление скалярного произведения векторов;
  • Вычисление [math] l_2 [/math]-нормы вектора;
  • Умножение вектора на скаляр;
  • Деление вектора на скаляр.

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

1.5.1 Схема реализации классического GMRES

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

Подготовка перед итерационным процессом:

  1. Выбрать начальное приближение [math] x_0 [/math];
  2. Посчитать невязку [math] r_0 = b - Ax_0 [/math];
  3. Вычислить [math] v_1 = \frac{r_0}{\|r_0\|_2} [/math].

k-ая итерация:

  1. [math] h_{ik} = (Av_k, v_i), \quad i=1,\ldots,k [/math];
  2. [math] \hat{v}_{k+1} = Av_k - \sum_{i=1}^k h_{ik}v_{i} [/math];
  3. [math] h_{k+1k} = \|\hat{v}_{k+1}\|_2 [/math];
  4. [math] v_{k+1} = \frac{\hat{v}_{k+1}}{h_{k+1k}} [/math];
  5. [math] x_k = x_0 + V_ky_k [/math], где [math]y_k[/math] минимизирует [math]\|r_0 - AV_ky_k\|_2[/math].

1.5.2 Схема реализации перезапускаемого GMRES

На практике используется перезапускаемая версия метода - GMRES(m):

1 Подготовка перед итерационным процессом:

1.1 Выбрать начальное приближение [math] x_0 [/math];
1.2 Посчитать невязку [math] r_0 = b - Ax_0 [/math];
1.3 Вычислить [math] v_1 = \frac{r_0}{\|r_0\|_2} [/math].

2. Построение ортонормированного базиса [math] K_m [/math]:

Для всех [math] j [/math] от 1 до m по нарастанию выполнять:
2.1 [math] h_{ij} := (Av_j, v_i), \quad i=1,\ldots,j [/math];
2.2 [math] \hat{v}_{j+1} := Av_j - \sum_{i=1}^j h_{ij}v_{i} [/math];
2.3 [math] h_{j+1j} = \|\hat{v}_{j+1}\|_2 [/math];
2.4 [math] v_{j+1} = \frac{\hat{v}_{j+1}}{h_{j+1j}} [/math].

3 Вычисление [math] x_m [/math]:

3.1 [math] x_m = x_0 + V_my_m [/math], где [math]y_m[/math] минимизирует [math]\|r_0 - AV_my_m\|_2[/math];
3.2 Вычислить [math] r_m [/math];
3.3 Если требуемая точность достигнута, остановиться.

4 Перезапуск:

4.1 [math]x_0 = x_m[/math];
4.2 [math]v_1 = \frac{r_m}{\|r_m\|_2}[/math];
4.3 Перейти к шагу 2.

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

Если пренебречь сложностью вычисления [math] y_k [/math], то общую сложность алгоритма GMRES можно разделить на две части:

1. Сложность вычисления ортонормированного базиса пространства [math] K_k [/math]:

a) Для вычисления [math] j [/math]-го вектора базиса [math] K_k, j \lt k[/math] требуется:
[math] NZ + (2j + 1)n [/math] мультипликативных операций, где [math] NZ [/math] - количество ненулевых элементов матрицы [math] A [/math].
b) Вычисление последнего вектора базиса требует:
[math] n(k + 1) [/math] мультипликативных операций.
c) Общая мультипликативная сложность вычисления ортонормированного базиса [math] K_k [/math]:
[math] nk(k + 1) + kNZ [/math].

2. Сложность вычисления приближённого решения [math] x_k [/math]:

Вычисление этой формулы требует [math] nk [/math] мультипликативных операций.

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

2 Литература

<references \>

  1. 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)
  2. C. C. Paige and M. A. Saunders. Solution of sparse indefinite systems of linear equations, SIAM J. Numerical Analysis 12, 617-629 (1975)