2.3 Lehmer's Algorithm

An alternate approach to speeding up Euclid's algorithm is due to Lehmer. One notices that when a and b have the same size, the integer part w of the quotient a/b is often single digit. Secondly, the process underlying Euclid's algorithm is the application of successive linear transformations (we will see this also in the Extended GCD computations later)

$$\begin{pmatrix} u v \end{pmatrix}$$ $$\mapsto$$ $$\begin{pmatrix} A & B C & D \end{pmatrix}$$ ^. $$\begin{pmatrix} u v \end{pmatrix}$$ = $$\begin{pmatrix} A u + B v C u + Dv \end{pmatrix}$$

Thus, we can repeatedly try to find the w as long as it is ``small'' and keep track of the operations involved by means of a matrix. When we cannot proceed further, we apply this matrix to the original data, and then try again. Occasionally, we many have to break out of a ``deadlock'' by performing a long division.

One way to check that w is small is to compare it with the quotient a₀/b₀ of the leading digits of a and b, when a and b have the same number of digits. More specifically if the integer part of a₀/(b₀ + 1) and the integer part of (a₀ + 1)/b₀ are equal, then they are equal to w. The detailed algorithm is given below (we assume a $\geq$ b).

If b is a single digit then we apply Euclid's algorithm to get the answer. Otherwise, we set x to be the leading digit (in base M) of a and y to the the leading digit of b at the same place (i. e. if a = (a₀...a_p) then b = (b₀...b_p) for the same p). We compute the invertible matrix $\left(\vphantom{ \begin{smallmatrix}A & B C & D \end{smallmatrix} }\right.$$\begin{smallmatrix}A & B C & D \end{smallmatrix}$$\left.\vphantom{ \begin{smallmatrix}A & B C & D \end{smallmatrix} }\right)$, by an iteration of the following steps after initialising it as the identity matrix.

We compute the integer quotient w₁ of (x + A)/(y + C) and the integer quotient w₂ of (x + B)/(y + D). We set w = w₁ if w₁ = w₂ and perform the matrix multiplication

$$\begin{pmatrix} 0 & 1 1 & -w \end{pmatrix}$$ ^. $$\begin{pmatrix} A & B C & D \end{pmatrix}$$ = $$\begin{pmatrix} C & D A-wC & B-wD \end{pmatrix}$$

We then replace our matrix by the resulting matrix. Similarly, we replace x by y and y by x - wy. Then we iterate.

Lemma 4   In the above situation at most one of y + C and y + D is 0. Moreover, we have the inequalities

$$\leq \leq \leq$$ 0 0

$$\leq \leq \leq$$ 0 0

This ensures that all operations involved in this part of the computation are digits. This procedure will exit when w₁ $\neq$ w₂ or when one of y + C or y + D is 0.

If the matrix computed has B = 0 then we need to perform a long division operation using a and b (this case occurs only when the first w₁ and w₂ and do not match). Otherwise, we replace a by Aa + Bb and b by Ca + Db and go back to the start.

The advantage of this procedure is that long division is only performed when absolutely essential; i. e. when the quotient w is larger than M. One can show that this case occurs sufficiently infrequently to account for the overhead of computing the matrix $\left(\vphantom{ \begin{smallmatrix}A & B C & D \end{smallmatrix} }\right.$$\begin{smallmatrix}A & B C & D \end{smallmatrix}$$\left.\vphantom{ \begin{smallmatrix}A & B C & D \end{smallmatrix} }\right)$.

Proof. (of the lemma) We only need to note that at any stage x and y are the result of applying Euclid's algorithm to the original pair x and y. Similarly, x + A and y + C at the same stage are the result of applying Euclid's algorithm to the pair x + 1 and y; and x + B and y + D are associated to the pair x and y + 1. The lemma follows from the consequence (of Euclid's algorithm) that these are all decreasing and initially lie between 0 and M. $\qedsymbol$