几种求逆矩阵的方法总结，以Matlab语言为例
- *0* 引言
- *1* 简单描述+函数实现
- *2* 方法调用计算对比
0 引言
最近在使用函数库求解逆矩阵的时候发现同一个矩阵使用不同的语言、不同的求解方法会产生不同精度的结果，特别是阶数很高的方阵，一些库中的算法为了使计算速度提升，用了非常规方法。这里整理了3种常用的矩阵求逆方法，伴随矩阵法、LU分解法和高斯消元法，并用Matlab进行了实现，一些过程参考了NET博文：
1 简单描述+函数实现
伴随矩阵法
对于一个矩阵A，如果它的伴随矩阵存在，并且A的行列式不为零，那么A的逆矩阵等于它的伴随矩阵除以A的行列式。
设A是一个n阶矩阵，其伴随矩阵为A*。那么A的逆矩阵为 A − 1 A^{-1} A−1，根据公式：
A − 1 = A ∗ / ∣ A ∣ A^{-1} = A* / |A| A−1=A∗/∣A∣
其中，|A|表示A的行列式。
因此，如果你已经求得了矩阵A的伴随矩阵A*，并且知道A的行列式|A|不为零，那么你可以通过上述公式来求得A的逆矩阵A^(-1)。
% 伴随矩阵法求逆矩阵
function value = inv_Adjoint(matrixA,tol)
As = size(matrixA);
N = As(1);
if(nargin < 2)
tol = 1e-15;end
if( abs(det(matrixA)) < tol )
msg ='the matrix is not full rank';
error(msg);
end
bb = adjoint(matrixA); % 伴随矩阵法
d = det(matrixA); % 矩阵的行列式子
value = bb/d;
end
LU分解法；
LU分解是一种将矩阵分解为一个下三角矩阵和一个上三角矩阵的方法。这种分解方法可以有效地将矩阵求逆的计算量减少一半。具体步骤如下：
- 将原矩阵表示为A = LU，其中L是一个下三角矩阵，U是一个上三角矩阵。
- 解下三角矩阵方程LY = I，其中Y是一个列向量，I是单位矩阵。
- 解上三角矩阵方程UX = Y，其中X是一个列向量。
- 求解出的X就是原矩阵A的逆矩阵。
这种方法的好处是可以避免直接使用矩阵求逆的方法，提高计算精度。由于LU分解方法是通过矩阵的分解来求逆，因此可以在一定程度上避免了矩阵求逆时可能出现的数值不稳定性。
然而，需要注意的是，当原矩阵的行列式为0时，即矩阵不可逆时，LU分解方法无法计算逆矩阵。在这种情况下，我们需要采用其他方法来求解。
总结起来，对于矩阵求逆，LU分解是一种精度较高的方法。它通过将矩阵分解为下三角矩阵和上三角矩阵来进行计算，从而提高了计算精度，并且避免了一些数值不稳定性的问题。
% LU分解求逆
function value = inv_lu(matrixA,tol)
As = size(matrixA);
N = As(1);
value = zeros(N,N);
if(nargin < 2)
tol = 1e-15;end
if( abs(det(matrixA)) < tol )
msg ='the matrix is not full rank';
error(msg);
end
MatrixB=eye(N);
Y = value;
array1 =1:N;
for i=1:N-1
[~, j]=max(abs(matrixA(i:N,i)));
num1 = matrixA(i,:);
matrixA(i,:) = matrixA(j+i-1,:);
matrixA(j+i-1,:) = num1;
num1 = array1(i);
array1(i) = array1(j+i-1);
array1(j+i-1) = num1;
if (matrixA(i,i)==0)
break
end
for j=i+1:N
num1=matrixA(j,i)/matrixA(i,i);
matrixA(j,i) = num1;
matrixA(j,i+1:N)=matrixA(j,i+1:N)-num1*matrixA(i,i+1:N);
end
end
for i=1:N
Y(1,i) = MatrixB(array1(1),i);
for j=2:N
Y(j,i)= MatrixB(array1(j),i)-matrixA(j,1:j-1)*Y(1:j-1,i);
end
value(N,i)=Y(N,i)/matrixA(N,N);
for j=N-1:-1:1
value(j,i)=(Y(j,i)-matrixA(j,j+1:N)*value(j+1:N,i))/matrixA(j,j);
end
end
end
高斯消元法;
要用高斯消元法求一个n*n矩阵的逆矩阵，可以按以下步骤进行：
- 将原矩阵和单位矩阵合并成一个增广矩阵，形成一个4*8的矩阵。
- 利用高斯消元法将矩阵的左半部分化为上三角矩阵。
- 利用高斯消元法将矩阵的右半部分也化为上三角矩阵。
- 将得到的上三角矩阵进行回代，将其化为单位矩阵。
- 得到的矩阵的右半部分即为所求的逆矩阵。
% 高斯消元法求逆矩阵
function value = inv_gaosi(matrixA,tol)
rc = size(matrixA);
N = rc(1);
value = zeros(N,N);
if(nargin < 2)
tol = 1e-15;end
if( abs(det(matrixA)) < tol )
msg ='the matrix is not full rank';
error(msg);
end
Meye = eye(N);
for i=1:N
maxV = matrixA(i,i);
index=i;
for j=i+1:N
if(abs(matrixA(j,i))>abs(maxV))
maxV = matrixA(j,i);
index=j;
end
end
for j=1:N
num1 = matrixA(i,j);
matrixA(i,j) = matrixA(index,j);
matrixA(index,j) = num1;
num1=Meye(i,j);
Meye(i,j) = Meye(index,j);
Meye(index,j) = num1;
end
for j=i+1:N
num1 = matrixA(j,i)/matrixA(i,i);
for k=1:N
matrixA(j,k) = matrixA(j,k) - num1*matrixA(i,k);
Meye(j,k) = Meye(j,k)-num1*Meye(i,k);
end
end
end
for i=N:-1:1
num1 =matrixA(i,i);
for j=i:N
matrixA(j,i)=matrixA(j,i)/num1;
end
for j=1:N
Meye(i,j)=Meye(i,j)/num1;
end
for j=i-1:-1:1
num1=matrixA(j,i);
matrixA(j,i)=0;
for k=1:N
Meye(j,k)=Meye(j,k)-num1*Meye(i,k);
end
end
end
value=Meye;
end
2 方法调用计算对比
clc;clear;
format long
a = [0.863926332645926,2.744822685285432E-002,0.358800923218055,7.587357258855056E-002;...
0.737486519704736,0.275982186877394,0.928772666000609,0.219795611073554;...
0.172452279207503,0.378156212690881,0.506403133293737,1.243293753676099E-002;...
0.938238414353813,0.119337792575225,0.952934573404228,0.859760719854225];
aa = inv(a); %
bb =inv_gaosi(a); % 高斯消元法
cc =inv_lu(a);
dd = inv_Adjoint(a);
%% .............................................................
%% ......... out ............
aa =
1.997894526750277 -1.386596126529541 0.814427151532846 0.166388882786055
1.987598379664842 -4.975548628572671 5.811811940191023 1.012537147110808
-2.163170299120444 4.249207887174017 -2.673735235716419 -0.856735351177312
-0.058547952499082 -2.505916255800319 1.268026177518930 1.790575344111539
>> bb
bb =
1.997894526750277 -1.386596126529541 0.814427151532846 0.166388882786054
1.987598379664842 -4.975548628572672 5.811811940191024 1.012537147110807
-2.163170299120444 4.249207887174018 -2.673735235716419 -0.856735351177311
-0.058547952499082 -2.505916255800319 1.268026177518930 1.790575344111538
>> cc
cc =
1.997894526750277 -1.386596126529541 0.814427151532846 0.166388882786054
1.987598379664842 -4.975548628572672 5.811811940191024 1.012537147110807
-2.163170299120444 4.249207887174018 -2.673735235716419 -0.856735351177311
-0.058547952499082 -2.505916255800319 1.268026177518930 1.790575344111538
>> dd
dd =
1.997894526750279 -1.386596126529543 0.814427151532847 0.166388882786054
1.987598379664844 -4.975548628572673 5.811811940191022 1.012537147110807
-2.163170299120446 4.249207887174019 -2.673735235716419 -0.856735351177311
-0.058547952499082 -2.505916255800320 1.268026177518929 1.790575344111539
仔细观察不难看出，在小数点后15位左右不同方法间会存在微笑差异，这个差异或者叫误差在一般计算中无所谓，可以被忽略，但对精度要求特高的过程，就需要斟酌一下，使用哪种求逆方法可以最大限度减小计算误差带来的影响。
还有哪些好的方法可以留言交流，希望对你有所帮助！
