MATLAB之高斯消元法

第一种方法：文件名gauss1.m

function x=gauss1(A,b)
%参量说明：A,系数矩阵；B，常数列向量；zg，增广矩阵
%将增广矩阵化为上三角，再回带求解x
%此方法较为常规，将zg(k,k)元素乘以-zg(i,k)/zg(k,k)加到第i行
%从1:n-1列，主对角元素的以下行，通过两层循环来遍历
zg=[A,b];
zj=rref(zg);
n=length(b);
ra=rank(A);
rz=rank(zg);
temp=rz-ra;
if temp>0
    disp('无解');
    return;
end
if ra==rz
    if ra==n
        x=zeros(n,1);
        for p=1:n-1
            for k=p+1:n
                m=zg(k,p)/zg(p,p);
                zg(k,p:n+1)=zg(k,p:n+1)-m*zg(p,p:n+1);
            end
        end
        b=zg(1:n,n+1);
        A=zg(1:n,1:n);
        x(n)=b(n)/A(n,n);
        for q=n-1:-1:1
            x(q)=(b(q)-sum(A(q,q+1:n)*x(q+1:n)))/A(q,q);
        end
    end
end

第二种方法：文件名gauss2.m

function x=gauss2(A,b)
%参量说明：A,系数矩阵；B，常数列向量；zg，增广矩阵
%将增广矩阵化为上三角，再回带求解x
%该方法为：zg(k,k)右下角的n-1阶方阵为第k列除zg(k,k)元素外的列向量除以zg(k,k)
%乘以第K行除了zg(k,k)元素外的行向量的值
n=length(b);
zg=[A,b];
for k=1:(n-1)
    zg((k+1):n,(k+1):(n+1))= zg((k+1):n,(k+1):(n+1))-zg((k+1):n,k)/zg(k,k)*zg(k,(k+1):(n+1));
    zg((k+1):n,k)=zeros(n-k,1);
end

x=zeros(n,1);
x(n)=zg(n,n+1)/zg(n,n);
for k=n-1:-1:1
    x(k,:)=(zg(k,n+1)-zg(k,(k+1):n)*x((k+1):n))/zg(k,k);
end

脚本文件：gaussscript.m

clc;
%test1
A=[2,4,-6;1,5,3;1,3,2];
b=[-4;10;5];
x=gauss1(A,b)
x=gauss2(A,b)
%test2
A=[2,3,4;3,5,2;4,3,30];
b=[6,5,32]';
x=gauss1(A,b)
x=gauss2(A,b)

运行结果：

 

gaussscript

x =

-3
 2
 1

x =

-3
 2
 1

x =

-13
 8
 2

x =

-13
 8
 2