gauss.m

来自「matlab具有强大的数值计算能力」· M 代码 · 共 60 行

function [x,det,index]=Gauss(A,b)
%求线性方程组的列主元素,高斯消去法,其中:
%A为系数矩阵
%b为方程组的右端项
%x为方程组的解
%index为指标变量
[n,m]=size(A);nb=length(b);
if n~=m
    error('the rows and columns of matrix A must be equal!');
    return;
end
%当程组右项维数匹时,停止计算
if m~=nb
    error('the columns of A must be equal the length of b!')
    return;
end 
%开始计算先初值
index=1;det=1;x=zeros(n,1);
for k=1:n-1
    %选主元
    a_max=0;
    for i=k:n
        if abs(A(i,k))>a_max
            a_max=abs(A(i,k));r=i;
        end
    end
    if a_max<1e-10
        index=0;return;
    end
    %交换两行
    if r>k
        for j=k:n
            z=A(k,j);A(k,j)=A(r,j);A(r,j)=z;
        end
        z=b(k);b(k)=b(r);b(r)=z;det=-det;
    end
    %消元计数
    
    for i=k+1:n
       m=A(i,k)/A(k,k);
   
            for j=k+1:n
                A(i,j)=A(i,j)-m*A(k,j);
            end
            b(i)=b(i)-m*b(k);
        end
        det=det*A(k,k);
    end
    det=det*A(n,n);
    
%输出X
if abs(A(n,n))<1e-10
    index=0;return;
end
for k=n:-1:1 
    for j=k+1:n
       b(k)=b(k)-A(k,j)*x(j);
    end 
    x(k)=b(k)/A(k,k);
end