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

• ${\displaystyle \det A^{-1}={\frac {1}{\det A}}}$ , где ${\displaystyle \ \det }$  обозначает определитель.
• ${\displaystyle \ (AB)^{-1}=B^{-1}A^{-1}}$  для двух квадратных обратимых матриц ${\displaystyle A}$  и ${\displaystyle B}$ .
• ${\displaystyle \ (A^{T})^{-1}=(A^{-1})^{T}}$ , где ${\displaystyle (...)^{T}}$  обозначает транспонированную матрицу.
• ${\displaystyle \ (kA)^{-1}=k^{-1}A^{-1}}$  для любого коэффициента ${\displaystyle k\not =0}$ .
• ${\displaystyle \ E^{-1}=E}$ .
• Если необходимо решить систему линейных уравнений ${\displaystyle Ax=b}$ , (b — ненулевой вектор) где ${\displaystyle x}$  — искомый вектор, и если ${\displaystyle A^{-1}}$  существует, то ${\displaystyle x=A^{-1}b}$ . В противном случае либо размерность пространства решений больше нуля, либо их нет вовсе.

Если матрица обратима, то для нахождения обратной матрицы можно воспользоваться одним из следующих способов:

Точные (прямые) методы

Метод Жордана—Гаусса

Возьмём две матрицы: саму A и единичную E. Приведём матрицу A к единичной матрице методом Гаусса—Жордана применяя преобразования по строкам (можно также применять преобразования и по столбцам). После применения каждой операции к первой матрице применим ту же операцию ко второй. Когда приведение первой матрицы к единичному виду будет завершено, вторая матрица окажется равной A⁻¹.

При использовании метода Гаусса первая матрица будет умножаться слева на одну из элементарных матриц ${\displaystyle \Lambda _{i}}$  (трансвекцию или диагональную матрицу с единицами на главной диагонали, кроме одной позиции):

${\displaystyle \Lambda _{1}\cdot \dots \cdot \Lambda _{n}\cdot A=\Lambda A=E\Rightarrow \Lambda =A^{-1}}$ .

${\displaystyle \Lambda _{m}={\begin{bmatrix}1&\dots &0&-a_{1m}/a_{mm}&0&\dots &0\\&&&\dots &&&\\0&\dots &1&-a_{m-1m}/a_{mm}&0&\dots &0\\0&\dots &0&1/a_{mm}&0&\dots &0\\0&\dots &0&-a_{m+1m}/a_{mm}&1&\dots &0\\&&&\dots &&&\\0&\dots &0&-a_{nm}/a_{mm}&0&\dots &1\end{bmatrix}}}$ .

Вторая матрица после применения всех операций станет равна ${\displaystyle \Lambda }$ , то есть будет искомой. Сложность алгоритма — ${\displaystyle O(n^{3})}$ .

С помощью матрицы алгебраических дополнений

Матрица, обратная матрице ${\displaystyle A}$ , представима в виде

${\displaystyle {A}^{-1}={{{\mbox{adj}}(A)} \over {\det(A)}}}$ 

где ${\displaystyle {\mbox{adj}}(A)}$  — присоединенная матрица (матрица, составленная из алгебраических дополнений для соответствующих элементов транспонированной матрицы).

Сложность алгоритма зависит от сложности алгоритма расчета определителя O_det и равна O(n²)·O_det.

Использование LU/LUP-разложения

Матричное уравнение ${\displaystyle AX=I_{n}}$  для обратной матрицы ${\displaystyle X}$  можно рассматривать как совокупность ${\displaystyle n}$  систем вида ${\displaystyle Ax=b}$ . Обозначим ${\displaystyle i}$ -й столбец матрицы ${\displaystyle X}$  через ${\displaystyle X_{i}}$ ; тогда ${\displaystyle AX_{i}=e_{i}}$ , ${\displaystyle i=1,\ldots ,n}$ , поскольку ${\displaystyle i}$ -м столбцом матрицы ${\displaystyle I_{n}}$  является единичный вектор ${\displaystyle e_{i}}$ . другими словами, нахождение обратной матрицы сводится к решению n уравнений с одной матрицей и разными правыми частями. После выполнения LUP-разложения (время O(n³)) на решение каждого из n уравнений нужно время O(n²), так что и эта часть работы требует времени O(n³)^[1].

Если матрица A невырождена, то для неё можно рассчитать LUP-разложение ${\displaystyle PA=LU}$ . Пусть ${\displaystyle PA=B}$ , ${\displaystyle B^{-1}=D}$ . Тогда из свойств обратной матрицы можно записать: ${\displaystyle D=U^{-1}L^{-1}}$ . Если умножить это равенство на U и L то можно получить два равенства вида ${\displaystyle UD=L^{-1}}$  и ${\displaystyle DL=U^{-1}}$ . Первое из этих равенств представляет собой систему из n² линейных уравнений для ${\displaystyle {\frac {n(n+1)}{2}}}$ , из которых известны правые части (из свойств треугольных матриц). Второе представляет также систему из n² линейных уравнений для ${\displaystyle {\frac {n(n-1)}{2}}}$ , из которых известны правые части (также из свойств треугольных матриц). Вместе они представляют собой систему из n² равенств. С помощью этих равенств можно рекуррентно определить все n² элементов матрицы D. Тогда из равенства (PA)⁻¹ = A⁻¹P⁻¹ = B⁻¹ = D получаем равенство ${\displaystyle A^{-1}=DP}$ .

В случае использования LU-разложения не требуется перестановки столбцов матрицы D, но решение может разойтись даже если матрица A невырождена.

Сложность алгоритма — O(n³).

Итерационные методы

Методы Шульца

${\displaystyle {\begin{cases}\Psi _{k}=E-AU_{k},\\U_{k+1}=U_{k}\sum _{i=0}^{n}\Psi _{k}^{i}\end{cases}}}$ 

Оценка погрешности

Выбор начального приближения

Проблема выбора начального приближения ${\displaystyle U_{0}}$  в рассматриваемых здесь процессах итерационного обращения матриц не позволяет относиться к ним как к самостоятельным универсальным методам, конкурирующими с прямыми методами обращения, основанными, например, на LU-разложении матриц. Имеются некоторые рекомендации по выбору ${\displaystyle U_{0}}$ , обеспечивающие выполнение условия ${\displaystyle \rho (\Psi _{0})<1}$  (спектральный радиус матрицы меньше единицы), являющегося необходимым и достаточным для сходимости процесса. Однако при этом, во-первых, требуется знать сверху оценку спектра обращаемой матрицы A либо матрицы ${\displaystyle AA^{T}}$  (а именно, если A — симметричная положительно определённая матрица и ${\displaystyle \rho (A)\leq \beta }$ , то можно взять ${\displaystyle U_{0}={\alpha }E}$ , где ${\displaystyle \alpha \in \left(0,{\frac {2}{\beta }}\right)}$ ; если же A — произвольная невырожденная матрица и ${\displaystyle \rho (AA^{T})\leq \beta }$ , то полагают ${\displaystyle U_{0}={\alpha }A^{T}}$ , где также ${\displaystyle \alpha \in \left(0,{\frac {2}{\beta }}\right)}$ ; можно конечно упростить ситуацию и, воспользовавшись тем, что ${\displaystyle \rho (AA^{T})\leq {\mathcal {k}}AA^{T}{\mathcal {k}}}$ , положить ${\displaystyle U_{0}={\frac {A^{T}}{\|AA^{T}\|}}}$ ). Во-вторых, при таком задании начальной матрицы нет гарантии, что ${\displaystyle \|\Psi _{0}\|}$  будет малой (возможно, даже окажется ${\displaystyle \|\Psi _{0}\|>1}$ ), и высокий порядок скорости сходимости обнаружится далеко не сразу.

Матрица 2 × 2

${\displaystyle \mathbf {A} ^{-1}={\begin{bmatrix}a&b\\c&d\\\end{bmatrix}}^{-1}={\frac {1}{\det \mathbf {A} }}{\begin{bmatrix}d&-b\\-c&a\\\end{bmatrix}}={\frac {1}{ad-bc}}{\begin{bmatrix}d&-b\\-c&a\\\end{bmatrix}}}$ ^[2]

Обращение матрицы 2 × 2 возможно только при условии, что ${\displaystyle ad-bc=\det A\neq 0}$ .

• Реализация с полным выбором ведущего элемента на C++