# How many Taylor series terms are needed to accurately approximate $\sqrt{a+x}-\sqrt{a}$?

Naive evaluation of $$\sqrt{a + x} - \sqrt{a}$$ when $$|a| >> |x|$$ suffers from catastrophic cancellation and loss of significance.

WolframAlpha gives the Taylor series for $$\sqrt{a+x}-\sqrt{a}$$ as: $$\frac{x}{2 \sqrt{a}} - \frac{x^2}{8 a^{3/2}} + \frac{x^3}{16 a^{5/2}} - \frac{5 x^4}{128 a^{7/2}} + \frac{7 x^5}{256 a^{9/2}} + O(x^6)$$ which (I think) equals: $$\sqrt{a} \left( \frac{1}{2} \left(\frac{x}{a}\right) - \frac{1}{8} \left(\frac{x}{a}\right)^2 + \frac{1}{16} \left(\frac{x}{a}\right)^3 - \frac{5}{128} \left(\frac{x}{a}\right)^4 + \frac{7}{256} \left(\frac{x}{a}\right)^5 + O\left(\left(\frac{x}{a}\right)^6\right) \right)$$

How quickly do the coeffients decrease?

How many terms are needed to reach $$53$$ bits of accuracy (IEEE double precision) in the result given that $$10^{-300} < \left|\frac{x}{a}\right| < 1$$ is known?

Alternatively, what are the threshold values of $$\left|\frac{x}{a}\right|$$ where the number of terms changes?

What about rounding errors, assuming each value is stored in double precision?

• @Winther shouldn't that be $x / \left(\sqrt{a + x} + \sqrt{a}\right)$? thanks, investigating... – Claude May 13 at 19:00

To avoid cancellation error the first thing to do is to write:

$$\sqrt{a+x}-\sqrt{a}=\frac{x}{\sqrt{a+x}+\sqrt{a}}=\sqrt{a}\frac{x}{a} \frac{1}{1+\sqrt{1+\frac{x}{a}}}$$

then with $$y=\frac{x}{a}$$ you must approximate this $$\sqrt{a}\frac{y}{1+\sqrt{1+y}}$$ fonction for $$y\in[10^{-300},1]$$. This function has nothing pathological and IMHO can be computed in a straightforward way.

If you really want to use Taylor series for $$y\sim 0$$ $$\sqrt{a}\frac{y}{1+\sqrt{1+y}}=\sqrt{a}(\frac{y}{2}-\frac{y^2}{8}+\frac{y^3}{16}-\frac{5 y^4}{128}+\frac{7 y^5}{256}+O\left(y^6\right))$$ I assume that the series is alternating, hence the error term $$e$$ is majored by $$|e|<\sqrt{a}\frac{7y^5}{256}$$. For instance if you want $$|e|<10^{-q}$$ you can use the Taylors series for $$0\le y \le y_*$$ where $$y_*$$ is such that $$\sqrt{a}\frac{7y_*^5}{256}<10^{-q}$$ which gives $$y_*<10^{-q/5}(\frac{256}{7\sqrt{a}})^{1/5}$$

Example: with $$q=5$$, $$a=3$$

We get $$y_*<0.184042$$.

That means that you can use the Taylors series for $$y_*=\frac{x_*}{a}<0.184042$$, hence $$x_*<3\times 0.184042 \approx 0.552125$$.

Let's try with $$x=0.55$$.

With the initial formula we find: $$\sqrt{a+x}-\sqrt{a}\approx 0.152094$$

With the Taylor series, with $$y=\frac{0.55}{3}$$ we get $$\sqrt{a}(\frac{y}{2}-\frac{y^2}{8}+\frac{y^3}{16}-\frac{5 y^4}{128}+\frac{7 y^5}{256})\approx 0.152095$$

We see that the error $$|e|=|0.152094-0.152095|\approx 1.17957\times 10^{-6}$$ is less that $$10^{-q}=10^{-5}$$ as expected

The Taylor series is

$$\sqrt{a+x} - \sqrt{a} = \sum_{k=1}^\infty (-1)^{k+1} \frac{(2k)!}{(k!)^2(2k-1)} 4^{-k} a^{1/2-k} x^k$$ If $$|x/a| < 1$$, the absolute values of the terms decrease, since if $$c_k = (2k)!/((k!)^2 (2k-1) 4^k)$$, $$\frac{c_{k+1}}{c_k} = \frac{2k-1}{2k+2} < 1$$ Thus if $$a > x > 0$$ the absolute value of the error is always less than that of the next term. However, if $$x/a$$ is close to $$1$$ the convergence is rather slow: $$c_k \sim \frac{1}{2 \sqrt{\pi} k^{3/2}}$$ so that won't be less than $$2^{-53}$$ unless $$k > 1.862 \times 10^{10}$$ approximately.

It was pointed out by Robert Israel that the series does badly when $$|x| \approx |a|$$, but in that case the loss of significance of the naive evaluation is small.

It was also suggested by Winther (and a since-deleted answer) to rewrite as $$\frac{x}{\sqrt{a+x}+\sqrt{a}}$$ The series for the denominator is similar to the series in the question, only with a leading constant term. This means that when $$\left|\frac{x}{a}\right|$$ is small enough, the addition of terms eventually becomes insignificant in double arithmetic.

If $$\left|\frac{x}{a}\right| < 2^{-52}$$, $$1$$ term is sufficient. Otherwise

If $$\left|\frac{x}{a}\right| < 2^{-25}$$, $$2$$ terms are sufficient. Otherwise

If $$\left|\frac{x}{a}\right| < 2^{-16}$$, $$3$$ terms are sufficient. Otherwise

If $$\left|\frac{x}{a}\right| < 2^{-11}$$, $$4$$ terms are sufficient. Otherwise

If $$\left|\frac{x}{a}\right| < 2^{-9}$$, $$5$$ terms are sufficient. Otherwise

If $$\left|\frac{x}{a}\right| > 2^{-9}$$, the loss of precision in the addition $$a + x$$ is relatively small.

But in fact, perhaps $$\frac{x}{\sqrt{a + x} + \sqrt{a}}$$ evaluated in double precision is good enough for all $$|x| << |a|$$ and series are unnecessary?