Open to feedback of working code (please review)

Question

Lavorizia Vaughn 2021-12-2

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1600920-open-to-feedback-of-working-code-please-review

评论： Lavorizia Vaughn 2021-12-2

Hey whats up folks, I was you could give my code a look .I have implements both the backward Euler method and Newtons method with f=f(t,y), dfdy=f'(t,y), maxiter=maximum number of iterations, and N=the number of steps. My code is below. I am open to feedback and possible changes. Any help would surely be appreciated.

my code:

function [t,w] = backeuler_four(f, dfdy, a, b, alpha, N, maxiter, tol)
h = (b-a)/N;
t = a:h:b;
w = t*0;
w(1) = alpha;
for i = 1:N
    w0=w(i);
    wj=w0;
    for j=1:maxiter
        wj=wj-(wj - w0 - h*f(t(i+1),wj)) / (1 - h*dfdy(t(i+1),wj));
        error=(wj - w0 - h*f(t(i+1),wj)) / (1 - h*dfdy(t(i+1),wj));
        fprintf('%d %g\n', j, abs(error));
        if abs(error)<=tol, break;end
    end
end
fprintf('\n');
if abs(error) > tol, error('No Newton convergence.'); end
w(i+1)=wj;

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Walter Roberson 2021-12-2

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1600920-open-to-feedback-of-working-code-please-review#answer_845405

在 MATLAB Online 中打开

if abs(error)<=tol, break;end

That only breaks out of one level of for loop.

Suppose you get convergence at i = 2, j = 7. Then you leave the for j loop. But you are still inside the for i loop, and you overwrite error and so on. So your final test is really testing whether you got convergence when i == N.

3 个评论
显示 1更早的评论隐藏 1更早的评论

Walter Roberson 2021-12-2

编辑：Walter Roberson 2021-12-2

在 MATLAB Online 中打开

If you need a solution for each i value, then do not break out of the inner loop -- but in such a case, you should possibly store the error for each i value.

If you only need to execute until you find one solution, then you should break out of the outer loop.

function [t,w] = backeuler_four(f, dfdy, a, b, alpha, N, maxiter, tol)
h = (b-a)/N;
t = a:h:b;
w = t*0;
w(1) = alpha;
converged = false;
for i = 1:N
    w0=w(i);
    wj=w0;
    for j=1:maxiter
        wj=wj-(wj - w0 - h*f(t(i+1),wj)) / (1 - h*dfdy(t(i+1),wj));
        error=(wj - w0 - h*f(t(i+1),wj)) / (1 - h*dfdy(t(i+1),wj));
        fprintf('%d %g\n', j, abs(error));
        if abs(error)<=tol
            converged = true;
            break;
        end
    end
    if converged; break; end
end
fprintf('\n');
if ~converged
    error('No Newton convergence.');
end
w(i+1)=wj;   %not sure what this is about