用Matlab对矩阵进行LU分解法
简介
- 01
利用矩阵分解来求先行方程组，可以节省内存，节省计算时间，因此在工程计算中最常用的技术。其中LU分解法是最基本也是最常用的方法。
操作方法
- 01
现将系数矩阵A进行LU分解，得到LU=PA；
- 02
然后解Ly=Pb;
- 03
再解Ux=y得到原方程组的解。因为矩阵L，U的特殊结构，是的上面两个方程组可以很容易的求出来。下面给出LU分解法求解线性方程组Ax=b的Matlab程序。
- 04
function x=solvebyLU(A,b)
- 05
% 该函数利用LU分解法求线性方程组Ax=b的解 flag=isexist(A,b); %调用第一小节中的isexist函数判断方程组解的情况 if flag==0 disp('该方程组无解！'); x=[]; return; else r=rank(A); [m,n]=size(A); [L,U,P]=lu(A); b=P*b; % 解Ly=b y(1)=b(1); if m>1 for i=2:m y(i)=b(i)-L(i,1:i-1)*y(1:i-1)'; end end y=y'; % 解Ux=y得原方程组的一个特解 x0(r)=y(r)/U(r,r); if r>1 for i=r-1:-1:1 x0(i)=(y(i)-U(i,i+1:r)*x0(i+1:r)')/U(i,i); end end x0=x0'; if flag==1 %若方程组有唯一解 x=x0; return; else %若方程组有无穷多解 format rat; Z=null(A,'r'); %求出对应齐次方程组的基础解系 [mZ,nZ]=size(Z); x0(r+1:n)=0; for i=1:nZ t=sym(char([107 48+i])); k(i)=t; %取k=[k1,k2...,]; end x=x0; for i=1:nZ x=x+k(i)*Z(:,i); %将方程组的通解表示为特解加对应齐次通解形式 end end end
- 06
应用举例：
- 07
Matlab命令窗口输入下面程序 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% >> A=[1 1 -3 -1;3 -1 -3 4;1 5 -9 -8]; %方程系数矩阵 >> b=[1 4 0]'; >> x=solvebyLU(A,b) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 输出结果为： x = (3*k1)/2 - (3*k2)/4 + 5/4 (3*k1)/2 + (7*k2)/4 - 1/4 k1 k2