最速下降求解具有希尔伯特矩阵的线性系统
Steepest descent to find the solution to a linear system with a Hilbert matrix
我正在使用 最速下降法 来计算具有 5x5 希尔伯特矩阵的线性系统的解。我相信代码很好,因为它给了我正确的答案。
我的问题是:
我认为迭代太多次才能得到正确的答案。我相信我可能错过了算法中的某些内容,但我现在不确定是什么。
我不确定这是否是实现该算法的最有效方法,此外,"tol"选择哪个有点令人困惑。
任何关于这些的见解都将不胜感激(尤其是 1.)。谢谢!
% Method of Steepest Descent with tol 10^-6
h = hilb(5); %Hilbert 5x5 matrix
b = [1;1;1;1;1]; %solution matrix
solution = zeros(d,1); %Initialization
residual = h*solution - b;
tol = 10^(-6)
count = 0;
while residual'*residual > tol;
roe = (residual'*residual)/(residual'*h*residual);
solution = solution - roe*residual;
residual = h*solution - b;
count = count + 1;
end
count
solution
%Method of Steepest Descent with tol 10^-12
solution = zeros(d,1);
residual = h*solution - b;
tol = 10^(-12)
count = 0;
while residual'*residual > tol;
roe = (residual'*residual)/(residual'*h*residual);
solution = solution - roe*residual;
residual = residual - roe*h*residual;
count = count + 1;
end
count
solution
%another_solution = invhilb(5)*b %Check for solution
我不具备从数学方面处理你的问题的知识,但从编程的角度来看,有一点我可以指出。
确实你是对的。在得到结果之前需要太多次迭代:
Elapsed time is 6.487824 seconds.
count =
292945
当您定义步长以接近最终结果时,必须优化步长。否则你要么接近答案的速度太慢,要么因为你的步长太大而多次绕过最佳答案并绕着它走。
为了优化步长,我先根据你的脚本形成一个函数(加上一些小的修改):
function [solution, count, T] = SDhilb(d, step)
h = hilb(d);
tic
b = ones(d, 1);
solution = zeros(d, 1);
residual = h * solution - b;
res2 = residual' * residual;
tol = 10^(-6);
count = 0;
while res2 > tol;
roe = res2/(residual' * h * residual);
solution = solution - step * roe * residual; % here the step size appears
residual = h * solution - b;
count = count + 1;
res2 = residual' * residual; % let's calculate this once per iteration
end
T = toc;
现在将此函数用于一系列步长 (0.1:0.1:1.3) 和几个希尔伯特矩阵 (d = 2, 3, 4, 5) 很明显 1
不是一个有效步长:
相反,step = 0.9
似乎更有效率。让我们看看在 hilb(5)
:
情况下它的效率如何
[result, count, T] = SDhilb(5, 0.9)
result =
3.1894
-85.7689
481.4906
-894.8742
519.5830
count =
1633
T =
0.0204
快两个数量级以上。
以类似的方式,您可以尝试 tol
的不同值,看看它对速度的影响有多大。在那种情况下没有最优值:公差越小,你需要等待的时间越长。
看起来你正确地实现了算法(最速 descent/gradient 下降与精确线搜索以最小化凸二次函数)。
收敛很慢,因为问题是病态的:您考虑的希尔伯特矩阵的条件数大于 400000。已知梯度下降在问题病态时很慢-有条件的。
改为考虑条件良好的问题,例如通过将恒等式添加到希尔伯特矩阵 (h = hilb(5)+eye(5)),相同的代码仅在 7 次迭代(和条件数)后终止因为该矩阵小于 3).
我正在使用 最速下降法 来计算具有 5x5 希尔伯特矩阵的线性系统的解。我相信代码很好,因为它给了我正确的答案。
我的问题是:
我认为迭代太多次才能得到正确的答案。我相信我可能错过了算法中的某些内容,但我现在不确定是什么。
我不确定这是否是实现该算法的最有效方法,此外,"tol"选择哪个有点令人困惑。
任何关于这些的见解都将不胜感激(尤其是 1.)。谢谢!
% Method of Steepest Descent with tol 10^-6
h = hilb(5); %Hilbert 5x5 matrix
b = [1;1;1;1;1]; %solution matrix
solution = zeros(d,1); %Initialization
residual = h*solution - b;
tol = 10^(-6)
count = 0;
while residual'*residual > tol;
roe = (residual'*residual)/(residual'*h*residual);
solution = solution - roe*residual;
residual = h*solution - b;
count = count + 1;
end
count
solution
%Method of Steepest Descent with tol 10^-12
solution = zeros(d,1);
residual = h*solution - b;
tol = 10^(-12)
count = 0;
while residual'*residual > tol;
roe = (residual'*residual)/(residual'*h*residual);
solution = solution - roe*residual;
residual = residual - roe*h*residual;
count = count + 1;
end
count
solution
%another_solution = invhilb(5)*b %Check for solution
我不具备从数学方面处理你的问题的知识,但从编程的角度来看,有一点我可以指出。
确实你是对的。在得到结果之前需要太多次迭代:
Elapsed time is 6.487824 seconds.
count =
292945
当您定义步长以接近最终结果时,必须优化步长。否则你要么接近答案的速度太慢,要么因为你的步长太大而多次绕过最佳答案并绕着它走。
为了优化步长,我先根据你的脚本形成一个函数(加上一些小的修改):
function [solution, count, T] = SDhilb(d, step)
h = hilb(d);
tic
b = ones(d, 1);
solution = zeros(d, 1);
residual = h * solution - b;
res2 = residual' * residual;
tol = 10^(-6);
count = 0;
while res2 > tol;
roe = res2/(residual' * h * residual);
solution = solution - step * roe * residual; % here the step size appears
residual = h * solution - b;
count = count + 1;
res2 = residual' * residual; % let's calculate this once per iteration
end
T = toc;
现在将此函数用于一系列步长 (0.1:0.1:1.3) 和几个希尔伯特矩阵 (d = 2, 3, 4, 5) 很明显 1
不是一个有效步长:
相反,step = 0.9
似乎更有效率。让我们看看在 hilb(5)
:
[result, count, T] = SDhilb(5, 0.9)
result =
3.1894
-85.7689
481.4906
-894.8742
519.5830
count =
1633
T =
0.0204
快两个数量级以上。
以类似的方式,您可以尝试 tol
的不同值,看看它对速度的影响有多大。在那种情况下没有最优值:公差越小,你需要等待的时间越长。
看起来你正确地实现了算法(最速 descent/gradient 下降与精确线搜索以最小化凸二次函数)。
收敛很慢,因为问题是病态的:您考虑的希尔伯特矩阵的条件数大于 400000。已知梯度下降在问题病态时很慢-有条件的。
改为考虑条件良好的问题,例如通过将恒等式添加到希尔伯特矩阵 (h = hilb(5)+eye(5)),相同的代码仅在 7 次迭代(和条件数)后终止因为该矩阵小于 3).