Итерационные методы решения СЛАУ: Якоби, Зейдель и релаксация

Метод простых итераций (метод Якоби)

Исторически первым и самым простым алгоритмом является метод Карла Якоби. Его логика предельно ясна: мы берем i-е уравнение системы и выражаем из него i-ю неизвестную через все остальные. В матричном виде это означает расщепление матрицы A на диагональную часть (D), нижнюю (L) и верхнюю (U). Метод Якоби использует значения неизвестных, вычисленные на предыдущем шаге (k), для нахождения новых значений на текущем шаге (k+1). На компьютерах этот алгоритм работает фантастически просто: он сводится к серии умножений разреженной матрицы на вектор. Главное достоинство метода Якоби — его идеальная распараллеливаемость. Каждая неизвестная на новом шаге может вычисляться на отдельном ядре процессора или потоке видеокарты совершенно независимо от других.

Метод Гаусса-Зейделя: ускорение сходимости

Метод Якоби надежен, но сходится к ответу мучительно медленно. Математики Карл Гаусс и Филипп Зейдель предложили элегантную модификацию. При расчете значений на шаге (k+1) мы вычисляем переменные по очереди (x1, затем x2, затем x3...). В методе Зейделя мы не ждем окончания полного цикла итерации, а сразу же используем только что найденное новое значение x1 для вычисления x2. Свежая информация мгновенно запускается в работу. Алгебраически это означает, что для обращения используется не только диагональная матрица D, но и нижнетреугольная (D+L). Это незначительное изменение алгоритма увеличивает скорость сходимости метода практически в два раза по сравнению с методом Якоби.

Условие сходимости: диагональное преобладание

Итерационные методы (в отличие от Гаусса) не гарантируют нахождения решения для любой произвольной невырожденной матрицы. Они могут разойтись, и тогда числа уйдут в бесконечность. Достаточным (но не необходимым) условием сходимости методов Якоби и Зейделя является наличие у матрицы коэффициентов свойства строгого диагонального преобладания по строкам или столбцам. Это означает, что модуль элемента, стоящего на главной диагонали в каждой строке, должен быть строго больше суммы модулей всех остальных элементов этой же строки. К счастью для инженеров, уравнения математической физики (например, дискретизация уравнений Пуассона или Навье-Стокса на расчетных сетках) естественным образом порождают матрицы, обладающие именно таким диагональным преобладанием.

Вершина классики: метод последовательной верхней релаксации (SOR)

В 1950-х годах Дэвид Янг разработал метод последовательной верхней релаксации (Successive Over-Relaxation, SOR), который стал прорывом в вычислительной математике. Идея SOR заключается в том, чтобы не просто брать новое значение из метода Зейделя, а экстраполировать его, «перешагивая» чуть дальше в предполагаемом направлении решения. В алгоритм вводится параметр релаксации omega (обычно в диапазоне от 1 до 2). При правильном (оптимальном) выборе этого параметра метод SOR способен сходиться в десятки и сотни раз быстрее обычного Гаусса-Зейделя. Сегодня классические стационарные методы Якоби и SOR редко используются в одиночку; чаще всего они выступают в роли «предобуславливателей» (preconditioners) для еще более мощных современных алгоритмов в подпространствах Крылова, таких как метод сопряженных градиентов (CG).