gaussquanzhuyuan.m

来自「这是线性议程组的直接解法」· M 代码 · 共 37 行

function [x,q]=Gaussquanzhuyuan(A,b)
B=[A,b];
if det(A)==0
    disp('方程组不存在唯一解！');
    return;
end
n=length(A);
t=[1:n];
for k=1:n
   m=max(max(abs(B(k:n,k:n))));
   [h,l]=find(B(k:n,k:n)==m);
   if h>1
       temp=B(k,:);
       B(k,:)=B(h+k-1,:);
       B(h+k-1,:)=temp;
   end
   if l>1
       temp=B(:,k);
       B(:,k)=B(:,l+k-1);
       B(:,l+k-1)=temp;
       temp=t(k);
       t(k)=t(l+k-1);
       t(l+k-1)=temp;
   end
   for kk=k+1:n
    B(kk,:)=(-B(kk,k)*B(k,:)/B(k,k)+B(kk,:));
   end
end
x(n,1)=B(n,n+1)/B(n,n);
for k=n-1:-1:1
    sx=B(k,k+1:n)*x(k+1:n,1);
    x(k,1)=(B(k,n+1)-sx)/B(k,k);
end
x=sortrows([t',x],1);
x=x(:,2)
q=(n^3)/3+n^2-1/3*n;
end