3.1 方程组的逆矩阵解法及其MATLAB 程序3.1.3 线性方程组有解的判定条件及其MATLAB 程序 判定线性方程组A n m ⨯b X =是否有解的MATLAB 程序function [RA,RB,n]=jiepb(A,b)B=[A b];n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0,disp('请注意:因为RA~=RB ,所以此方程组无解.') return endif RA==RB if RA==ndisp('请注意:因为RA=RB=n ,所以此方程组有唯一解.') elsedisp('请注意:因为RA=RB<n ,所以此方程组有无穷多解.') end end例3.1.4 判断下列线性方程组解的情况.如果有唯一解,则用表 3-2方法求解.(1) ⎪⎪⎩⎪⎪⎨⎧=-+-=+-+=-++=+-+;0742,0634,0723,05324321432143214321x x x x x x x x x x x x x x x x (2) ⎪⎪⎩⎪⎪⎨⎧=++-=+-+=-+-=+-+;0327,01613114,02332,075434321432143214321x x x x x x x x x x x x x x x x (3) ⎪⎩⎪⎨⎧=+=+-=-+;8311,1023,22421321321x x x x x x x x (4) ⎪⎩⎪⎨⎧=--+=+-+=+-+.12,2224,12w z y x w z y x w z y x解 在MATLAB 工作窗口输入程序>> A=[2 3 -1 5;3 1 2 -7;4 1 -3 6;1 -2 4 -7]; b=[ 0; 0; 0; 0]; [RA,RB,n]=jiepb(A,b)运行后输出结果为请注意:因为RA=RB=n ,所以此方程组有唯一解. RA = 4,RB =4,n =4 在MATLAB 工作窗口输入>>X=A\b,运行后输出结果为 X =(0 0 0 0)’.(2) 在MATLAB 工作窗口输入程序>> A=[3 4 -5 7;2 -3 3 -2;4 11 -13 16;7 -2 1 3];b=[ 0; 0; 0; 0];[RA,RB,n]=jiepb(A,b)运行后输出结果请注意:因为RA=RB<n ,所以此方程组有无穷多解. RA =2,RB =2,n =4(3) 在MATLAB 工作窗口输入程序>> A=[4 2 -1;3 -1 2;11 3 0]; b=[2;10;8]; [RA,RB,n]=jiepb(A,B)运行后输出结果请注意:因为RA~=RB ,所以此方程组无解. RA =2,RB =3,n =3(4)在MATLAB 工作窗口输入程序>> A=[2 1 -1 1;4 2 -2 1;2 1 -1 -1]; b=[1; 2; 1]; [RA,RB,n]=jiepb(A,b)运行后输出结果请注意:因为RA=RB<n ,所以此方程组有无穷多解. RA =2,RB =2,n =33.2 三角形方程组的解法及其MATLAB 程序3.2.2 解三角形方程组的MATLAB 程序 解上三角形线性方程组b AX =的MATLAB 程序function [RA,RB,n,X]=shangsan(A,b)B=[A b]; n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0,disp('请注意:因为RA~=RB ,所以此方程组无解.') return endif RA==RB if RA==ndisp('请注意:因为RA=RB=n ,所以此方程组有唯一解.') X=zeros(n,1); X(n)=b(n)/A(n,n); for k=n-1:-1:1X(k)=(b(k)-sum(A(k,k+1:n)*X(k+1:n)))/A(k,k);end elsedisp('请注意:因为RA=RB<n ,所以此方程组有无穷多解.')end end例3.2.2 用解上三角形线性方程组的MATLAB 程序解方程组⎪⎪⎩⎪⎪⎨⎧==+-=-+-=++-.63,456,7472,203254434324321x x x x x x x x x x . 解 在MATLAB 工作窗口输入程序>>A=[5 -1 2 3;0 -2 7 -4;0 0 6 5;0 0 0 3]; b=[20; -7; 4;6];[RA,RB,n,X]=shangsan(A,b)运行后输出结果请注意:因为RA=RB=n ,所以此方程组有唯一解. RA = RB =4, 4, n =4,X =[2.4 -4.0 -1.0 2.0]’3.3 高斯(Gauss )消元法和列主元消元法及其MATLAB 程序3.3.1 高斯消元法及其MATLAB 程序用高斯消元法解线性方程组b AX =的MATLAB 程序f unction [RA,RB,n,X]=gaus(A,b)B=[A b]; n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0,disp('请注意:因为RA~=RB ,所以此方程组无解.') return endif RA==RB if RA==ndisp('请注意:因为RA=RB=n ,所以此方程组有唯一解.') X=zeros(n,1); C=zeros(1,n+1); for p= 1:n-1for k=p+1:nm= B(k,p)/ B(p,p);B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1); end endb=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n); for q=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q); endelsedisp('请注意:因为RA=RB<n ,所以此方程组有无穷多解.') end end例3.3.2 用高斯消元法和MATLAB 程序求解下面的非齐次线性方程组,并且用逆矩阵解方程组的方法验证.⎪⎪⎩⎪⎪⎨⎧-=+---=+--=+--=-+-.142,16422,0,13432143214324321x x x x x x x x x x x x x x x 解 在MATLAB 工作窗口输入程序>> A=[1 -1 1 -3; 0 -1 -1 1;2 -2 -4 6;1 -2 -4 1]; b=[1;0; -1;-1]; [RA,RB,n,X] =gaus (A,b)运行后输出结果请注意:因为RA=RB=n ,所以此方程组有唯一解. RA =4RB =4n =43.3.2 列主元消元法及其MATLAB 程序用列主元消元法解线性方程组b AX =的MATLAB 程序function [RA,RB,n,X]=liezhu(A,b)B=[A b]; n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0,disp('请注意:因为RA~=RB ,所以此方程组无解.') return endif RA==RB if RA==ndisp('请注意:因为RA=RB=n ,所以此方程组有唯一解.') X=zeros(n,1); C=zeros(1,n+1); for p= 1:n-1X = 0 -0.5000 0.5000 0[Y,j]=max(abs(B(p:n,p))); C=B(p,:); B(p,:)= B(j+p-1,:); B(j+p-1,:)=C; for k=p+1:nm= B(k,p)/ B(p,p);B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1); end endb=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n); for q=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q); endelsedisp('请注意:因为RA=RB<n ,所以此方程组有无穷多解.') end end例3.3.3 用列主元消元法解线性方程组的MATLAB 程序解方程组⎪⎪⎩⎪⎪⎨⎧-=+---=+--=-+-=+--.142,16422,13,0432143214321432x x x x x x x x x x x x x x x . 解 在MATLAB 工作窗口输入程序>> A=[0 -1 -1 1;1 -1 1 -3;2 -2 -4 6;1 -2 -4 1]; b=[0;1;-1;-1]; [RA,RB,n,X]=liezhu(A,b)运行后输出结果请注意:因为RA=RB=n ,所以此方程组有唯一解. RA = 4,RB = 4,n = 4,X =[0 -0.5 0.5 0]’3.4 LU 分解法及其MATLAB 程序3.4.1判断矩阵LU 分解的充要条件及其MATLAB 程序 判断矩阵A 能否进行LU 分解的MATLAB 程序function hl=pdLUfj(A)[n n] =size(A); RA=rank(A); if RA~=ndisp('请注意:因为A 的n 阶行列式hl 等于零,所以A 不能进行LU 分解.A 的秩RA 如下:'), RA,hl=det(A); returnendif RA==nfor p=1:n,h(p)=det(A(1:p, 1:p));, endhl=h(1:n); for i=1:nif h(1,i)==0disp('请注意:因为A 的r 阶主子式等于零,所以A 不能进行LU 分解.A的秩RA 和各阶顺序主子式值hl 依次如下:'),hl;RA,returnend endif h(1,i)~=0disp('请注意:因为A 的各阶主子式都不等于零,所以A 能进行LU 分解.A的秩RA 和各阶顺序主子式值hl 依次如下:')hl;RA end end例3.4.1 判断下列矩阵能否进行LU 分解,并求矩阵的秩.(1)⎪⎪⎪⎭⎫ ⎝⎛6547121321;(2)⎪⎪⎪⎭⎫ ⎝⎛654721321;(3)⎪⎪⎪⎭⎫⎝⎛654321321. 解 (1)在MATLAB 工作窗口输入程序>> A=[1 2 3;1 12 7;4 5 6];hl=pdLUfj(A)运行后输出结果为请注意:因为A 的各阶主子式都不等于零,所以A 能进行LU 分解.A 的秩RA 和各阶顺序主子式值hl 依次如下:RA = 3, hl = 1 10 -48(2)在MATLAB 工作窗口输入程序>> A=[1 2 3;1 2 7;4 5 6];hl=pdLUfj(A)运行后输出结果为请注意:因为A 的r 阶主子式等于零,所以A 不能进行LU 分解.A 的秩RA 和各阶顺序主子式值hl 依次如下:RA = 3, hl =1 0 12(3)在MATLAB 工作窗口输入程序>> A=[1 2 3;1 2 3;4 5 6];hl=pdLUfj(A)运行后输出结果为请注意:因为A 的n 阶行列式hl 等于零,所以A 不能进行LU 分解.A 的秩RA 如下RA = 2, hl = 03.4.2 直接LU 分解法及其MATLAB 程序 将矩阵A 进行直接LU 分解的MATLAB 程序function hl=zhjLU(A)[n n] =size(A); RA=rank(A); if RA~=ndisp('请注意:因为A 的n 阶行列式hl 等于零,所以A 不能进行LU 分解.A的秩RA 如下:'), RA,hl=det(A);return endif RA==n for p=1:nh(p)=det(A(1:p, 1:p)); endhl=h(1:n); for i=1:nif h(1,i)==0disp('请注意:因为A 的r 阶主子式等于零,所以A 不能进行LU 分解.A的秩RA 和各阶顺序主子式值hl 依次如下:'), hl;RAreturn end endif h(1,i)~=0disp('请注意:因为A 的各阶主子式都不等于零,所以A 能进行LU 分解.A 的秩RA 和各阶顺序主子式值hl 依次如下:')for j=1:nU(1,j)=A(1,j); endfor k=2:n for i=2:n for j=2:nL(1,1)=1;L(i,i)=1; if i>jL(1,1)=1;L(2,1)=A(2,1)/U(1,1);L(i,1)=A(i,1)/U(1,1);L(i,k)=(A(i,k)- L(i,1:k-1)*U(1:k-1,k))/U(k,k); elseU(k,j)=A(k,j)-L(k,1:k-1)*U(1:k-1,j); end end end endhl;RA,U,L end end例3.4.3 用矩阵进行直接LU 分解的MA TLAB 程序分解矩阵⎪⎪⎪⎪⎪⎭⎫⎝⎛=3010342110100201A . 解 在MATLAB 工作窗口输入程序>> A=[1 0 2 0;0 1 0 1;1 2 4 3;0 1 0 3]; hl=zhjLU(A)运行后输出结果请注意:因为A 的各阶主子式都不等于零,所以A 能进行LU 分解.A 的秩RA和各阶顺序主子式值hl 依次如下:RA = 4U = 1 0 2 00 1 0 10 0 2 10 0 0 2 3.4.4 判断正定对称矩阵的方法及其MATLAB 程序 判断矩阵A 是否是正定对称矩阵的MATLAB 程序function hl=zddc(A) [n n] =size(A); for p=1:nh(p)=det(A(1:p, 1:p)); endhl=h(1:n);zA=A'; for i=1:nif h(1,i)<=0disp('请注意:因为A 的各阶顺序主子式hl 不全大于零,所以A 不是正定的.A 的转置矩阵zA 和各阶顺序主子式值hl 依次如下:'), hl;zA,returnend endif h(1,i)>0disp('请注意:因为A 的各阶顺序主子式hl 都大于零,所以A 是正定的.A 的转置矩阵zA 和各阶顺序主子式值hl 依次如下:')hl;zA end例3.4.5 判断下列矩阵是否是正定对称矩阵:L = 1 0 0 0 0 1 0 0 1 2 1 0 0 1 0 1 hl = 1 1 2 4(1)⎪⎪⎪⎪⎪⎭⎫ ⎝⎛--98754113211143214321.0;(2) ⎪⎪⎪⎪⎪⎭⎫ ⎝⎛------19631690230311211; (3) ⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛----212100212100002121002121;(4)⎪⎪⎪⎭⎫⎝⎛---401061112. 解 (1)在MATLAB 工作窗口输入程序>> A=[0.1 2 3 4;-1 2 -3 4;11 21 13 41;5 7 8 9];hl=zddc (A)运行后输出结果请注意: A 不是对称矩阵请注意:因为A 的各阶顺序主子式hl 不全大于零,所以A 不是正定的.A 的转置矩阵zA 和各阶顺序主子式值hl 依次如下:zA = 1/10 -1 11 5 2 2 21 7 3 -3 13 8 4 4 41 9 hl = 1/10 11/5 -1601/10 3696/5因此,A 即不是正定矩阵,也不是对称矩阵.(2)在MATLAB 工作窗口输入程序>> A=[1 -1 2 1;-1 3 0 -3;2 0 9 -6;1 -3 -6 19],hl=zddc(A)运行后输出结果A = 1 -1 2 1 -1 3 0 -3 2 0 9 -6 1 -3 -6 19 请注意: A 是对称矩阵请注意:因为A 的各阶顺序主子式hl 都大于零,所以A 是正定的.A 的转置矩阵zA 和各阶顺序主子式值hl 依次如下:zA = 1 -1 2 1 -1 3 0 -3 2 0 9 -6 1 -3 -6 19 hl = 1 2 6 24 (3)在MATLAB 工作窗口输入程序>> A=[1/sqrt(2) -1/sqrt(2) 0 0; -1/sqrt(2) 1/sqrt(2) 0 0; 0 01/sqrt(2) -1/sqrt(2); 0 0 -1/sqrt(2) 1/sqrt(2)], hl=zddc (A) 运行后输出结果A= 985/1393 -985/1393 0 0 -985/1393 985/1393 0 0 0 0 985/1393 -985/1393 0 0 -985/1393 985/1393 请注意: A 是对称矩阵请注意:因为A 的各阶顺序主子式hl 不全大于零,所以A 不是正定的.A 的转置矩阵zA 和各阶顺序主子式值hl 依次如下:zA = 985/1393 -985/1393 0 0 -985/1393 985/1393 0 0 0 0 985/1393 -985/1393 0 0 -985/1393 985/1393 hl = 985/1393 0 0 0可见,A 不是正定矩阵,是半正定矩阵;因为A = A T因此,A 是对称矩阵.(4)在MATLAB 工作窗口输入程序>> A=[-2 1 1;1 -6 0;1 0 -4];hl=zddc (A)运行后输出结果A = -2 1 11 -6 01 0 -4请注意:A是对称矩阵请注意:因为A的各阶顺序主子式hl不全大于零,所以A不是正定的.A的转置矩阵zA和各阶顺序主子式值hl依次如下:zA = -2 1 1 hl = -2 11 -381 -6 01 0 -4可见A不是正定矩阵,是负定矩阵;因为A= A T因此,A是对称矩阵.。