Open to feedback of working code (please review)
2 次查看(过去 30 天)
显示 更早的评论
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 个评论
采纳的回答
Walter Roberson
2021-12-2
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 个评论
Walter Roberson
2021-12-2
编辑:Walter Roberson
2021-12-2
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
更多回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Loops and Conditional Statements 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!