数值计算解线性方程组

高斯消去法

function result = Gauss(n,A,b)
%GAUSS 此处显示有关此函数的摘要
%   此处显示详细说明


result=[];
 
for j = 1:n-1
    for i = j+1:n
        m = A(i,j)/A(j,j);
        A(i,:)= A(i,:)-m*A(j,:);
        b(i,1)= b(i,1)-m*b(j,1);
    end
end
 
 
for i=n:-1:1
    sum=0;
    for j=n:-1:i+1
        sum=sum+result(j,1)*A(i,j);
    end
    result(i,1)=(b(i,1)-sum)/A(i,i);
end


disp(result)
end

//main
C=zeros(100,100);
for i=1:99
     C(i,i+1)=-1;
     C(i+1,i)=-1;
     C(i,i)=3; 
end
C(100,100)=3;
 
d=ones(100,1);
d(1,1)=2;
d(100,1)=2;
Gauss(100,C,d);

jacobi迭代

function x=jacobi(A,b,n,x0,eps,N)
D = diag(diag(A));
L = -tril(A,-1);
U = -triu(A,1);
B = D\(L+U);
if(max(abs(eig(B))>1))
    disp('不收敛')
    return
end
x=zeros(n,1);
k=0;
while k<N
    for i = 1:n
        s = b(i);
        for j=1:i-1
            s = s-A(i,j)*x0(j);
        end
        for j=i+1:n
            s = s-A(i,j)*x0(j);
        end
        x(i) = s/A(i,i);
    end
    if norm(x-x0)<eps
        break;
    end
    x0=x;
    k=k+1;  
end
if k>N
    disp('达到最大循环次数');
end
end

//main
A = [2 -2 -1;
    4 1 -2;
    -2 1 -1];
b = [-2;1;-3];

Gauss(3,A,b);

C=zeros(100,100);
for i=1:99
     C(i,i+1)=-1;
     C(i+1,i)=-1;
     C(i,i)=3; 
end
C(100,100)=3;
 
d=ones(100,1);
d(1,1)=2;
d(100,1)=2;

Gauss(100,C,d);

disp('jacobi what')
x0=[0;0;0];
jacobi(A,b,3,x0,10^-3,100)
disp('jacobi what')
x1=zeros(100,1);
jacobi(C,d,100,x1,10^-3,100)

disp('gswhat')
gs(A,b,3,x0,10^-3,100)
disp('gswhat')
gs(C,d,100,x1,10^-3,100)

G-S 迭代

function x=gs(A,b,n,x0,eps,N)
D = diag(diag(A));
L = tril(A,-1);
U = triu(A,1);
B = -(D+L)\U;
f=(diag(diag(A))+tril(A,-1))\b;

if(max(abs(eig(B))>1))
    disp('不收敛')
    return
end
x=zeros(n,1);
k=0;
while k<N
    for i = 1:n
        s = b(i);
        for j=1:i-1
            s = s-A(i,j)*x(j);
        end
        for j=i+1:n
            s = s-A(i,j)*x0(j);
        end
        x(i) = s/A(i,i);
    end
    if norm(x-x0)<eps
        break;
    end
    x0=x;
    k=k+1;  
end
if k>N
    disp('达到最大循环次数');
disp(x);
end

高斯列主元

%列主元消去法解线性方程组
clear;
clc;
%fprintf('列主元消去法解线性方程组:\n')
n=3;
A = [2 -2 -1;
    4 1 -2;
    -2 1 -1];
b = [-2;1;-3];

      %注意，这里我用的S是行向量，计算时注意转置
%n=input('请输入系数矩阵A的阶数:n=');
%A=input('请输入系数矩阵A:A=');
%b=input('请输入结果向量b:b=');%b为行向量，计算时需要转置
%消元计算---->产生上三角矩阵：
A1=A;                         %为了保留原始参数，以下用A1存储产生的上三角矩阵
b1=b;
for j=1:n                     %从第一列到第n列
    for i=j+1:n               %每一行从主元素的下侧第一个元素开始
        for k=i:n%这个循环选出第j列最大的元素做主元
            if abs(A1(k,j))>abs(A1(j,j))%if语句，别忘了绝对值，ps：冒泡排序，真滴好用
                for k1=j:n              %用来交换第k行和第j行
                    a1=A1(j,k1);
                    a2=A1(k,k1);
                    A1(j,k1)=a2;
                    A1(k,k1)=a1;
                end
                t=b1(k);
                b1(k)=b1(j);
                b1(j)=t;%别忘了结果向量，也要相应的交换行
            end
        end
%上述增补的部分就是列主元消去的精髓，下面的部分和高斯消去一样
        if A1(j,j)~=0         %每次的主元素都不可为零，这是能使用高斯消去法的前提
            M(i,j)=A1(i,j)/A1(j,j);
            for k=j:n         %第i行j列开始计算到第i行n列
                A1(i,k)=A1(i,k)-M(i,j)*A1(j,k);
            end
                b1(i)=b1(i)-M(i,j)*b1(j);   
        else 
            error('求解出错，高斯消去法不适用于此方程组')
        end
    end
end
%回代求解：
 x(n)=b1(n)/A1(n,n);    %先求出解向量的最后一个元素，用于回代
 i=n-1;                 %从n-1行开始，回代求解
while i>0               %从第n-1行到第一行，依次回代
     s=0;               %对s赋初值0，
     for j=i+1:n
         s=A1(i,j)*x(j)+s;%计算：（已知解向量×系数矩阵对应元素）的和
     end
     x(i)=(b1(i)-s)/A1(i,i);%回代求解，解向量中已知元素+1.用于回代求x(i-1)
     i=i-1;
end
fprintf('经过消元回代，此方程组的解为:\n')
disp(x);