Help with badly working loop
7 views (last 30 days)
Show older comments
Hi, hello everyone, I've been working on a method that implements the Backwartd Euler method with the Newton method, amd I have the following code:
function [t,w] = backeuler(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
Now, when I terminate:
f = @(t,y) y^2 * (1-y);
dfdy = @(t,y) 2*y - 3*y^2;
a = 0; b = 2000; alpha = 0.9; maxiter = 20; tol = 1e-12;
[t3,w3] = backeuler(f, dfdy, a, b, alpha, 1, maxiter, tol);
[t4,w4] = backeuler(f, dfdy, a, b, alpha, 5, maxiter, tol);
[t5,w5] = backeuler(f, dfdy, a, b, alpha, 10, maxiter, tol);
I've been getting:
%1st solution
1 0.0270214
2 0.00149354
3 4.46627e-06
4 3.98803e-11
5 7.88716e-18
%2nd solution
1 0.0268375
2 0.00147092
3 4.32538e-06
4 3.7348e-11
5 4.40305e-17
%3rd soltuion
1 0.0266097
2 0.00144316
3 4.15576e-06
4 3.44115e-11
5 1.74336e-17
Which are the three solutions for the givens. But I should be getting more solutions, for example, I only have one solution for [t4,w4] = backeuler(f, dfdy, a, b, alpha, 5, maxiter, tol), but I should have four.
%1st solution for [t3,w3] = backeuler(f, dfdy, a, b, alpha, 1, maxiter, tol);
1 0.128469
2 0.0270214
3 0.00149354
4 4.46627e-006
5 3.98803e-011
6 7.88716e-018
%2nd(1) solution for [t4,w4] = backeuler(f, dfdy, a, b, alpha, 5, maxiter, tol);
1 0.128063
2 0.0268375
3 0.00147092
4 4.32538e-006
5 3.7348e-011
6 4.40305e-017
%2nd(2) soution for [t4,w4] = backeuler(f, dfdy, a, b, alpha, 5, maxiter, tol);
1 0.000249002
2 1.23664e-007
3 3.05022e-014
%2nd(3) solution for [t4,w4] = backeuler(f, dfdy, a, b, alpha, 5, maxiter, tol);
1 6.20646e-007
2 7.68531e-013
%2nd(4) solution for [t4,w4] = backeuler(f, dfdy, a, b, alpha, 5, maxiter, tol);
1 1.54774e-009
2 1.16283e-017
%2nd soltuoin(5)
1 3.8597e-012
2 9.69023e-018
Any help would be appreciated.
14 Comments
Voss
on 4 Dec 2021
I think you need to use the value wj that has been converged on as the initial wj for the next i-iteration. So, right before or after setting converged = true, add this line: w(i+1) = wj;
Answers (0)
See Also
Categories
Find more on Loops and Conditional Statements in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!