MATLAB Answers

Newton's method iterations

5 views (last 30 days)
James Holden
James Holden on 1 Oct 2021
Commented: Alan Stevens on 4 Oct 2021
I have a code written for Netwon's method below. However, not all initial starting guesses would converge. How would I find out which initial values for x would converge for the function in, say 20 or 30 or 50 iterations? For example, I'll start with x0 = 0. Does this converge? Then I'll try x0 = 0.1. Does this converge? The x0 = 0.2, then 0.3, etc. Which ones would converge?
I am trying to figure out how to summarize this information in a table. Below is the code that I have so far for Newton's Method.
%Newton's Method x0 = 7
f = @(x) 2*exp(-2*x) + 4*sin(x) - 2*cos(2*x);
fp = @(x) 4*(-exp(-2*x) + sin(2*x) + cos(x));
x0 = 7;
N = 10;
tol = 1E-6;
x(1) = x0;
n = 2;
nfinal = N + 1;
while (n <= N + 1)
fe = f(x(n - 1));
fpe = fp(x(n - 1));
x(n) = x(n - 1) - fe/fpe;
if (abs(fe) <= tol)
nfinal = n;
break;
end
n = n + 1;
end
plot(0:nfinal - 1,x(1:nfinal),'o-')
title('Solution:')
xlabel('Iterations')
ylabel('X')
table([1:length(x)]',x')
x(n) = x(n - 1) - fe/fpe;
fprintf('%3d: %20g %20g\n', n, x(n), abs(fe));
if (abs(fe) <= tol)

Answers (2)

Alan Stevens
Alan Stevens on 1 Oct 2021
I suggest you plot a graph of your function, then you can see where good initial estimates would be. For example
f = @(x) 2*exp(-2*x) + 4*sin(x) - 2*cos(2*x);
x = 0:0.1:20;
y = f(x);
plot(x,y),grid
  6 Comments
Alan Stevens
Alan Stevens on 4 Oct 2021
Try this
f = @(x) 2*exp(-2*x) + 4*sin(x) - 2*cos(2*x);
fp = @(x) 4*(-exp(-2*x) + sin(2*x) + cos(x));
x0 = 0.1:0.1:8;
cv = 0; fv = 0;
N = 10;
tol = 1e-6;
for j = 1:numel(x0)
x = x0(j);
n = 0;
fe = 1;
while (abs(fe) > tol) && (n <= N)
n = n+1;
fe = f(x);
fpe = fp(x);
x = x - fe/fpe;
end
if abs(fe)<=tol
cv=cv+1;
feval1(cv,1) = fe;
C(cv,:) = [x0(j), n];
else
fv = fv+1;
feval2(fv,1) = fe;
F(fv,1) = x0(j);
end
end
format shorte
if cv>0
disp('Converged')
disp('[x0, n, fe]')
disp([C, feval1])
end
disp(' ')
if fv>0
disp(['Failed to converge in ' num2str(N) ' iterations'])
disp('x0, fe')
disp([F, feval2])
end

Sign in to comment.


Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!