Numerical Difficulty in Evaluating a Function.

Question

I was asked, given the function $$f(x)=\ln{\left(x-\sqrt{x^{2}-1}\right)},\ x \gg 1,$$ to discuss the numerical difficulty in evaluating for the given values of x.

Would it be correct to say that, since $$\lim\limits_{x\to\infty}\left(x-\sqrt{x^{2}-1}\right)=\lim\limits_{x\to\infty}\frac{1}{x+\sqrt{x^{2}-1}}=0,$$

and due to the finite accuracy of floating point arithmetic, that there will eventually be a loss of significance that'll make the subtraction zero, and thus make the function undefined?

Is there anything else that could be said?

Answer 1

If you consider $$A=x-\sqrt{x^{2}-1}$$ there will be already problems because of what you properly wrote. This is why it is better to write $$A=x-\sqrt{x^{2}-1}\times \frac{x+\sqrt{x^{2}+1} } {x+\sqrt{x^{2}-1} }=\frac 1{ x+\sqrt{x^{2}-1}} \sim\frac{1}{2 x}$$ which will not make any problem; the same for $\log(A)$.

Answer 2

The error of evaluating $x-\sqrt{x^2-1}$ is about the same as the error of evaluating $|x|+\sqrt{x^2+1}$ which is somewhere around $(1+2|x|)\mu$. While the exact value in the second case has the same growth for large $x$, in the original expression it goes to zero like $\frac1{2x}$, so that the relative error will grow unbounded like $\sim4x^2\mu$.

At the point where $x^2\mu>1$, that is, $x>10^8$, the subtraction of $1$ under the root has no longer any effect and the value of the original expression is zero, giving a relative error of $1$.

Answer 3

I can be useful to write down the relative error in the computed result in terms of the condition number (independent of the algorithm/ expression used to compute $f$) and of the truncation errors committed in each step of the algorithm (let's say $k$ steps). You can obtain, as a first order approximation an expression of the type $$ \varepsilon_f \approx Q_0(x) \varepsilon_x + \sum_{i=1}^k Q_k(x) e_{a_k} $$

The numerical stability is only guaranteed if all coefficients $Q_0, \cdots, Q_k$ are bounded.

In this case, if you assume that $f$ is computed using the algorithm $$ z_1 = x^2-1, \quad z_2 = \sqrt{z_1}, \quad z_3 = x-z_2, \quad z_4= \log z_3 $$

You get that

\begin{align*} \varepsilon_{z_4} = & -\dfrac{x}{\sqrt{x^2-1} \cdot \log(x-\sqrt{x^2-1})} \varepsilon_x -\dfrac{\sqrt{x^2-1}}{2 \log(x -\sqrt{x^2-1})\cdot (x-\sqrt{x^2-1})} e_{a_1}\\ & -\dfrac{\sqrt{x^2-1}}{\log(x -\sqrt{x^2-1})\cdot (x-\sqrt{x^2-1})} e_{a_2}+\dfrac{1}{\log(x-\sqrt {x^2-1})} e_{a_3} + e_{a_4} \end{align*}

Looking at these coefficients, you'll see that $Q_2$ and $Q_3$ are unbounded ($x\to +\infty$) and so the algorithm is unstable.