Matlab Newton-Raphson for x-tan(x)

Question

For the function $f(x)=x-tan(x)$ I am trying to find the 3 smallest roots to six decimal places.

Attempt:

We need to use an iterative method to solve this problem, and I choose to use Newton-Raphson for this equation which can also be written as $x=tan(x)$.

Clearly the roots of $tan(x)-x$ are given by the intersections of $y = x$ and $y = tan (x)$. So the first 3 roots are around $3 \pi/2$, $5\pi /2$, and $7 \pi/2$.

Here is my iteration scheme so far with $x_0 = 3 \pi/2$:

$$x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)} = x_n - \frac{x_n- tan(x_n)}{1- \frac{1}{cos^2(x_n)}}$$

Here's the code:

for i=1:n_itn
x(i)=x(i)-((x(i)-tan(x(i)))/(1-(1/cos(x(i)))^2))
end

But my code doesn't seem to work, I keep getting $4.7123889 (=x_0)$. There is also the following error:

Attempted to access x(2); index out of bounds because numel(x)=1.

So, what is wrong with my code, and what would be the correct iteration formula?

Any help is greatly appreciated.

P.S. I have tried $x(i+1)=...$ in the code, but all the outputs are still equal to the initial value of $x_0$.

Answer 1

Beside what has been said in comments, you selected as starting point $x_0=\frac{3\pi}2$ where the function is not defined.

Let us suppose that you start using $x_0=4.7$; the successive iterations will then be $$x_1=4.688331848$$ $$x_2=4.666984472$$ $$x_3=4.631183287$$ $$x_4=4.580473096$$ $$x_5=4.528429052$$ and so on for a solution which is close to $4.49341$. The process will converge but quite slowly.

The function is poorly conditionned because of the $\tan(x)$ term. Consider instead solving $$g(x)=x \cos(x)-\sin(x)$$ and do the same. The iterates will be $$x_1=4.499623657$$ $$x_2=4.493417956$$ $$x_3=4.493409458$$ which is the solution for ten significant figures.

Using this transform, the problem does not show any more discontinuities and solves much faster.

Answer 2

Try:

function [x]=MathEx1(x0,NN)
% x0 initial guess at root
% NN number of NR iterations to perform
% x array of length NN+1 with first element equal to x0
%    and the following elements the NN iterates

  x=zeros(1,NN+1);
  x(1)=x0;
  for i=2:NN+1
    x(i)=x(i-1)-(x(i-1)-tan(x(i-1))) / (1-(sec(x(i-1)))^2);
  end

endfunction

with x0=4.5this gives:

>> rr=MathEx1(4.5,10)
rr =

 Columns 1 through 9:

   4.5000   4.4936   4.4934   4.4934   4.4934   4.4934   4.4934   4.4934   4.4934

 Columns 10 and 11:

   4.4934   4.4934

>>
>> rr(11)-tan(rr(11))
ans =  -8.8818e-016

Note the positive roots near odd half multiple of $\pi$ are all slightly less than the odd half multiple and the iteration is ill-defined at such a point, so don't use an odd half multiple of $\pi$ as a starting value.