Алгоритм Бартельса – Стюарта - Bartels–Stewart algorithm

В числовая линейная алгебра, то Алгоритм Бартелса – Стюарта используется для численного решения Матричное уравнение Сильвестра ${ Displaystyle AX-XB = C}$. Разработано R.H. Bartels и G.W. Стюарт в 1971 году,^[1] это был первый численно стабильный метод, который можно систематически применять для решения таких уравнений. Алгоритм работает с использованием вещественные разложения Шура из ${ displaystyle A}$ и ${ displaystyle B}$ преобразовывать ${ Displaystyle AX-XB = C}$ в треугольную систему, которую затем можно решить с помощью прямой или обратной замены. В 1979 г. Г. Голуб, К. Ван Лоан и С. Нэш представили улучшенную версию алгоритма,^[2] известный как алгоритм Хессенберга – Шура. Остается стандартным подходом к решению Уравнения Сильвестра когда ${ displaystyle X}$ от малого до среднего размера.

Алгоритм

Позволять ${ displaystyle X, C in mathbb {R} ^ {m times n}}$, и предположим, что собственные значения ${ displaystyle A}$ отличны от собственных значений ${ displaystyle B}$. Тогда матричное уравнение ${ Displaystyle AX-XB = C}$ имеет уникальное решение. Алгоритм Бартелса – Стюарта вычисляет ${ displaystyle X}$ выполнив следующие действия:^[2]

1. вычислить вещественные разложения Шура

${ displaystyle R = U ^ {T} AU,}$

${ Displaystyle S = V ^ {T} B ^ {T} V.}$

Матрицы ${ displaystyle R}$ и ${ displaystyle S}$ - блочно-верхние треугольные матрицы, с диагональными блоками размером ${ Displaystyle 1 раз 1}$ или ${ displaystyle 2 times 2}$.

2. Установить ${ displaystyle F = U ^ {T} CV.}$

3. Решите упрощенную систему ${ displaystyle RY-YS ^ {T} = F}$, где ${ Displaystyle Y = U ^ {T} XV}$. Это можно сделать, используя прямую подстановку блоков. В частности, если ${ displaystyle s_ {k-1, k} = 0}$, тогда

${ Displaystyle (R-s_ {kk} I) y_ {k} = f_ {k} + sum _ {j = k + 1} ^ {n} s_ {kj} y_ {j},}$

где ${ displaystyle y_ {k}}$это ${ displaystyle k}$-й столбец ${ displaystyle Y}$. Когда ${ displaystyle s_ {k-1, k} neq 0}$, столбцы ${ Displaystyle [y_ {k-1} mid y_ {k}]}$ должны быть объединены и решены одновременно.

4. Установите ${ displaystyle X = UYV ^ {T}.}$

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

С использованием QR-алгоритм, то вещественные разложения Шура на шаге 1 требуется примерно ${ displaystyle 10 (m ^ {3} + n ^ {3})}$ проваливается, так что общие вычислительные затраты составляют ${ displaystyle 10 (m ^ {3} + n ^ {3}) + 2,5 (mn ^ {2} + nm ^ {2})}$.^[2]

Упрощения и особые случаи

В частном случае, когда ${ displaystyle B = -A ^ {T}}$ и ${ displaystyle C}$ симметрично, решение ${ displaystyle X}$ также будет симметричным. Эту симметрию можно использовать так, чтобы ${ displaystyle Y}$ находится более эффективно на шаге 3 алгоритма.^[1]

Алгоритм Хессенберга – Шура

Алгоритм Хессенберга – Шура^[2] заменяет разложение ${ Displaystyle R = U ^ {T} AU}$ на шаге 1 с разложением ${ Displaystyle H = Q ^ {T} AQ}$, где ${ displaystyle H}$ является матрица верхнего Гессенберга. Это приводит к системе вида ${ displaystyle HY-YS ^ {T} = F}$ это можно решить с помощью прямой замены. Преимущество такого подхода в том, что ${ Displaystyle H = Q ^ {T} AQ}$ можно найти с помощью Размышления домохозяина по цене ${ Displaystyle (5/3) м ^ {3}}$ провалы, по сравнению с ${ displaystyle 10m ^ {3}}$ флопов, необходимых для вычисления реального разложения Шура ${ displaystyle A}$.

Программное обеспечение и реализация

Подпрограммы, необходимые для варианта Хессенберга-Шура алгоритма Бартельса-Стюарта, реализованы в библиотеке SLICOT. Они используются в наборе инструментов системы управления MATLAB.

Альтернативные подходы

Для больших систем ${ Displaystyle { mathcal {O}} (м ^ {3} + п ^ {3})}$ Стоимость алгоритма Бартелса – Стюарта может быть непомерно высокой. Когда ${ displaystyle A}$ и ${ displaystyle B}$ являются разреженными или структурированными, поэтому линейные решения и матричные векторные умножения с их участием эффективны, итерационные алгоритмы потенциально могут работать лучше. К ним относятся методы на основе прогнозов, которые используют Крыловское подпространство итерации, методы на основе переменное направление неявное (ADI) итерации и гибридизации, включающие как проекцию, так и ADI.^[3] Итерационные методы также могут использоваться для непосредственного построения приближения низкого ранга к ${ displaystyle X}$ при решении ${ Displaystyle AX-XB = C}$.

использованная литература