雅可比求解器进入无限循环

Jacobi solver going into an infinite loop

我似乎找不到解决我的无限循环的方法。我编写了一个 Jacobi 求解器来求解线性方程组。

这是我的代码:

function [x, i] = Jacobi(A, b, x0, TOL)
[m n] = size(A);
i = 0;
x = [0;0;0];
while (true)
    i =1;

     for r=1:m
        sum = 0;
        for c=1:n
            if r~=c
                sum = sum + A(r,c)*x(c);
            else
                x(r) = (-sum + b(r))/A(r,c);
            end    
            x(r) = (-sum + b(r))/A(r,c);
 xxx   end                                       xxx
    end
    if abs(norm(x) - norm(x0)) < TOL;
        break
    end
    x0 = x;
    i = i + 1;
end

当我终止代码时,它在带有 xxx

的行结束

您的代码无法正常工作的原因是您的 if 语句在 for 循环中的逻辑。具体来说,您需要先累加 属于该行对角线的特定行的所有值。完成后,您 然后 执行除法。您还需要确保除以您正在关注的那一行的 A 的对角线系数,这对应于您要解决的 x 的组件。您还需要删除循环开头的 i=1 语句。您每次都在重置 i

换句话说:

function [x, i] = Jacobi(A, b, x0, TOL)
[m n] = size(A);
i = 0;
x = [0;0;0];
while (true)
     for r=1:m
        sum = 0;
        for c=1:n
            if r==c %// NEW
                continue;
             end
            sum = sum + A(r,c)*x(c); %// NEW
        end                                      
       x(r) = (-sum + b(r))/A(r,r); %// CHANGE
    end
    if abs(norm(x) - norm(x0)) < TOL;
        break
    end
    x0 = x;
    i = i + 1;
end

使用示例:

A = [6 1 1; 1 5 3; 0 2 4]
b = [1 2 3].';
[x,i] = Jacobi(A, b, [0;0;0], 1e-10)

x =

   0.048780487792648
  -0.085365853612062
   0.792682926806031


i =

    20

这意味着需要 20 次迭代才能获得具有公差 1e-10 的解决方案。将其与 MATLAB 的内置逆函数进行比较:

x2 = A \ b

x2 =

   0.048780487804878
  -0.085365853658537
   0.792682926829268

如您所见,我指定了 1e-10 的公差,这意味着我们保证有 10 位小数的准确性。我们当然可以看到 Jacobi 给我们的和 MATLAB 给我们内置的之间有 10 位小数的精度。