(The documentation of function matrix_inverse_givenLU)
Last Update: 2021/12/6
◆機能・用途(Purpose)
LU分解を用いて正方行列の逆行列を計算する。
行列そのものではなく、そのLU分解の結果を用いるバージョン。
Calculate the inverse of a square matrix using an LU decomposition.
The LU decomposition is given instead of the matrix itself.
◆形式(Format)
#include <matrix/inverse.h>
inline struct matrix matrix_inverse_givenLU(const struct LU A_LU)
◆引数(Arguments)
A_LU
逆行列を計算したい正方行列\(\myvector{A}\)のLU分解。
The LU decomposition of a square matrix \(\myvector{A}\)
for which the inverse is to be calculated.
◆戻り値(Return value)
行列\(\myvector{A}\)の逆行列。
The inverse of a matrix \(\myvector{A}\).
◆使用例(Example)
struct matrix A;
struct LU A_LU=LU_decomposition(A);
struct matrix Ainv=matrix_inverse_givenLU(A_LU);
◆計算方法(Computation method)
\(\myvector{A}\)を\(N\times N\)の正方行列、
\(\myvector{A}^{-1}\)をその逆行列、
\(\myvector{I}\)を\(N\times N\)の単位行列とするとき、
これらの間には
\[\begin{equation}
\myvector{A}\myvector{A}^{-1}=\myvector{I}
\label{eq.inverse.definition}
\end{equation}\]
の関係が成り立つ。
\(\myvector{I}\)は第\(j\)行のみが1で
それ以外の成分が全て0の列ベクトル\(\myvector{e_j}\)
を用いて
\[\begin{equation}
\myvector{I}=
\begin{pmatrix}
\myvector{e_1} & \cdots & \myvector{e_N}
\end{pmatrix}
\label{eq.I.decomposition}
\end{equation}\]
のように分解できる。
\(\myvector{A}^{-1}\)の第\(j\)列のみを取り出した列ベクトルを
\(\myvector{x_j}\)とおけば
\[\begin{equation}
\myvector{A}^{-1}=
\begin{pmatrix}
\myvector{x_1} & \cdots & \myvector{x_N}
\end{pmatrix}
\label{eq.Ainv.decomposition}
\end{equation}\]
である。(\ref{eq.I.decomposition})(\ref{eq.Ainv.decomposition})を
(\ref{eq.inverse.definition})に代入すると
\[\begin{equation}
\myvector{A}
\begin{pmatrix}
\myvector{x_1} & \cdots & \myvector{x_N}
\end{pmatrix}
=
\begin{pmatrix}
\myvector{e_1} & \cdots & \myvector{e_N}
\end{pmatrix}
\label{eq.problem.use_columnvectors}
\end{equation}\]
となるので\(N\)個の連立1次方程式
\[\begin{equation}
\myvector{A}\myvector{x_j}=\myvector{e_j}
\hspace{2em}(j=1,\cdots,N)
\label{eq.problem.component}
\end{equation}\]
に分解できる。これらを
関数solution_LU_givenLUdecomposition (matrix/inverse.h)
を用いて解くことにより
\(\myvector{x_1},\cdots,\myvector{x_N}\)を求めれば
(\ref{eq.Ainv.decomposition})式によって逆行列が求まる。
「Numerical Recipes in C」(Press et al., 1988, Cambridge University Press)
によれば逆行列の計算に当たって
このアプローチよりも効率の良い方法は無いとのことである。
Let \(\myvector{A}\) be an \(N\times N\) square matrix,
\(\myvector{A}^{-1}\) be the inverse matrix of it,
and \(\myvector{I}\) be an \(N\times N\) unit matrix.
They satisfy eq. (\ref{eq.inverse.definition}).
The matrix \(\myvector{I}\) can be decomposed into column vectors
as (\ref{eq.I.decomposition}),
where \(\myvector{e_j}\) is a column vector whose \(j\)th row is 1
and the other components are 0.
The matrix \(\myvector{A}^{-1}\) can also be decomposed
as (\ref{eq.Ainv.decomposition}),
where the column vector \(\myvector{x_j}\) represents
the \(j\)th column of \(\myvector{A}^{-1}\).
Substituting eqs (\ref{eq.I.decomposition}) and (\ref{eq.Ainv.decomposition})
into (\ref{eq.inverse.definition}) results in
(\ref{eq.problem.use_columnvectors}).
Thus the problem can be decomposed into
\(N\) sets of simultaneous equations given by (\ref{eq.problem.component}).
Solving these equations for \(\myvector{x_1},\cdots,\myvector{x_N}\)
using function solution_LU_givenLUdecomposition
(matrix/inverse.h),
\(\myvector{A}^{-1}\) can be constructed by eq. (\ref{eq.Ainv.decomposition}).
According to “Numerical Recipes in C”
(Press et al., 1988, Cambridge University Press),
there is no more efficient method than this approach
for the computation of an inverse matrix.