# fixed point Iterative method for finding root of an equation

595 views (last 30 days)

Show older comments

The code below gives the root and the iteration at which it occur. The code goes into an infinite loop when the function contains any logarithmic or exponential function. But works fine with algebraic and trignometric function.

Can anyone point out the error?

Thank you in advance.

% let the equation whose root we want to find be x^3-5*x-7 = 0

% simplified eqation example:- f = @(x) (5x+7)^(1/3)

function [root,iteration] = fixedpoint(a,f) %input intial approiximation and simplified form of function

if nargin<1 % check no of input arguments and if input arguments is less than one then puts an error message

fprintf('Error! Atleast one input argument is required.');

return;

end

if nargin<2 % check no of input arguments and if input arguments is less than two then puts a=0

a=0;

end

x(1) = f(a) ;

i =1 ;

temp = false;

while temp == false

x(i+1) = f(x(i));

if x(i+1) == x(i)

temp = true;

root = x(i-1);

iteration = i-1;

return;

end

i = i +1 ;

end

x

end

##### 8 Comments

John D'Errico
on 30 Sep 2020

@Milind Amga - you misunderstood my point, which was to Alan's question.

You would rarely want to use fixed point iteration in the real world. (I have, but only on rare occasions.) Yes, you are solving this as requested because it was homework. It teaches you both about MATLAB as well as teaching you one method for solving a problem, even though that method is not really a very good one in practice, because it has some serious flaws as outlined. (In the end, I hope your teacher spends the time for you to understand why the method works, as well as the limitations on the method and why you probably want to use better methods in the end. It is the latter parts that often get skipped over.)

I made the statement about fzero because I see far too many scientists and engineers using crude methods they learned in school to solve a problem, not knowing that good code already exists to do the same thing. For example, long ago, I was asked to help a scientist to do a computation using determinants to solve a problem. They were solving 2x2 systems of equations, and they wanted to change it so they could solve 5x5 or 6x6 systems. But they did not know how to form the necessary higher order determinants.

My answer was no, I would not help them to do what they asked me to do, but that I would show them how to solve the problem using better methods.

### Accepted Answer

John D'Errico
on 30 Sep 2020

Edited: John D'Errico
on 30 Sep 2020

To understand fixed point iteration, we need to know why and when it will diverge. When does it converge? For some problems, the iterations might never be absolutely convergent, at least not without getting creative.

You comment about "better convergence". It is a question of divergence. Consider the trivial problem,

f = @(x) 1.1*x - 2.2;

Yes, we know the solution is x == 2. But does MATLAB?

function [xfinal,fval,ferr,itercount] = myfp(F,x0,ftol,maxit,verbosity)

% [xfinal,fval,itercount] = myfp(F,x0,ftol,maxit)

%

% myfp: solve for a root of the function F(x) == 0 using fixed point

% iteration. This will be achieved using the iteration as

% x = F(x) + x

%

% Arguments: (input)

% F - function of one variable

% x0 - initial value

% ftol - tolerance on the function value being zero

% maxit - maximum number of iterations allowed

% verbosity - (optional) flag to define output behavior, 0 = none, 1 = stuff on screen

%

% xfinal - final estimate of the root

% fval - function value at that point

% itercount - the number of iterations taken

if (nargin < 5) || isempty(verbosity)

% the default is to be comatose

verbosity = 1;

end

% initialize for iteration

itercount = 0;

xfinal = x0;

fval = F(x0);

ferr = abs(fval);

if verbosity

disp('Current x, Fval, Ferr')

end

% just a while loop that goes until it runs out of time, or it gets lucky

% and "converges"

while (ferr > ftol) && (itercount < maxit)

itercount = itercount + 1;

xprev = xfinal;

% the iteration itself

fval = F(xfinal);

xfinal = fval + xfinal;

% how far from 0 is the current iteration?

ferr = abs(fval);

if verbosity

disp([xfinal,fval,ferr])

end

end

That should work. As I set it up, I am solving the problem as

F(x) - x + x == 0

and therefore we have the iterations as

x = F(x) + x

Now, let me try it in MATLAB.

[xfinal,fval,ferr,itercount] = myfp(F,1.5,0.001,5,1);

Current x, Fval, Ferr

0.95 -0.55 0.55

-0.205 -1.155 1.155

-2.6305 -2.4255 2.4255

-7.72405 -5.09355 5.09355

-18.420505 -10.696455 10.696455

So it clearly diverges. The longer I allow it to run, the further it will go, even if I start quite close to the solution (unless of course it is within the tolerance at the start point.)

So let me fix it by just dividing by 2.

F = @(x) (1.1*x - 2.2)/2;

[xfinal,fval,ferr,itercount] = myfp(F,1.5,0.001,5,1);

Current x, Fval, Ferr

1.225 -0.275 0.275

0.79875 -0.42625 0.42625

0.1380625 -0.6606875 0.6606875

-0.886003125 -1.024065625 1.024065625

-2.47330484375 -1.58730171875 1.58730171875

Hmm. Still diverges. Why?

Because the derivative of the right hand side as I created the fixed point iteration must have the property

abs(diff(F(x) + x)) < 1

in the vicinity of the root. If that derivative is larger than 1, divergence occurs. The scheme will wander away.

But mow let me change it to be:

F = @(x) (1.1*x - 2.2)/(-2);

>> [xfinal,fval,ferr,itercount] = myfp(F,1.5,0.001,5,1);

Current x, Fval, Ferr

1.775 0.275 0.275

1.89875 0.12375 0.12375

1.9544375 0.0556875000000001 0.0556875000000001

1.979496875 0.0250593749999999 0.0250593749999999

1.99077359375 0.01127671875 0.01127671875

[xfinal,fval,ferr,itercount] = myfp(F,1.5,0.001,100,0)

xfinal =

1.99962165967871

fval =

0.000462415948242256

ferr =

0.000462415948242256

itercount =

9

We can see it is clearly converging now, and in fact, required only 9 iterations.

Now, suppose we wish to find a root of the function f(x) == exp(x) - 2. My code does that by formulating the iteration as

x = exp(x) - 2 + x

And of course, since the derivative of the right hand side is just

exp(x) + 1

then fixed point iteratiion must always diverge. The starting value will not matter, unless it is EXACTLY at log(2). and even then, even the tiniest difference in the least significant bits will start to push it away from the root. The value of ftol would save you there though.

F = @(x) exp(x) - 2;

[xfinal,fval,ferr,itercount] = myfp(F,1.5,0.001,5,1);

Current x, Fval, Ferr

3.98168907033806 2.48168907033806 2.48168907033806

55.5891937168657 51.6075046465276 51.6075046465276

1.38701157272672e+24 1.38701157272672e+24 1.38701157272672e+24

Inf Inf Inf

Inf Inf Inf

So things got bad very fast. Modify the function simply as

F2 = @(x) -0.5*F(x)

[xfinal,fval,ferr,itercount] = myfp(F2,1.5,0.001,5,1);

Current x, Fval, Ferr

0.259155464830968 -1.24084453516903 1.24084453516903

0.611237841842661 0.352082377011693 0.352082377011693

0.68988235566251 0.0786445138198487 0.0786445138198487

0.693141856814415 0.00325950115190499 0.00325950115190499

0.693147180545774 5.32373135941899e-06 5.32373135941899e-06

As you can see, that final iteration will have converged.

But now, try to use the same trick to solve the related problem

F = @(x) exp(x) - 10;

>> F2 = @(x) -0.5*F(x);

>> [xfinal,fval,ferr,itercount] = myfp(F2,1.5,0.001,10,1);

Current x, Fval, Ferr

4.25915546483097 2.75915546483097 2.75915546483097

-26.1159481242028 -30.3751035890338 30.3751035890338

-21.1159481242051 4.99999999999773 4.99999999999773

-16.1159481245427 4.99999999966238 4.99999999966238

-11.1159481746502 4.99999994989251 4.99999994989251

-6.11595561126095 4.99999256338923 4.99999256338923

-1.11705929394996 4.998896317311 4.998896317311

3.71932035616684 4.8363796501168 4.8363796501168

-11.898858916735 -15.6181792729018 15.6181792729018

-6.89886231581379 4.99999660092118 4.99999660092118

And you should see it bounces around. It is easily fixed of course, as long as I choose a proper transformation of f.

F2 = @(x) -0.1*F(x);

[xfinal,fval,ferr,itercount] = myfp(F2,1.5,0.001,10,1);

Current x, Fval, Ferr

2.05183109296619 0.551831092966194 0.551831092966194

2.27361730438218 0.221786211415983 0.221786211415983

2.30216954873883 0.0285522443566585 0.0285522443566585

2.30258500666749 0.000415457928655094 0.000415457928655094

log(10)

ans =

2.30258509299405

As you see, after 4 iterations it was quite happy.

The problem with an exponential like that however, is exp(x) gets arbitrarily large, so any constant multiplier will fail for some variation of the function. We could get tricky though, if we can find some transformation of F that will always have a small derivative of the fixed point iterant. since the derivative of exp(x) will always be large for some values of x, we need to be tricky.

### More Answers (2)

Alan Stevens
on 30 Sep 2020

I defer to John on the guarantee. Certainly it's possible with x - exp(-x) = 0:

% x - exp(-x) = 0; rearrange as x = exp(-x)

x0 = 0.5; % Initial guess

f = @(x) exp(-x);

tol = 10^-10;

flag = true;

its = 0;

while flag == true && its<100

x = f(x0);

if abs(x-x0)<tol

flag = false;

end

x0 = x;

its = its+1;

end

disp([x, x-exp(-x), its])

Akash G M
on 19 Dec 2020

##### 0 Comments

### See Also

### Community Treasure Hunt

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

Start Hunting!