Самый крутой спуск, чтобы найти решение линейной системы с гильбертовой матрицей

Question

Я использую метод наискорейшего спуска, чтобы найти решение линейной системы с матрицей Гильберта 5x5. Я считаю, что код хорош в том смысле, что он дает мне правильный ответ.

Моя проблема заключается в том, что:

1. Я думаю, что это занимает слишком много итераций, чтобы получить к правильному ответу. Я считаю, что, возможно, что-то пропустил в алгоритме, но я не уверен, что на данный момент.

2. Я не уверен, что это самый эффективный способ реализации алгоритма, и, кроме того, он немного запутан, на что «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 

Answer 1

У меня нет знания, чтобы справиться с вашей проблемой из математической стороны, но с точки зрения программирования есть точка я мог бы отметить.

Действительно, вы правы. Это занимает слишком много итераций, пока не дойдет до результата:

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) и пара матриц Гильберта (д = 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, чтобы увидеть, насколько драматически это может повлиять на скорость. В этом случае нет оптимального значения: чем меньше допуск, тем больше времени вам нужно подождать.

Answer 2

Похоже, что вы правильно реализовали алгоритм (крутой спуск/градиентный спуск с точным линейным поиском для минимизации выпуклой квадратичной функции).

Конвергенция медленная, потому что проблема плохо обусловленная: матрица Гильберта, которую вы считаете, имеет состояние выше 400000. Градиентный спуск, как известно, медленный, когда проблема плохо обусловлена.

Рассмотрение проблемы с условной условностью, например, путем добавления идентичности к матрице Гильберта (h = hilb (5) + eye (5)), тот же код завершается только после 7 итераций (и номер условия для эта матрица меньше 3).