Участник:Alexander34396/Обобщенный метод минимальных невязок
Эта работа ждет рассмотрения преподавателем Дата последней правки страницы: 18.10.2016 Авторы этой статьи считают, что задание выполнено. |
Основные авторы описания: А.А.Ульянов (разделы 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.9, 2.7), К.Петров (разделы 1.7, 1.8, 1.10, 2.4).
Обобщенный метод минимальных невязок | |
Последовательный алгоритм | |
Последовательная сложность | [math] O(nm^2) [/math] |
Объём входных данных | [math] (n + 1)^2 [/math] |
Объём выходных данных | [math] n [/math] |
Содержание
- 1 Свойства и структура алгоритма
- 1.1 Общее описание алгоритма
- 1.2 Математическое описание алгоритма
- 1.3 Вычислительное ядро алгоритма
- 1.4 Макроструктура алгоритма
- 1.5 Схема реализации последовательного алгоритма
- 1.6 Последовательная сложность алгоритма
- 1.7 Информационный граф
- 1.8 Ресурс параллелизма алгоритма
- 1.9 Входные и выходные данные алгоритма
- 1.10 Свойства алгоритма
- 2 Программная реализация алгоритма
- 2.1 Особенности реализации последовательного алгоритма
- 2.2 Локальность данных и вычислений
- 2.3 Возможные способы и особенности параллельной реализации алгоритма
- 2.4 Масштабируемость алгоритма и его реализации
- 2.5 Динамические характеристики и эффективность реализации алгоритма
- 2.6 Выводы для классов архитектур
- 2.7 Существующие реализации алгоритма
- 3 Литература
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];
- матрица [math] A [/math] размером [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]:
- [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];
- [math] h_{j+1j} = \|\hat{v}_{j+1}\|_2 [/math];
- [math] v_{j+1} = \frac{\hat{v}_{j+1}}{h_{j+1j}} [/math].
При введении матричных обозначений можно записать:
- [math] z_k = V_ky_k [/math], где [math] y_k \in \mathbb{R}^k [/math] - вектор коэффициентов;
- [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 может быть записан следующим образом:
- найти ортонормированный базис [math] V_k [/math] подпространства [math] K_k [/math] с помощью ортогонализации Арнольди;
- найти [math] y_k [/math], минимизирующий [math] \|r_k\|_2 [/math];
- вычислить [math] x_k = x_0 + V_ky_k [/math], где [math] x_0 [/math] - начальное приближение;
- вычислить [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 Схема реализации последовательного алгоритма
На практике используется перезапускаемая версия метода - 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_m [/math], то общую сложность алгоритма GMRES можно разделить на две части:
1. Сложность вычисления ортонормированного базиса пространства [math] K_m [/math]:
- a) Для вычисления [math] j [/math]-го вектора базиса [math] K_m, j \lt m [/math] требуется:
- [math] NZ + (2j + 1)n [/math] мультипликативных операций, где [math] NZ [/math] - количество ненулевых элементов матрицы [math] A [/math].
- b) Вычисление последнего вектора базиса требует:
- [math] n(m + 1) [/math] мультипликативных операций.
- c) Общая мультипликативная сложность вычисления ортонормированного базиса [math] K_m [/math]:
- [math] O(nm^2) [/math].
2. Сложность вычисления приближённого решения [math] x_m [/math]:
- a) Вычисление этой формулы требует [math] nm [/math] мультипликативных операций.
1.7 Информационный граф
Здесь input - входные данные алгоритма (квадратная матрица, вектор правой части, начальное приближение, точность решения);
[math] A_i [/math] - шаги алгоритма Арнольди, [math] x_i [/math] - нахождение невязки на i-той итерации;
if - вычисление [math] r_k [/math] и в случае достижения нужной точности - остановка.
При этом, шаги алгоритма Арнольди выполняются строго последовательно, но вычисление i-той невязки может происходить параллельно с (i+1)-м шагом алгоритма Арнольди.
1.8 Ресурс параллелизма алгоритма
1.9 Входные и выходные данные алгоритма
Входные данные:
- [math] A [/math] – квадратная матрица размера [math] n [/math]-на-[math] n [/math];
- [math] b [/math] – вектор правой части размера [math] n [/math];
- [math] x_0 [/math] – начальное приближение размера [math] n [/math];
- [math] \epsilon [/math] – погрешность решения.
Для перезапускаемой версии алгоритма дополнительно:
- [math] m [/math] – число итераций.
Объём входных данных: [math] n^2 + 2n + 1 [/math].
Выходные данные:
- [math] x_k [/math] - вектор приближённого решения размера [math] n [/math].
Объём выходных данных: [math] n [/math].
1.10 Свойства алгоритма
2 Программная реализация алгоритма
2.1 Особенности реализации последовательного алгоритма
2.2 Локальность данных и вычислений
2.3 Возможные способы и особенности параллельной реализации алгоритма
2.4 Масштабируемость алгоритма и его реализации
2.5 Динамические характеристики и эффективность реализации алгоритма
2.6 Выводы для классов архитектур
2.7 Существующие реализации алгоритма
Реализации алгоритма GMRES присутствуют во многих программных пакетах:
3 Литература
<references \>
- ↑ 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)
- ↑ C. C. Paige and M. A. Saunders. Solution of sparse indefinite systems of linear equations, SIAM J. Numerical Analysis 12, 617-629 (1975)