Алгоритм Томаса

Алгоритм Томаса (также известный как метод прогонки для трёхдиагональных матриц) — это специализированный прямой метод решения систем линейных алгебраических уравнений (СЛАУ) с трёхдиагональной матрицей. Предложен британским математиком Ллевеллином Томасом в 1949 году. Алгоритм является вариантом метода исключения Гаусса, адаптированным для матриц ленточной структуры, и отличается высокой вычислительной эффективностью и низкими требованиями к памяти.

История

Метод был разработан Ллевеллином Томасом (Llewellyn Thomas, 1903–1992) в контексте решения задач численного анализа, возникающих при моделировании физических процессов. Первоначально алгоритм применялся для решения дифференциальных уравнений в частных производных, в частности, при расчётах теплопроводности и гидродинамики. В 1949 году Томас опубликовал работу, в которой формализовал метод, однако широкое распространение алгоритм получил в 1950–1960-е годы с развитием вычислительной техники. В СССР и России метод часто называют «методом прогонки» (или «алгоритмом прогонки»), что связано с характерной последовательностью вычислений — «прогонкой» коэффициентов.

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

Постановка задачи

Рассматривается система линейных уравнений вида:

\[ a_i x_{i-1} + b_i x_i + c_i x_{i+1} = d_i, \quad i = 1, 2, \dots, n \]

где:

• \( a_i \) — поддиагональные элементы (для \( i = 1 \) полагается \( a_1 = 0 \)),
• \( b_i \) — диагональные элементы,
• \( c_i \) — наддиагональные элементы (для \( i = n \) полагается \( c_n = 0 \)),
• \( d_i \) — элементы правой части,
• \( x_i \) — искомые неизвестные.

Матрица системы является трёхдиагональной, то есть ненулевые элементы расположены только на главной диагонали и двух соседних с ней.

Прямой ход (прямая прогонка)

Цель прямого хода — исключить поддиагональные элементы \( a_i \), преобразовав систему к виду с двухдиагональной матрицей. Для этого вычисляются прогоночные коэффициенты \( P_i \) и \( Q_i \):

\[ P_1 = -\frac{c_1}{b_1}, \quad Q_1 = \frac{d_1}{b_1} \]

Для \( i = 2, 3, \dots, n \):

\[ P_i = -\frac{c_i}{b_i + a_i P_{i-1}}, \quad Q_i = \frac{d_i - a_i Q_{i-1}}{b_i + a_i P_{i-1}} \]

Обратный ход (обратная прогонка)

После завершения прямого хода неизвестные находятся последовательно с конца:

\[ x_n = Q_n \]

Для \( i = n-1, n-2, \dots, 1 \):

\[ x_i = P_i x_{i+1} + Q_i \]

Условия применимости

Алгоритм Томаса корректен, если в процессе вычислений не происходит деления на ноль. Достаточным условием устойчивости является диагональное преобладание матрицы:

\[ |b_i| \geq |a_i| + |c_i|, \quad i = 1, 2, \dots, n \]

причём хотя бы для одного \( i \) неравенство строгое. В этом случае алгоритм численно устойчив. Если условие не выполняется, метод может давать неверные результаты или требовать перестановки строк (перестановки не нарушают трёхдиагональную структуру, но усложняют реализацию).

Вычислительная сложность

Алгоритм Томаса имеет линейную сложность \( O(n) \) по числу арифметических операций, что значительно эффективнее общего метода Гаусса (\( O(n^3) \)) для полных матриц. Для решения системы из \( n \) уравнений требуется примерно \( 8n \) операций (умножений и сложений). Объём памяти — \( O(n) \), так как хранятся только три диагонали и вектор правой части.

Применение

Численное решение дифференциальных уравнений

Алгоритм Томаса широко применяется при решении краевых задач для обыкновенных дифференциальных уравнений (ОДУ) и дифференциальных уравнений в частных производных (ДУЧП) методом конечных разностей. Например, при дискретизации одномерного уравнения теплопроводности или уравнения Пуассона на равномерной сетке возникает трёхдиагональная система.

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

В задачах газовой динамики и гидродинамики (например, при расчёте течений в трубах или каналах) метод прогонки используется для решения систем, возникающих при аппроксимации уравнений Навье — Стокса.

Интерполяция кубическими сплайнами

При построении кубических сплайнов для интерполяции данных возникает система с трёхдиагональной матрицей, которая решается методом Томаса.

Обработка сигналов и изображений

В задачах фильтрации и сглаживания (например, при реализации фильтра Гаусса или медианного фильтра) могут возникать трёхдиагональные системы, решаемые алгоритмом Томаса.

Пример

Рассмотрим систему из 4 уравнений:

\[ \begin{cases} 2x_1 + 3x_2 = 8 \\ x_1 + 4x_2 + x_3 = 10 \\ x_2 + 5x_3 + 2x_4 = 15 \\ 3x_3 + 6x_4 = 12 \end{cases} \]

Матрица коэффициентов: \[ A = \begin{pmatrix} 2 & 3 & 0 & 0 \\ 1 & 4 & 1 & 0 \\ 0 & 1 & 5 & 2 \\ 0 & 0 & 3 & 6 \end{pmatrix}, \quad d = \begin{pmatrix} 8 \\ 10 \\ 15 \\ 12 \end{pmatrix} \]

Прямой ход:

• \( P_1 = -3/2 = -1.5 \), \( Q_1 = 8/2 = 4 \)
• \( i=2 \): \( P_2 = -1 / (4 + 1 \cdot (-1.5)) = -1 / 2.5 = -0.4 \), \( Q_2 = (10 - 1 \cdot 4) / 2.5 = 6 / 2.5 = 2.4 \)
• \( i=3 \): \( P_3 = -2 / (5 + 1 \cdot (-0.4)) = -2 / 4.6 \approx -0.43478 \), \( Q_3 = (15 - 1 \cdot 2.4) / 4.6 = 12.6 / 4.6 \approx 2.73913 \)
• \( i=4 \): \( P_4 = 0 \) (так как \( c_4 = 0 \)), \( Q_4 = (12 - 3 \cdot 2.73913) / (6 + 3 \cdot (-0.43478)) = (12 - 8.21739) / (6 - 1.30434) = 3.78261 / 4.69566 \approx 0.80556 \)

Обратный ход:

• \( x_4 = Q_4 \approx 0.80556 \)
• \( x_3 = P_3 \cdot x_4 + Q_3 \approx (-0.43478) \cdot 0.80556 + 2.73913 \approx -0.35000 + 2.73913 = 2.38913 \)
• \( x_2 = P_2 \cdot x_3 + Q_2 \approx (-0.4) \cdot 2.38913 + 2.4 = -0.95565 + 2.4 = 1.44435 \)
• \( x_1 = P_1 \cdot x_2 + Q_1 \approx (-1.5) \cdot 1.44435 + 4 = -2.16653 + 4 = 1.83347 \)

Проверка подстановкой в исходную систему подтверждает корректность решения.

Варианты и обобщения

Метод встречных прогонок

Используется для решения систем с трёхдиагональной матрицей, когда требуется повышенная устойчивость. Вычисления ведутся одновременно с двух сторон — от начала и конца системы.

Метод циклической прогонки

Применяется для систем с циклической трёхдиагональной матрицей (когда \( a_1 \neq 0 \) и \( c_n \neq 0 \)). Требует дополнительных вычислений для учёта замыкания.

Метод матричной прогонки

Обобщение на случай, когда элементы матрицы являются не числами, а блоками (например, матрицами меньшего размера). Используется при решении многомерных задач методом переменных направлений.

Критика и ограничения

• Алгоритм не применим к системам с произвольной разреженной структурой — только к трёхдиагональным.
• При нарушении диагонального преобладания возможна численная неустойчивость.
• Для систем с плохо обусловленной матрицей (например, с большим разбросом коэффициентов) метод может давать большую погрешность.
• В параллельных вычислениях алгоритм плохо масштабируется из-за последовательной природы прогонки.

Интересные факты

• Алгоритм Томаса является частным случаем метода исключения Гаусса для ленточных матриц, но был разработан независимо.
• В русскоязычной литературе метод часто называют «методом прогонки», а сам термин «прогонка» ввёл академик А. Н. Тихонов.
• Алгоритм используется в библиотеках численного анализа, таких как LAPACK, SciPy, MATLAB, а также в специализированных пакетах для моделирования физических процессов.

Источники

• Thomas, L. H. (1949). "Elliptic Problems in Linear Difference Equations over a Network". Watson Scientific Computing Laboratory, Columbia University.
• Самарский А. А., Гулин А. В. «Численные методы». — М.: Наука, 1989.
• Бахвалов Н. С., Жидков Н. П., Кобельков Г. М. «Численные методы». — М.: Бином, 2008.
• Press, W. H., Teukolsky, S. A., Vetterling, W. T., Flannery, B. P. "Numerical Recipes: The Art of Scientific Computing". — Cambridge University Press, 2007.