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);
end










end

To avoid $a_{kk}$ is too small lead to big error.