解线性方程组的列主元素高斯消去法和LU分解法实验报告.docVIP

解线性方程组的列主元素高斯消去法和LU分解法 一、实验目的：通过数值实验，从中体会解线性方程组选主元的必要性和LU分解法的优点，以及方程组系数矩阵和右端向量的微小变化对解向量的影响。 二、实验内容：解下列两个线性方程组 （1） （2） 三、实验要求： （1） 用你熟悉的算法语言编写程序用列主元高斯消去法和LU分解求解上述两个方程组，输出Ax=b中矩阵A及向量b, A=LU分解的L及U，detA及解向量x. （2） 将方程组（1）中系数3.01改为3.00，0.987改为0.990，用列主元高斯消去法求解变换后的方程组，输出列主元行交换次序，解向量x及detA，并与（1）中结果比较。 （3） 将方程组（2）中的2.099999改为2.1，5.900001改为5.9，用列主元高斯消去法求解变换后的方程组，输出解向量x及detA，并与（1）中的结果比较。 （4）用MATLAB的内部函数inv求出系数矩阵的逆矩阵，再输入命令x=inv(A)*b,即可求出上述各个方程组的解，并与列主元高斯消去法和LU分解法求出的解进行比较，体会选主元的方法具有良好的数值稳定性。用MATLAB的内部函数det求出系数行列式的值，并与（1）、（2）、（3）中输出的系数行列式的值进行比较。 四、实验过程： （1）列主元高斯消去法的主程序为 function [RA,RB,n,X]=liezhuY(A,b) B=[A b]; n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; D=det(A) if zhica0, disp(请注意：因为RA~=RB，所以此方程组无解.) return end if RA==RB if RA==n disp(请注意：因为RA=RB=n，所以此方程组有唯一解.) X=zeros(n,1); C=zeros(1,n+1); for p= 1:n-1 [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:n m= B(k,p)/ B(p,p); B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1); end end b=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n); for q=n-1:-1:1 X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q); end else disp(请注意：因为RA=RBn，所以此方程组有无穷多解.) end end 解方程组（1） 在MATLAB工作窗口输入 A=[3.01 6.03 1.999;1.27 4.16 -1.23;0.987 -4.81 9.34];b=[1;1;1];[RA,RB,n,X]=liezhuY(A,b) 运行后输出结果为 请注意：因为RA=RB=n，所以此方程组有唯一解. D=-0.1225 RA =3 RB =3 n =3 X = 397.8654 -157.6242 -123.1120 解方程组（2） 在MATLAB工作窗口输入 A=[10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 0 2];b=[8;5.900001;5;1];[RA,RB,n,X]=liezhu(A,b) 运行后输出结果为 请注意：因为RA=RB=n，所以此方程组有唯一解. D=-762.0000 RA =4 RB =4 n =4 X =0.0000 -1.0000 1.0000 1.0000 LU分解法及MATLAB主程序为 function hl=zhjLU(A) [n n] =size(A); RA=rank(A); D=det（A） if RA~=n disp(请注意：因为A的n阶行列式hl等于零，所以A不能进行LU分解.A的秩RA如下:), RA,hl=det(A); return end if RA==n for p=1:n h(p)=det(A(1:p, 1:p)); end hl=h(1:n); for i=1:n if h(1,i)==0 disp(请注意：因为A的r阶主子式等于零，所以A不能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下：), hl;RA return end end if h(1,i)~=0 disp(请注