Selection of the numerical method

Question

Given the transcendental equation: $$ \frac{\tan x}{x} + c = 0, $$ where $c$ is any real number.

I tried Newton's method, but it is very bad fit. Which numerical method will be the smartest solution in this case?

Answer 1

If you can find an interval which contains a single root with a sign change, and no singularities, then the secant method or a variant on it should work. This method avoids the problems associated with picking a bad starting point.

Answer 2

Avoid the singularities by solving

$$\sin x+cx\cos x=0$$ and try starting values $\dfrac\pi2+k\pi$ (near the vertical asymptotes of the original function).

Answer 3

There will be 2 solutions in $(-\frac {\pi}{2}, \frac {\pi}{2})$ if $c> 1$

and none if $c< 1$

For each interval $(\frac {k}{2} \pi , (\frac {k}{2} + 1) \pi)$ there will be a solution.

I would choose a point an interval. Choose a point in that interval, and apply Newtons method.

If the suggested next point is inside the targeted interval, use it. If it is outside of the interval, use the midpoint of $x_n$ and the boundary point of the interval as $x_{n+1}$

Answer 4

Usually it is better to have a compressing function like $\arctan$ instead of an expanding one like $\tan$. Thus solve $$x+\arctan(cx)=kπ$$ which has a contracting fixed-point form $$x=g(x)=kπ−\arctan(cx)$$ for any $k∈\Bbb Z$ and $|c|<1$.

For any $k\ne 0$ this iteration is contractive after the first step.

• After the first step the iteration is restricted to $[(k-0.5)\pi,(k+0.5)\pi]$.
• On this interval $$|g'(x)|=\frac{|c|}{1+(cx)^2}\le\frac1{2|x|}<1$$ which is true for all $|x|>1$, thus especially for $|x|>(|k|-0.5)\pi\ge\frac\pi2$.

Use Aitken's delta-square process or directly Newton for faster convergence.

def aitken(g,x0,eps):
    x0 = g(x0)
    while True:
        x1 = g(x0);
        if abs(x1-x0)<eps: break
        x2 = g(x1)
        x0 = x0 - (x1-x0)**2/(x2-2*x1+x0)
    return x0

for k in range(1,n+1):
    print k, aitken(lambda x: k*pi-atan(c*x), k*pi, eps);  

As a graph one might visualize this plotting $-\arctan(cx)$ against $x-kπ$. This shows that one will find roots close to $(k-0.5\text{sign}(c))\pi$ for $k\ne 0$.

For $c<-1$ this has two non-zero solutions for $k=0$ close to $\pm\frac\pi2$. Using $\tan x\approx x+\frac13x^3$ finds these roots close to $\pm\min(\frac\pi2,\sqrt{3(|c|-1)}$. Using this as initial value, the iteration of $g$ still converges, but a strict proof might be more involved.

if c<-1: print 0, aitken(lambda x: -atan(c*x), min(pi/2, (-3*(c+1))**0.5), eps)

Answer 5

Using Yves Daoust's answer, you will find a root in each interval $(k \pi,(k+1)\pi)$ so starting Newton method with $x_0=\frac \pi 2+k \pi$ will make the method converging quite well.

We can improve the initial guess using the simplest $[1,1]$ Padé approximant of $$f(x)=\sin(x)+c x \cos(x)$$ and get $$\color{blue}{x_0=a+\frac{2 a c}{2 c \left(a^2 c+1\right)+1}}\qquad \text{where}\qquad a=\frac \pi 2+k \pi$$

For illustration purposes, let us use $c=1.234$. The table below reports the estimate and the solution for a few values of $k$ between $0$ and $10$. $$\left( \begin{array}{ccc} 0 & 1.923787490 & 1.962413795 \\ 1 & 4.875967347 & 4.877045230 \\ 2 & 7.955291292 & 7.955494773 \\ 3 & 11.06858652 & 11.06865718 \\ 4 & 14.19416434 & 14.19419687 \\ 5 & 17.32548133 & 17.32549895 \\ 6 & 20.45992874 & 20.45993934 \\ 7 & 23.59626779 & 23.59627466 \\ 8 & 26.73383619 & 26.73384090 \\ 9 & 29.87224814 & 29.87225150 \\ 10 & 33.01126382 & 33.01126630 \end{array} \right)$$ Repeating the calculation with $c=12.345$ $$\left( \begin{array}{ccc} k & \text{estimate} & \text{solution} \\ 0 & 1.620661976 & 1.620734854 \\ 1 & 4.729513663 & 4.729514741 \\ 2 & 7.864281368 & 7.864281569 \\ 3 & 11.00293616 & 11.00293623 \\ 4 & 14.14289442 & 14.14289445 \\ 5 & 17.28344637 & 17.28344638 \\ 6 & 20.42431830 & 20.42431831 \\ 7 & 23.56538232 & 23.56538232 \\ 8 & 26.70657067 & 26.70657067 \\ 9 & 29.84784411 & 29.84784412 \\ 10 & 32.98917834 & 32.98917834 \end{array} \right)$$

Repeating the calculation with $c=123.456$ $$\left( \begin{array}{ccc} k & \text{estimate} & \text{solution} \\ 0 & 1.575936039 & 1.575936117 \\ 1 & 4.714107235 & 4.714107235 \\ 2 & 7.855012829 & 7.855012829 \\ 3 & 10.99631090 & 10.99631090 \\ 4 & 14.13773988 & 14.13773988 \\ 5 & 17.27922837 & 17.27922837 \\ 6 & 20.42074891 & 20.42074891 \\ 7 & 23.56228867 & 23.56228867 \\ 8 & 26.70384088 & 26.70384088 \\ 9 & 29.84540161 & 29.84540161 \\ 10 & 32.98696842 & 32.98696842 \end{array} \right)$$

Since the successive derivatives are simple, you could also use Halley or Householder methods. Starting with the proposed guess, I bet that a couple of iterations will give the solution to quite high accuracy.

Edit

Trying to improve, using the $[2,1]$ Padé approximant of $$f(x)=\sin(x)+c x \cos(x)$$ we can get $$\color{blue}{x_0=a+\frac{6 c \left(a^2 c+1\right)+3}{a c \left(6 c \left(a^2 c+2\right)+5\right)}}\qquad \text{where}\qquad a=\frac \pi 2+k \pi$$

For the first root $(k=0)$ in the worked examples, we would get $1.972142041$, $1.620737162$ and $1.575936117$.