Полный метод циклической редукции

Циклическая редукция для трёхдиагональной матрицы, точечный вариант
Последовательный алгоритм
Последовательная сложность	[math]\displaystyle{ O(n) }[/math]
Объём входных данных	[math]\displaystyle{ 4n-2 }[/math]
Объём выходных данных	[math]\displaystyle{ n }[/math]
Параллельный алгоритм
Высота ярусно-параллельной формы	[math]\displaystyle{ O(log n) }[/math]
Ширина ярусно-параллельной формы	[math]\displaystyle{ O(n) }[/math]

Основные авторы описания: А.В.Фролов.

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

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

Циклическая редукция - один из вариантов метода исключения неизвестных в приложении к решению СЛАУ вида [math]\displaystyle{ Ax = b }[/math], где

[math]\displaystyle{ A = \begin{bmatrix} a_{11} & a_{12} & 0 & \cdots & \cdots & 0 \\ a_{21} & a_{22} & a_{23}& \cdots & \cdots & 0 \\ 0 & a_{32} & a_{33} & \cdots & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots & \ddots & 0 \\ 0 & \cdots & \cdots & a_{n-1 n-2} & a_{n-1 n-1} & a_{n-1 n} \\ 0 & \cdots & \cdots & 0 & a_{n n-1} & a_{n n} \\ \end{bmatrix}, x = \begin{bmatrix} x_{1} \\ x_{2} \\ \vdots \\ x_{n} \\ \end{bmatrix}, b = \begin{bmatrix} b_{1} \\ b_{2} \\ \vdots \\ b_{n} \\ \end{bmatrix} }[/math]

Бывает, однако, что при изложении сути методов решения трёхдиагональных СЛАУ элементы правой части и матрицы системы обозначают и нумеруют по-другому, например СЛАУ может иметь вид

[math]\displaystyle{ A = \begin{bmatrix} b_{1} & a_{1} & 0 & \cdots & \cdots & 0 \\ c_{2} & b_{2} & a_{2} & \cdots & \cdots & 0 \\ 0 & c_{3} & b_{3} & \cdots & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots & \ddots & 0 \\ 0 & \cdots & \cdots & c_{n-1} & b_{n-1} & a_{n-1} \\ 0 & \cdots & \cdots & 0 & c_{n} & b_{n} \\ \end{bmatrix}\begin{bmatrix} x_{1} \\ x_{2} \\ \vdots \\ x_{n} \\ \end{bmatrix} = \begin{bmatrix} f_{1} \\ f_{2} \\ \vdots \\ f_{n} \\ \end{bmatrix} }[/math]

или, если записывать отдельно по уравнениям, то

[math]\displaystyle{ b_{1} x_{1} + a_{1} x_{2} = f_{1} }[/math],

[math]\displaystyle{ c_{i} x_{i-1} + b_{i} x_{i} + a_{i} x_{i+1} = f_{i}, 2 \le i \le n-1 }[/math],

[math]\displaystyle{ c_{n} x_{n-1} + b_{n} x_{n} = f_{n} }[/math].

Циклическая редукция, как и все варианты прогонки, заключается в исключении из уравнений неизвестных, однако, в отличие от них, в ней исключение ведут одновременно по всей СЛАУ. В принципе, её можно считать вариантом метода редукции, выполняемого максимально возможное для данной СЛАУ число раз.

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

Лучше всего схема циклической редукции разработана для случая [math]\displaystyle{ n = 2^{m}-1 }[/math]. Эта схема состоит из прямого и обратного ходов. Прямой ход состоит из последовательного уменьшения в СЛАУ количества уравнений почти в 2 раза (за счёт подстановки из уравнений с нечётными номерами заменяются уравнения с чётными), пока не останется одно уравнение, обратный - в получении всё большего количества компонент решения исходной СЛАУ. Оба хода - как прямой, так и обратный - разбиты на шаги. Здесь мы приведём тот вариант алгоритма, в котором операции экономятся за счёт предварительной нормировки уравнений, используемых для исключения неизвестных.

Прямой ход редукции

В начале считается, что все [math]\displaystyle{ c^{(0)}_{i} = c_{i}, b^{(0)}_{i} = b_{i}, a^{(0)}_{i} = a_{i}, f^{(0)}_{i} = f_{i}, x^{(0)}_{i} = x_{i} }[/math]

На k-м шаге теперь выполняем процедуру редукции системы уравнений размерности n.

Для каждого из уравнений

[math]\displaystyle{ c^{(k)}_{i} x^{(k)}_{i-1} + b^{(k)}_{i} x^{(k)}_{i} + a^{(k)}_{i} x^{(k)}_{i+1} = f^{(k)}_{i} }[/math]

с нечётными [math]\displaystyle{ i }[/math] с помощью деления уравнения на [math]\displaystyle{ b^{(k)}_{i} }[/math] выполняется его замена на уравнение

[math]\displaystyle{ (c^{(k)}_{i}/b^{(k)}_{i}) x^{(k)}_{i-1} + x^{(k)}_{i} + (a^{(k)}_{i}/b^{(k)}_{i}) x^{(k)}_{i+1} = f^{(k)}_{i}/b^{(k)}_{i} }[/math]

Рисунок 1. Микрограф "узла" подготовительного шага прямого хода алгоритма циклической редукции

Уравнение

[math]\displaystyle{ b^{(k)}_{1} x^{(k)}_{1} + a^{(k)}_{1} x^{(k)}_{2} = f^{(k)}_{1} }[/math]

аналогично меняется на уравнение

[math]\displaystyle{ x^{(k)}_{1} + (a^{(k)}_{1}/b^{(k)}_{1}) x^{(k)}_{2} = f^{(k)}_{1}/b^{(k)}_{1} }[/math]:

а уравнение

[math]\displaystyle{ c^{(k)}_{n} x^{(k)}_{n-1} + b^{(k)}_{n} x^{(k)}_{n} = f^{(k)}_{n} }[/math]

меняется на уравнение

[math]\displaystyle{ (c^{(k)}_{n}/b^{(k)}_{n}) x^{(k)}_{n-1} + x^{(k)}_{n} = f^{(k)}_{n}/b^{(k)}_{n} }[/math]:

Рисунок 2a. Микрограф "узла" с нечётным номером прямого хода алгоритма циклической редукции

Для каждого же из уравнений

[math]\displaystyle{ c^{(k)}_{i} x^{(k)}_{i-1} + b^{(k)}_{i} x^{(k)}_{i} + a^{(k)}_{i} x^{(k)}_{i+1} = f^{(k)}_{i} }[/math]

с чётными [math]\displaystyle{ i }[/math] (кроме [math]\displaystyle{ 2 }[/math] и [math]\displaystyle{ n-2 }[/math]) выполняется, с учётом [math]\displaystyle{ x^{(k+1)}_{i/2} = x^{(k)}_{i} }[/math] его замена на уравнение

[math]\displaystyle{ c^{(k+1)}_{i/2} x^{(k+1)}_{(i-2)/2} + b^{(k+1)}_{i/2} x^{(k+1)}_{i/2} + a^{(k+1)}_{i/2} x^{(k+1)}_{(i+2)/2} = f^{(k+1)}_{i/2} }[/math]

Рисунок 2b. Микрограф "узла" с чётным номером прямого хода алгоритма циклической редукции

при этом

[math]\displaystyle{ c^{(k+1)}_{i/2} = - c^{(k)}_{i}(c^{(k)}_{i-1}/b^{(k)}_{i-1}) }[/math],

[math]\displaystyle{ a^{(k+1)}_{i/2} = - a^{(k)}_{i}(a^{(k)}_{i+1}/b^{(k)}_{i+1}) }[/math],

[math]\displaystyle{ b^{(k+1)}_{i/2} = b^{(k)}_{i} - c^{(k)}_{i}(a^{(k)}_{i-1}/b^{(k)}_{i-1}) - a^{(k)}_{i}(c^{(k)}_{i+1}/b^{(k)}_{i+1}) }[/math],

[math]\displaystyle{ f^{(k+1)}_{i/2} = f^{(k)}_{i} - c^{(k)}_{i}f^{(k)}_{i-1}/b^{(k)}_{i-1} - a^{(k)}_{i}f^{(k)}_{i-1}/b^{(k)}_{i-1} }[/math].

Для 2го уравнения выполняется его замена на уравнение

[math]\displaystyle{ b^{(k+1)}_{1} x^{(k+1)}_{1} + a^{(k+1)}_{1} x^{(k+1)}_{(2} = f^{(k+1)}_{1} }[/math]

при этом

[math]\displaystyle{ a^{(k+1)}_{1} = - a^{(k)}_{2}(a^{(k)}_{3}/b^{(k)}_{3}) }[/math],

[math]\displaystyle{ b^{(k+1)}_{1} = b^{(k)}_{2} - c^{(k)}_{2}(a^{(k)}_{1}/b^{(k)}_{1}) - a^{(k)}_{2}(c^{(k)}_{3}/b^{(k)}_{3}) }[/math],

[math]\displaystyle{ f^{(k+1)}_{1} = f^{(k)}_{2} - c^{(k)}_{2}f^{(k)}_{1}/b^{(k)}_{1} - a^{(k)}_{2}f^{(k)}_{1}/b^{(k)}_{1} }[/math]

[math]\displaystyle{ n-1 }[/math]-е уравнение заменяется на

[math]\displaystyle{ c^{(k+1)}_{(n-1)/2} x^{(k+1)}_{(n-3)/2} + b^{(k+1)}_{(n-1)/2} x^{(k+1)}_{(n-1)/2} = f^{(k+1)}_{(n-1)/2} }[/math]

при этом

[math]\displaystyle{ c^{(k+1)}_{(n-1)/2} = - c^{(k)}_{n-1}(c^{(k)}_{n-2}/b^{(k)}_{n-2}) }[/math],

[math]\displaystyle{ b^{(k+1)}_{(n-1)/2} = b^{(k)}_{n-1} - c^{(k)}_{n-1}(a^{(k)}_{n-2}/b^{(k)}_{n-2}) - a^{(k)}_{n-1}(c^{(k)}_{n}/b^{(k)}_{n}) }[/math],

[math]\displaystyle{ f^{(k+1)}_{(n-1)/2} = f^{(k)}_{n-1} - c^{(k)}_{n-1}f^{(k)}_{n-2}/b^{(k)}_{n-2} - a^{(k)}_{n-1}f^{(k)}_{n-2}/b^{(k)}_{n-2} }[/math].

По окончании всех этих манипуляций размерность k+1-й СЛАУ оказывается равной [math]\displaystyle{ (n-1)/2 }[/math].

Шаги повторяются до тех пор, пока после [math]\displaystyle{ m-1 }[/math] шагов редукции размерность СЛАУ не становится равной 1 и остаётся одно уравнение

[math]\displaystyle{ b^{(m-1)}_{1} x^{(m-1)}_{1} = f^{(m-1)}_{1} }[/math]

Обратный ход редукции

Рисунок 3. Микрограф "узла" обратного хода алгоритма циклической редукции

Из последнего уравнения, полученного прямым ходом, вычисляется

[math]\displaystyle{ x^{(m-1)}_{1} = f^{(m-1)}_{1}/b^{(m-1)}_{1} }[/math]

Теперь, последовательно уменьшая верхние индексы неизвестных, используется нечётные уравнения каждого шага для вычисления неизвестных с соотвествующими нечётными номерами. Чётные неизвестные получаются из тождеств [math]\displaystyle{ x^{(k)}_{i} = x^{(k+1)}_{i/2} }[/math], а для нечётных [math]\displaystyle{ i }[/math]

[math]\displaystyle{ x^{(k)}_{i} = f^{(k)}_{i}/b^{(k)}_{i} - (c^{(k)}_{i}/b^{(k)}_{i}) x^{(k)}_{i-1} - (a^{(k)}_{i}/b^{(k)}_{i}) x^{(k)}_{i+1} }[/math]

c "левого края" системы будет

[math]\displaystyle{ x^{(k)}_{1} = f^{(k)}_{1}/b^{(k)}_{1} - (a^{(k)}_{1}/b^{(k)}_{1}) x^{(k)}_{2} }[/math]

а с "правого"

[math]\displaystyle{ x^{(k)}_{n} = f^{(k)}_{n}/b^{(k)}_{n} - (c^{(k)}_{n}/b^{(k)}_{n}) x^{(k)}_{n-1} }[/math]

После вычисления всех [math]\displaystyle{ x^{(0)}_{i} }[/math] значения искомых неизвестных [math]\displaystyle{ x_{i} = x^{(k)}_{i} }[/math].

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

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

[math]\displaystyle{ c^{(k)}_{i}/b^{(k)}_{i} , a^{(k)}_{i}/b^{(k)}_{i} , f^{(k)}_{i}/b^{(k)}_{i} }[/math]

почти все используются для преобразований двух чётных уравнений, то при выделении "микровычислений", из которых следует составить шаги редукции и которые составляют его ядро, лучше отнести вычисления этих отношений к предыдущему шагу редукции. Таким образом, на "подготовительном шаге" микроядро будет для каждого нечётного [math]\displaystyle{ i }[/math] состоять только из трёх делений (кроме [math]\displaystyle{ 2 }[/math] и [math]\displaystyle{ n }[/math] - там будет по 2 деления), а затем на каждом последующем шаге редукции для каждого нечётного [math]\displaystyle{ i }[/math] - из двух умножений и двух вычислений выражений типа [math]\displaystyle{ a-bc-de }[/math], с последующими тремя делениями (на "краях" часть этих операций отсутствует или урезана, но общую картину это не очень меняет). Для чётных [math]\displaystyle{ i }[/math] деления в микроядре будут отсутствовать.

Что касается шагов обратного хода, то там для каждого [math]\displaystyle{ i }[/math] рано или поздно выполняется одна операция типа [math]\displaystyle{ a-bc-de }[/math] (на "краях" - типа [math]\displaystyle{ a-bc }[/math]).

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

Рисунок 4. Граф алгоритма прямого хода циклической редукции при n=15. Светло-зелёным обозначены микроядра "подготовительного шага" с тремя (или меньше) делениями, синим - микроядра нечётных узов, сине-зелёным - ядра чётных узлов (без делений).

Рисунок 5. Граф алгоритма обратного хода циклической редукции при n=15. В вершинах - операции вида a-bc-de, в вершинах с 1 входящей дугой - вида a-bc. Показаны только дуги, передающие значения найденных неизвестных.

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

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