piveli.m

fit programme for a serials of data

function X=piveli(A,b,delta)
% input A is an N*N matrix;
%       b is an N*1 matrix;
% output X is the solution of AX=b
n=length(b);
augA=[A,b];
t=zeros(1,n+1);
for k=1:n-1
   [N,j]=max(abs(augA(k:n,k)))
   if N<delta
      disp(sprintf('detA=0'));
      return;
   end
   if j~1
      t=augA(k,:);
      augA(k,:)=augA(k+j-1,:);
      augA(k+j-1,:)=t;
      augA
   end
   for j=k+1:n
      augA(j,:)=augA(j,:)*augA(k,k)/augA(j,k)-augA(k,:); 
   end
   augA
end
if augA(n,n)==0
   disp(sprintf('Matrix A is singular'));
   return;
end
for k=n:-1:2
   for j=1:k-1
       augA(j,:)=augA(j,:)*augA(k,k)/augA(j,k)-augA(k,:);
   end
   augA
end
for k=1:n
   X(k)=augA(k,n+1)/A(k,k) 
end