Обратная матрица

Алгоритм

В пакете ALGLIB реализован рекурсивный алгоритм обращения матриц. Исходная матрица A, заданная одной из своих треугольных факторизаций (LU-разложением или разложением Холецкого) разбивается на четыре части, после чего для матрицы, заданной LU-разложением, используется следующая формула:

Для вычисления матричных произведений общего вида и для умножения на U^-1/L^-1 используются кэш-эффективные подпрограммы из ALGLIB BLAS. Выражения вида (LU)^-1 вычисляются рекурсивно. Базовые случае (размер матрицы достаточно мал, чтобы она поместилась в кэш процессора) решаются с использованием другого алгоритма, основанного на явном обращении треугольных матриц, входящих в треугольную факторизацию A. Использование рекурсии имеет следующие преимущества:

Становится возможным эффективно использовать кэш процессора без явного указания его размера. С каждой итерацией обрабатываемая матрица становится меньше, пока не поместится полностью в кэш. Такие алгоритмы называются cache oblivious.
Рекурсивный алгоритм тесно использует ALGLIB BLAS, главным образом - операции третьего уровня (произведения матрица-матрица), обладающие максимальной эффективностью. При этом перемножаемые матрицы обычно квадратные (в противоположность перемножению узких длинных матриц), что позволяет заметно повысить быстродействие.
Реализованный таким образом алгоритм обладает следующим свойством - его быстродействие, выраженное в числе операций в секунду, растет с ростом размера матрицы. Чем крупнее матрица, тем эффективнее используется процессор.

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

Подпрограммы

Ниже приведен список подпрограмм для обращения матриц, доступных в ALGLIB:

rmatrixinverse
spdmatrixinverse
cmatrixinverse
hpdmatrixinverse
rmatrixluinverse
spdmatrixcholeskyinverse
cmatrixluinverse
hpdmatrixcholeskyinverse

Первые четыре подпрограммы служат для обращения вещественных (префикс R), симметричных положительно определенных (префикс SPD), комплексных (префикс C) и Эрмитовых положительно определенных (префикс HPD) матриц, заданных в явной форме. Следующие четыре подпрограммы решают ту же задачу, но для матриц, заданных соответствующей треугольной факторизацией.

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