Jacobi迭代不会结束 [英] Jacobi iteration doesn't end
问题描述
我正在尝试在MATLAB中实现Jacobi迭代，但是无法使其收敛.我已经在网上和其他地方查看了可以比较的工作代码，但是找不到与我的代码相似的东西，并且仍然可以使用.这是我所拥有的:
I'm trying to implement the Jacobi iteration in MATLAB but am unable to get it to converge. I have looked online and elsewhere for working code for comparison but am unable to find any that is something similar to my code and still works. Here is what I have:
function x = Jacobi(A,b,tol,maxiter)
n = size(A,1);
xp = zeros(n,1);
x = zeros(n,1);
k=0; % number of steps
while(k<=maxiter)
k=k+1;
for i=1:n
xp(i) = 1/A(i,i)*(b(i) - A(i,1:i-1)*x(1:i-1) - A(i,i+1:n)*x(i+1:n));
end
err = norm(A*xp-b);
if(err<tol)
x=xp;
break;
end
x=xp;
end
不管我用什么A和b，它都会炸毁.这可能是我忽略的一个小错误，但是如果有人能解释什么地方出了问题，我将不胜感激，因为这应该是正确的，但实际上并非如此.
This just blows up no matter what A and b I use. It's probably a small error I'm overlooking but I would be very grateful if anyone could explain what's wrong because this should be correct but is not so in practice.
推荐答案
您的代码正确.它似乎不起作用的原因是，您指定的系统在使用Jacobi迭代时可能不会收敛.
Your code is correct. The reason why it may not seem to work is because you are specifying systems that may not converge when you are using Jacobi iterations.
具体来说(由于@Saraubh)，如果您的矩阵A严格 对角线占主导地位.换句话说，对于矩阵中的每一行i，行i上所有列j的绝对总和必须小于 ，而对角线系数必须位于i.对角线本身.换句话说:
To be specific (thanks to @Saraubh), this method will converge if your matrix A is strictly diagonally dominant. In other words, for each row i in your matrix, the absolute summation of all of the columns j at row i without the diagonal coefficient at i must be less than the diagonal itself. In other words:
但是，即使不满足此条件，也有一些系统可以与Jacobi融合，但是在尝试将Jacobi用于系统之前，应将其用作一般规则.如果您使用Gauss-Seidel，它实际上会更稳定.唯一的不同是，您将重新使用x的解决方案，并随着行的向下移动将其输入到其他变量中.要制作此高斯-赛德尔(Gauss-Seidel)，您要做的就是在for循环中更改一个字符.对此进行更改:
However, there are some systems that will converge with Jacobi, even if this condition isn't satisfied, but you should use this as a general rule before trying to use Jacobi for your system. It's actually more stable if you use Gauss-Seidel. The only difference is that you are re-using the solution of x and feeding it into the other variables as you progress down the rows. To make this Gauss-Seidel, all you have to do is change one character within your for loop. Change it from this:
xp(i) = 1/A(i,i)*(b(i) - A(i,1:i-1)*x(1:i-1) - A(i,i+1:n)*x(i+1:n));
对此:
xp(i) = 1/A(i,i)*(b(i) - A(i,1:i-1)*xp(1:i-1) - A(i,i+1:n)*x(i+1:n));
**HERE**
以下是两个示例，我将向您展示:
Here are two examples that I will show you:
- 我们在这里指定了Jacobi不会收敛的系统，但是有解决方案.该系统不是对角线占主导地位.
- 我们在其中指定确实由Jacobi收敛的系统.具体来说，该系统在对角线占主导地位.
- Where we specify a system that does not converge by Jacobi, but there is a solution. This system is not diagonally dominant.
- Where we specify a system that does converge by Jacobi. Specifically, this system is diagonally dominant.
示例1
A = [1 2 2 3; -1 4 2 7; 3 1 6 0; 1 0 3 4];
b = [0;1;-1;2];
x = Jacobi(A, b, 0.001, 40)
xtrue = A \ b
x =
1.0e+09 *
4.1567
0.8382
1.2380
1.0983
xtrue =
-0.1979
-0.7187
0.0521
0.5104
现在，如果我使用Gauss-Seidel解决方案，这就是我得到的:
Now, if I used the Gauss-Seidel solution, this is what I get:
x =
-0.1988
-0.7190
0.0526
0.5103
哇！它收敛于高斯-塞德尔(Gauss-Seidel)，而不是雅可比(Jacobi)，即使该系统不是对角线主导系统，我也可能对此做出解释，我将在稍后提供.
Woah! It converged for Gauss-Seidel and not Jacobi, even though the system isn't diagonally dominant, I may have an explanation for that, and I'll provide later.
A = [10 -1 2 0; -1 -11 -1 3; 2 -1 10 -1; 0 3 -1 8];
b = [6;25;-11;15];
x = Jacobi(A, b, 0.001, 40);
xtrue = A \ b
x =
0.6729
-1.5936
-1.1612
2.3275
xtrue =
0.6729
-1.5936
-1.1612
2.3274
这就是我对高斯-赛德尔的看法:
This is what I get with Gauss-Seidel:
x =
0.6729
-1.5936
-1.1612
2.3274
这肯定会融合在一起，并且系统在对角线占主导地位.
This certainly converged for both, and the system is diagonally dominant.
因此，您的代码没有任何问题.您只是在指定一个系统，无法使用Jacobi解决.最好将Gauss-Seidel用于围绕这种解决方案的迭代方法.原因是因为您将立即使用当前迭代中的信息，并将其扩展到其他变量. Jacobi不这样做，这就是它更快分散的原因.对于Jacobi，您可以看到示例1未能收敛，而示例2却收敛了.高斯-塞德尔同时收敛.实际上，当它们都收敛时，它们就非常接近于真正的解决方案.
As such, there is nothing wrong with your code. You are just specifying a system that can't be solved using Jacobi. It's better to use Gauss-Seidel for iterative methods that revolve around this kind of solving. The reason why is because you are immediately using information from the current iteration and spreading this to the rest of the variables. Jacobi does not do this, which is the reason why it diverges more quickly. For Jacobi, you can see that Example #1 failed to converge, while Example #2 did. Gauss-Seidel converged for both. In fact, when they both converge, they're quite close to the true solution.
同样，您需要确保您的系统在对角线占主导地位，以便确保您具有收敛性.不执行此规则...好吧...您可能会冒险，因为它可能会收敛，也可能不会收敛.
Again, you need to make sure that your systems are diagonally dominant so you are guaranteed to have convergence. Not enforcing this rule... well... you'll be taking a risk as it may or may not converge.
祝你好运！
这篇关于Jacobi迭代不会结束的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
