gauss 2

function    [A,b,x]=gauss2(A,b)%--------------------------------------------------------------n=size(A,1);for k=1:(n-1)  max1=max(abs(A(k:n,k)));%find the max of kth to nth  index=find(A(k:n,k)==max1);  index=index+k-1;         %find the index of the row of the max   if k~=index  t1=A(k,:);  A(k,:)=A(index,:);   %interchange the kth and the index row  A(index,:)=t1;  t2=b(k);  b(k)=b(index);  b(index)=t2;  end%------------------------------------------------------------for i=k+1:n             m(i,k)=A(i,k)/A(k,k);%the number a_{ik}/a_{kk}           for j=1:n                 A(i,j)=A(i,j)-A(k,j)*m(i,k);            end                b(i)=b(i)-b(k)*m(i,k);        end end x=zeros(n,1);x(n)=b(n)/A(n,n);for k=n-1:-1:1    x(k)=(b(k)-(A(k,k+1:n)*x(k+1:n)))/A(k,k);endend

To avoid akk is too small lead to big error.