MATLAB实现Jacobi迭代

线性方程组迭代方法

用Jacobi迭代计算如下方程组

推导过程

编程思路

在求解线性方程组时， 程序中输入的是线性方程组的增广矩阵，并设定送
代结果的精度要求， 在开始迭代之前则可以进行以下判断:

1. 检测该用户输入矩阵是否有误，并判断该线性方程组是否有唯一解;
2. 判断系数矩阵的对角元是否有零元素。
   由于迭代过程的本质是矩阵乘法，因此程序中直接调用 MATLAB 的矩阵乘
   法来完成送代，选代均以零向量为初始向量，并且在送代之前判断该送代过程是
   否收敛。 倘若收敛，则输出结果并打印迭代次数。

MATLAB代码

function函数：

function solution = Jacobi(Ab,epsilon)
%% Ab为输入的增广矩阵
%% epsilon为输入的精度要求

%% 输入参数检查
if nargin ==1
    disp('请输入精度要求epsilon')
    return 
end

row = size(Ab,1);
col = size(Ab,2);
if ndims(Ab)~= 2 | col - row ~= 1
    disp('矩阵的大小有误，不能采用Jacobi迭代法')
    return
end

A = Ab(:,1:row);
b = Ab(:,col);

ddet = abs(det(A));
ddiag = abs(det(diag(diag(A))));
if ddet < eps | ddiag < eps
    disp('该方程的系数矩阵行列式不为0，方程组无解或有无穷解，或系数矩阵的对角线有零元，不能采用Jacobi迭代法')
    return
end

%%提取上下三角矩阵及对角矩阵
U = -triu(A,1);   %提取对角元素为0的上三角矩阵
L = -tril(A,-1);  %提取对角元素为0的下三角矩阵
Dinv = diag(diag(A).^(-1));

B = Dinv * (L + U);
f = Dinv * b;
if max(abs(eig(B))) > 1
    disp('迭代法不收敛！！！')
    return
end
%% 迭代过程
error = 10;
n = 0;
start = zeros(row,1);
xk = start;
xknext = B * xk +f;
error = norm(xk - xknext);
while error > epsilon
    xk = xknext;
    xknext = B * xk + f;
    error = norm(xk - xknext);
    n = n+1;
end
fprintf('Jacobi迭代次数：%d',n)
solution = xknext;

主函数：

Ab = [10 -2 -1 3;-2 10 -1 15;-1 -2 5 10]
Jacobi(Ab,1e-8)

**注意：**function函数和主函数是在同一个文件夹路径下的两个文件！