How to calculate error bounds for the function behind the following algorithm?

Question

D.E. Knuth, in his infamous The Art of Computer Programming, section 1.2.2, presents the following algorithm to efficiently calculate logarithms based on the method used by Henry Briggs:

Suppose that we have a binary computer and a number $x$, $1\leqslant x < 2$. Show that the following algorithm, which uses only shifting, addition, and subtraction operations proportional to the number of places of accuracy desired, may be used to calculate an approximation to $y=\log_bx$:

L1. [Initialize.] Set $y\leftarrow 0$, $z\leftarrow x$ shifted right $1$, $k\leftarrow 1$.

L2. [Test for end.] If $x=1$, stop.

L3. [Compare.] If $x-z < 1$, set $z\leftarrow z$ shifted right $1$, $k\leftarrow k+1$, and repeat this step.

L4. [Reduce values.] Set $x\leftarrow x-z, z\leftarrow x$ shifted right $k$, $y\leftarrow y+\log_b(2^k/(2^k-1))$, and go to L2.

My question arose from the exercise right after the aforementioned algorithm:

Determine the upper bounds on the accuracy of the algorithm in the previous exercise, based on the precision used in the arithmetic operations.

Knuth stays silent, giving no hints on its solution in his book.

I think I've seen something similar in the BKM Algorithm published in 1994. It uses a similar method to get logarithms and claims to have an error of $2^{-N}$ after $N$ iterations.

I kindly ask you to help me understand what does Knuth mean exactly by basing the determination of the upper bound on the precision used in the arithmetic operations, and how does one arrive at that error bound purely from this algorithm?

Answer 1

$y$ is a sum of the constants $\log_b\frac{2^k}{2^k-1}$, precomputed to a certain accuracy. These terms are summed until the value of $x$ is exhausted, i.e. (if I am right) the number of bits of $x$. So the total error is the error on the sum, itself bounded by $m$ times the accuracy, for $m$ terms.

Answer 2

Error of floating point operations is usually proportional to the size of the result.

If you add n terms starting with the largest, each intermediate result is about log n and the total rounding error eps * n * log n. If you start with the smallest terms, the intermediate results are all quite small and the total rounding error about eps * log n. N times smaller. (That’s useful whenever you use Taylor polynomials as well).