1、高斯消元法解线性方程组实验一、高斯消元法解线性方程组流程图:K= abs(A(z,j) c=i; i=z; z=c; else; end end if z=j ; else q=A(z,:); A(z,:)=A(j,:); A(j,:)=q; end if abs(A(j,j)= abs(A(z,j) c=i; i=z; z=c; else; end end if z=j ; else q=A(z,:); A(z,:)=A(j,:); A(j,:)=q; pc=p(z,:); p(z,:)=p(j,:); p(j,:)=pc; end if abs(A(j,j)=e error(方程系数矩阵为
2、非满秩矩阵); else; end endp;for j=1:Al-1 for i=j+1:Ah d=-A(i,j)/A(j,j); A(i,:)=d*A(j,:)+A(i,:); endendU=A(1:Ah,1:Al-1);A=Ae;L=(p*A)/U;pALUB=p*BW=size(L);%求解Ly=bAh=W(1);Al=W(2);Y=zeros(Al,1);for i =1:Ah S=0; for j=1:i-1 S=S+L(i,j)*Y(j,1); end Y(i,1)=(B(i,1)-S)/L(i,i);endYW=size(U);%求解Ux=yAh=W(1);Al=W(2);X
3、=zeros(Al,1);for i =Ah:-1:1 S=0; for j=Al-1:-1:i S=S+U(i,j+1)*X(j+1,1); end X(i,1)=(Y(i,1)-S)/U(i,i);endX运行结果: 实验二、直接三角法解线性方程组Li=input(请输入你要解的方程的增广矩阵的列数:);e=input(请输入你要求的非满秩矩阵的最小元绝对值:);fid=fopen(Untitled.m,r);A=fscanf(fid,%f,Li,inf);A=A;W=size(A);H=W(1);AZ=A(1:H,1:H);B=A(1:H,H+1);A=AZ;for j=1:H %求取矩
4、阵U的第一行 if abs(A(j,j)=e error(方程系数矩阵为非满秩矩阵,请选择杨雁勇的下一个可以选择列主元的程序); else; endendfor j=1:H %求取矩阵U的第一行 A(1,j)=A(1,j);endfor i=2:H %求取矩阵L的第1列 A(i,1)=A(i,1)/A(1,1);endfor r=2:H %求取矩阵U的第2-H列和第2-H行 for i=r:H A(r,i)= A(r,i)-A(r,1:r-1)*A(1:r-1,i); end for i=r+1:H A(i,r)=(A(i,r)-A(i,1:r-1)*A(1:r-1,r)/A(r,r); en
5、d if abs(A(r,r)j L(i,j)= A(i,j); else L(i,j)= 0; end end endendLfor i=1:H %将A拆成L与U for j=1:H if i=j U(i,j)= A(i,j); else U(i,j)= 0; end endendU %回带求取解W=size(L);%求解Ly=bAh=W(1);Al=W(2);Y=zeros(Al,1);for i =1:Ah S=0; for j=1:i-1 S=S+L(i,j)*Y(j,1); end Y(i,1)=(B(i,1)-S)/L(i,i);endYW=size(U);%求解Ux=yAh=W(
6、1);Al=W(2);X=zeros(Al,1);for i =Ah:-1:1 S=0; for j=Al-1:-1:i S=S+U(i,j+1)*X(j+1,1); end X(i,1)=(Y(i,1)-S)/U(i,i);endX运行结果: 实验三、列主元三角法解线性方程组Li=input(请输入你要解的方程的增广矩阵的列数:);e=input(请输入你要求的非满秩矩阵的最小元绝对值:);fid=fopen(Untitled.m,r);A=fscanf(fid,%f,Li,inf);A=A;W=size(A);H=W(1);AZ=A(1:H,1:H);B=A(1:H,H+1);A=AZ;K
7、=1;for i=1:H K=K*A(i,i);endif K=0 error(方程中有对角元上元素=0); else;endfor j=1:H %求取矩阵U的第一行 A(1,j)=A(1,j);endfor i=2:H %求取矩阵L的第1列 A(i,1)=A(i,1)/A(1,1);endfor r=2:H %求取矩阵U的第2-H列和第2-H行%检验A(r,r)是否为一个特别小的数 if abs(A(r,r)= abs(S(z,1) c=j; j=z; z=c; else; end end if z=r ; else q=A(z,r:H); A(z,r:H)=A(r,r:H); A(r,r:
8、H)=q; end else; end for i=r:H A(r,i)= A(r,i)-A(r,1:r-1)*A(1:r-1,i); end for i=r+1:H A(i,r)=(A(i,r)-A(i,1:r-1)*A(1:r-1,r)/A(r,r); end if abs(A(r,r)j L(i,j)= A(i,j); else L(i,j)= 0; end end endendLfor i=1:H %将A拆成L与U for j=1:H if ij L(i,j)= -A(i,j); else L(i,j)= 0; end endendLfor i=1:H %将A拆成L与U for j=1
9、:H if ij U(i,j)= -A(i,j); else U(i,j)= 0; end endendU X=zeros(H,1); for j=1:1000 S=X; for i=1:H a=X(i,1); X(i,1)=(B(i,1)-A(i,1:i-1)*X(1:i-1,1)-A(i,i+1:H)*X(i+1:H,1)/A(i,i); X(i,1)=a+w*( X(i,1)-a); end disp(X) R=S-X; c=abs(R(1,1); for i=1:H if c=abs(R(i,1) c=abs(R(i,1); else ; end end if c=e break; else; end end运行结果:实验体会: 通过实验我掌握了消元法解方程的一些基本算法以及用matlab实现矩阵的几种基本计算。对MATLAB软件有了更深的了解。 同时,在实验中,在输入向量时错误也很多,比如将:“(”写成“”等,致使程序不能运行,无法得到预期的实验结果,有时候还会在结束循环体是落掉“end”等,通过此次练习,我基本上熟悉了Matlab编程的基技能,希望在以后的学习中还能有更多的机会运用它。
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1