演算法解釋見此：https://blog.csdn.net/zengxyuyu/article/details/53056453

原始碼在此：

main.m

clear
clc
A = [8 -3 2;4 11 -1;6 3 12];
b = [20;33;36];
[x, n] = jacobi(A,b,[0,0,0]',1.0e-7,30)

guaseidel.m

​
function [x,n] = guaseidel(A,b,x0,eps,it_max)
% 求線性方程組的Gauss-Seidel迭代法，呼叫格式為

%  [x, k] = guaseidel(A,b,x0,eps,it_max)

%  其中, A 為線性方程組的係數矩陣，b 為常數項，eps 為精度要求，預設為1e-5,

%  it_max 為最大迭代次數，預設為100

%  x 為線性方程組的解，k迭代次數
if nargin == 3
    eps = 1.0e-6
    it_max= 200
elseif nargin == 4
    it_amx = 200
elseif nargin <3
    disp('輸入引數個數不足3個');
    return;
end
D = diag(diag(A));%求A的對角矩陣
L = -tril(A,-1);%求A的下三角矩陣,不帶對角線
U = -triu(A,1);%求A的上三角矩陣
G = (D-L)\U;
f = (D-L)\b;
x = G*x0+f;
n=1;  %迭代次數
while norm(x-x0)>=eps
    x0 = x;
    x = G*x0+f;
    n = n+1;
    if(n>=it_max)
        disp('Warning:迭代次數太多,可能不收斂');
        return;
    end
end

​