资源描述：
《线性方程组求解Matlab程序.doc》由会员上传分享，免费在线阅读，更多相关内容在应用文档-天天文库。
1、线性方程组求解1.直接法Gauss消元法：functionx=DelGauss(a,b)%Gauss消去法[n,m]=size(a);nb=length(b);det=1;%存储行列式值x=zeros(n,1);fork=1:n-1fori=k+1:nifa(k,k)==0returnendm=a(i,k)/a(k,k);forj=k+1:na(i,j)=a(i,j)-m*a(k,j);endb(i)=b(i)-m*b(k);enddet=det*a(k,k);enddet=det*a(n,n);fork=n:-1:1%回代forj=k+1:nb(k)=b(k)-a
2、(k,j)*x(j);endx(k)=b(k)/a(k,k);endExample:>>A=[1.0170-0.00920.0095;-0.00920.99030.0136;0.00950.01360.9898];>>b=[101]';>>x=DelGauss(A,b)x=0.9739-0.00471.0010列主元Gauss消去法：functionx=detGauss(a,b)%Gauss列主元消去法[n,m]=size(a);nb=length(b);det=1;%存储行列式值x=zeros(n,1);fork=1:n-1amax=0;%选主元fori=k:ni
3、fabs(a(i,k))>amaxamax=abs(a(i,k));r=i;endendifamax<1e-10return;endifr>k%交换两行forj=k:nz=a(k,j);a(k,j)=a(r,j);a(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;det=-det;endfori=k+1:n%进行消元m=a(i,k)/a(k,k);forj=k+1:na(i,j)=a(i,j)-m*a(k,j);endb(i)=b(i)-m*b(k);enddet=det*a(k,k);enddet=det*a(n,n);fork=n:-1:1
4、%回代forj=k+1:nb(k)=b(k)-a(k,j)*x(j);endx(k)=b(k)/a(k,k);endExample:>>x=detGauss(A,b)x=0.9739-0.00471.0010Gauss-Jordan消去法：functionx=GaussJacobi(a,b)%Gauss-Jacobi消去法[n,m]=size(a);nb=length(b);x=zeros(n,1);fork=1:namax=0;%选主元fori=k:nifabs(a(i,k))>amaxamax=abs(a(i,k));r=i;endendifamax<1e-10
5、return;endifr>k%交换两行forj=k:nz=a(k,j);a(k,j)=a(r,j);a(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;end%进行消元b(k)=b(k)/a(k,k);forj=k+1:na(k,j)=a(k,j)/a(k,k);endfori=1:nifi~=kforj=k+1:na(i,j)=a(i,j)-a(i,k)*a(k,j);endb(i)=b(i)-a(i,k)*b(k);endendendfori=1:nx(i)=b(i);endExample:>>x=GaussJacobi(A,b)x=0.9
6、739-0.00471.0010LU分解法：function[l,u]=lu(a)%LU分解n=length(a);l=eye(n);u=zeros(n);fori=1:nu(1,i)=a(1,i);endfori=2:nl(i,1)=a(i,1)/u(1,1);endforr=2:n%%%%fori=r:nuu=0;fork=1:r-1uu=uu+l(r,k)*u(k,i);endu(r,i)=a(r,i)-uu;end%%%%fori=r+1:nll=0;fork=1:r-1ll=ll+l(i,k)*u(k,r);endl(i,r)=(a(i,r)-ll)/u(
7、r,r);end%%%%Endfunctionx=lusolv(a,b)%LU分解求解线性方程组aX=biflength(a)~=length(b)error('Errorininputing!')return;endn=length(a);[l,u]=lu(a);y(1)=b(1);fori=2:nz=0;fork=1:i-1z=z+l(i,k)*y(k);endy(i)=b(i)-z;endx(n)=y(n)/u(n,n);fori=n-1:-1:1z=0;fork=i+1:nz=z+u(i,k)*x(k);endx(i)=(y(i)-z)/u(i,i);e
