clear all
%a1=[1,-1,2,-1;2,-2,3,-3;1,1,1,0;1,-1,4,3];%系數矩陣a1
%b=[-8,-20,-2,4]';
%a1=[0.012,0.01,0.167;1,0.8334,5.91;3200,1200,4.2];
%b=[0.6781,12.1,981]';
a1=[2,-2,5;2,3,1;-1,4,1];
b=[6,13,3]';
n=length(b);%方程組大小為n
a=zeros(n,n+1);
a(:,1:n)=a1;
a(1:n,n+1)=b;%增廣矩陣
%%%%%%%%%%%%%%%%%%%%%%%%%%%%消元過程%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:n-1
    if(a(i:n,i)==0)
       error('方程組沒有唯一解'); 
    end
    for p=i:n
        if(a(p,i)~=0)
            p;%找到a的第i列第一個不是0的行號
            break;
        end
    end 
    if(p~=i)  %如果p不等于i，則交換第p和第i行
      t=a(i,:);
      a(i,:)=a(p,:);
      a(p,:)=t;
    end
    for j=i+1:n
        a(j,i)=a(j,i)/a(i,i);%為節省存儲空間，將比例因子存于a的下三角處
        a(j,:)=a(j,:)-a(j,i)*a(i,:);
    end
end
if (a(n,n)==0)
     error('方程組沒有唯一解');
end
b=a(:,n+1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%回代過程%%%%%%%%%%%%%%%%%%%%%%%%%%%%
x(n)=b(n)/a(n,n);
for i=n-1:-1:1
    sum=0;
    for j=i+1:n
        sum=sum+a(i,j)*x(j);
    end
    x(i)=(b(i)-sum)/a(i,i);
end
jie=x'