Алгоритм сдваивания Стоуна для LU-разложения трёхдиагональной матрицы

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

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

Алгоритм сдваивания Стоуна для LU-разложения трёхдиагональной матрицы - часть метода сдваивания Стоуна для решения СЛАУ вида [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]

Метод сдваивания Стоуна впервые предложен в начале 70-х гг. 20го века в качестве альтернативы другим параллельным алгоритмам решения трёхдиагональных СЛАУ, например, методу циклической редукции.

Здесь рассматривается его первая часть - [math]\displaystyle{ LU }[/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} }[/math]

в виде произведения матриц

[math]\displaystyle{ L = \begin{bmatrix} 1 & 0 & 0 & \cdots & \cdots & 0 \\ l_{21} & 1 & 0 & \cdots & \cdots & 0 \\ 0 & l_{32} & 1 & \cdots & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots & \ddots & 0 \\ 0 & \cdots & \cdots & l_{n-1 n-2} & 1 & 0 \\ 0 & \cdots & \cdots & 0 & l_{n n-1} & 1 \\ \end{bmatrix} }[/math]

и

[math]\displaystyle{ U = \begin{bmatrix} u_{11} & u_{12} & 0 & \cdots & \cdots & 0 \\ 0 & u_{22} & u_{23}& \cdots & \cdots & 0 \\ 0 & 0 & u_{33} & \cdots & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots & \ddots & 0 \\ 0 & \cdots & \cdots & 0 & u_{n-1 n-1} & u_{n-1 n} \\ 0 & \cdots & \cdots & 0 & 0 & u_{n n} \\ \end{bmatrix} }[/math]

Важным моментом является то, что в условиях точных вычислений алгоритм сдваивания Стоуна вычисляет то же самое разложение, что и компактная схема метода Гаусса.

Теоретически метод Стоуна основан на том, что, если ввести величины [math]\displaystyle{ q_i }[/math] так, что [math]\displaystyle{ q_0 = 1 }[/math], [math]\displaystyle{ q_1 = a_{11} }[/math] и [math]\displaystyle{ q_i = a_{ii} q_{i-1} - a_{ii-1}a_{i-1i} q_{i-2} }[/math], то после вычисления всех [math]\displaystyle{ q_i }[/math] легко вычислить все [math]\displaystyle{ u_{ii} = q_i/q_{i-1} }[/math], [math]\displaystyle{ l_{ii-1} = a_{ii-1}/u_{i-1i-1} }[/math]. При этом не нужно вычислять [math]\displaystyle{ u_{ii+1} = a_{ii+1} }[/math].

Вычисление всех [math]\displaystyle{ q_i }[/math] производится с использованием ассоциативности матричного умножения. Из рекуррентной связи [math]\displaystyle{ q_i = a_{ii} q_{i-1} - a_{ii-1}a_{i-1i} q_{i-2} }[/math] следует матричное равенство

[math]\displaystyle{ \begin{bmatrix} q_i \\ q_{i-1} \\ \end{bmatrix} = \begin{bmatrix} a_{ii} & -a_{ii-1}a_{i-1i} \\ 0 & 1 \\ \end{bmatrix} \begin{bmatrix} q_{i-1} \\ q_{i-2} \\ \end{bmatrix} = A_i \begin{bmatrix} q_{i-1} \\ q_{i-2} \\ \end{bmatrix} = A_i A_{i-1}...A_2 \begin{bmatrix} q_{1} \\ q_{0} \\ \end{bmatrix} = A_i A_{i-1}...A_2 \begin{bmatrix} a_{11} \\ 1 \\ \end{bmatrix} }[/math]

Произведения матриц [math]\displaystyle{ A_i A_{i-1}...A_2 }[/math] вычисляются параллельной схемой сдваивания.

Программная реализация алгоритма

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

Из-за большой избыточности вычислений алгоритм сдваивания Стоуна никогда не предназначался для последовательной реализации.

Существующие реализации алгоритма

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