解线性方程组的高斯消去法

clc;
clear all;
disp('此为列主元消去法');
A=input('矩阵A=');
b=input('右端项b=');
eps=input('控制精度eps=');
[n,n]=size(A);
B=[A, b]; %%矩阵增广
ra=rank(A); %求系数矩阵的秩
rb=rank(B); %求增广矩阵的秩
%判断方程解的情况
if ra~=rb
    disp('系数矩阵的秩与增广矩阵的秩不一致,此方程无解')
    return;
elseif (ra==rb)&&(ra==n)
    disp('此方程组有唯一解');
else disp('此方程组有无穷多解');
end
% end
disp('矩阵A的行列式为:')
det(A)
t=zeros(n+1)
for k=1:1:n-1
    %找出每一列的最大值
    max=k;
    for j=k+1:1:n
        if abs(A(j,k))>abs(A(max,k))
            max=j;
        end
    end
    %执行换行
    t(1,1)=b(k,1);
    b(k,1)=b(max,1);
    b(max,1)=t(1,1);
    w(1,:)=A(k,:);
    A(k,:)=A(max,:);
    A(max,:)=w(1,:);
    %消元过程
    for i=k+1:1:n
        m(i,k)=A(i,k)/A(k,k);
        for j=k+1:1:n
            A(i,j)=A(i,j)-A(k,j)*m(i,k);
        end
        b(i)=b(i)-m(i,k)*b(k);
    end
end
%回代过程
x=zeros(n,1);
x(n)=b(n)/A(n,n);
for i=n-1:-1:1
    x(i)=(b(i)-A(i,i+1:n)*x(i+1:n))/A(i,i);
end
x

完全主元高斯消去法

function  x = GuassFull( A,b)
%说明A是方阵,b是列向量
%example:
%A=[1,2,3;4,8,9;10,11,12];
%B = [2,3,5];
%x =GuassFull(A,B');
% Created By wangxiaole 14.9.18
% This function use Guass methods to calculate Linear Equations
% Choose the maximum element in all Matrix A;
[m,n] = size(A);
if rank(A) == m
    if  m == 1 %只有一个方程
        if A ~=0
            x = b/A;
        else
            x = inf;
        end
    else
        
        MaxEle = max(max(abs(A)));
        [m1,n1] = find(abs(A) == MaxEle);
        Row = A(:,n1(1));
        A(:,n1(1))=[];
        A = [Row,A];
        Column  =A(m1(1),:);
        A(m1(1),:) =[];
        A = [Column;A];
        
        B_Column = b(m1(1));
        b(m1(1)) = [];
        b =[B_Column ;b];
        
        for i = 2:m
            temp  = A(i,1);
            
            A(i,:) = A(i,:) + -1*temp/A(1,1)*A(1,:);
            b(i)= -1*temp/A(1,1)*b(1) +b(i);
        end
        A1 = A;
        B1 = b;
        
        A1(1,:)=[];
        A1(:,1)=[];
        B1(1)= [];
        x = GuassFull(A1,B1);
        tempx = (sum(-1*A(1,2:end).*x') +b(1))/A(1,1);
        if(x == inf)
            x = inf;
        else
            if Row == 1
                x = [tempx;x];
            else
                x = [x(1:n1(1)-1);tempx;x(n1(1):end)];
            end
        end
    end
else
    x = inf;
end
end

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值