现在时间是:
当前位置:首 页 >> 矩阵计算>> 教学区>> 文章列表

130213025+王颖洁+作业3+Gauss消去法列选全选主元

作者:130213025王颖洁   发布时间:2015-06-28 17:42:53   浏览次数:150
全选主元Gauss消去法
function x=gauss(A,b)
[m,n]=size(A);
p=[];
q=[];
tol=1e-14;
x=zeros(n,1);
if m==n
   for k=1:n-1
       [r,c]=find(A==max(max(abs(A(k:n,k:n)))));
       p(k)=r(1);
       q(k)=c(1);
       A([k,p(k)],:)=A([p(k),k],:);
       A(:,[k,q(k)])=A(:,[q(k),k]);
       if A(k,k)>tol
           A(k+1:n,k)=A(k+1:n,k)/A(k,k);
           A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n);
       else
           break;
       end
   end
L=tril(A,-1)+eye(n)
U=triu(A,0)
for k=1:length(p)
   b([k,p(k)])=b([p(k),k]);
end
b(1)=b(1)/L(1,1);
for j=2:n
  b(j)=(b(j)-L(j,1:j-1)*b(1:j-1))/L(j,j);
end
x(n)=b(n)/U(n,n);
for j=n-1:-1:1
  x(j)=(b(j)-x(j+1:n)*U(j,j+1:n))/U(j,j);
end
for k=length(q):-1:1
   x([k,q(k)])=x([q(k),k]);
end
else
   'A is not a square!'
end
 
运行结果:
A=[0.001 1.00;1.00 2.00];
b=[1.00;3.00];
x=gauss(A,b)
err=A*x-b
 
L =
   1.0000         0
   0.5000    1.0000
U =
   2.0000    1.0000
        0   -0.4990
x =
   1.0020
   0.9990
err =
    0
    0
 
 
列选主元Gauss消去法
function x=gauss(A,b)
[m,n]=size(A);
p=[];
tol=1e-14;
x=zeros(n,1);
if m==n
   for k=1:n-1
       [r,c]=find(A==max(max(abs(A(k:n,k:n)))));
       p(k)=r(1);
       A([k,p(k)],:)=A([p(k),k],:);
       if A(k,k)>tol
           A(k+1:n,k)=A(k+1:n,k)/A(k,k);
           A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n);
       else
           break;
       end
   end
L=tril(A,-1)+eye(n)
U=triu(A,0)
for k=1:length(p)
   b([k,p(k)])=b([p(k),k]);
end
b(1)=b(1)/L(1,1);
for j=2:n
  b(j)=(b(j)-L(j,1:j-1)*b(1:j-1))/L(j,j);
end
x(n)=b(n)/U(n,n);
for j=n-1:-1:1
  x(j)=(b(j)-x(j+1:n)*U(j,j+1:n))/U(j,j);
end
else
   'A is not a square!'
end
 
运行结果:
A=[0.001 1.00;1.00 2.00];
b=[1.00;3.00];
x=gauss(A,b)
err=A*x-b
 
L =
   1.0000         0
   0.0010    1.0000
U =
   1.0000    2.0000
        0    0.9980
x =
   1.0020
   0.9990
err =
    0

    0







上一篇:没有了    下一篇:没有了

Copyright ©2018    计算数学达人 All Right Reserved.

技术支持:自助建站 | 领地网站建设 |短信接口 版权所有 © 2005-2018 lingw.net.粤ICP备16125321号 -5