TAOCP 4.7 Exercise 14

Let U(z)=z+U_k z^k+\cdots,\qquad k\ge 2, and define the iterates of $U$ by

Section 4.7: Manipulation of Power Series

Exercise 14. ▶ [HM30] Fill in the details of Brent and Traub's method for calculating $U^{[n]}(z)$ when $U(z) = z + U_k z^k + \cdots$, using (27) and (28).

Verified: yes
Solve time: 6m57s


Setup

Let

$$ U(z)=z+U_k z^k+\cdots,\qquad k\ge 2, $$

and define the iterates of $U$ by

$$ U^{[0]}(z)=z,\qquad U^{[m+1]}(z)=U!\left(U^{[m]}(z)\right). $$

The problem is to make explicit Brent and Traub's procedure for computing the first $N$ coefficients of $U^{[n]}(z)$, using equations (27) and (28).

Equation (27) gives the semigroup law

$$ U^{[m+n]}=U^{[m]}\circ U^{[n]}, $$

and equation (28) is the special case

$$ U^{[2m]}=U^{[m]}\circ U^{[m]}. $$

The task is therefore to compute $U^{[n]}$ by repeated composition, exactly as ordinary exponentiation is computed by repeated squaring.

Solution

Write the binary expansion of $n$ as

$$ n=\sum_{j=0}^{r} \varepsilon_j 2^j, \qquad \varepsilon_j\in{0,1}, \qquad \varepsilon_r=1 . $$

Define power series

$$ A_0(z)=U(z), $$

and recursively

$$ A_{j+1}(z)=A_j(A_j(z)), \qquad 0\le j<r. $$

By (28),

$$ A_j(z)=U^{[2^j]}(z) $$

for every $j$. Indeed, $A_0=U^{[1]}$; if

$A_j=U^{[2^j]}$, then

$$ A_{j+1} =A_j\circ A_j =U^{[2^j]}\circ U^{[2^j]} =U^{[2^{j+1}]} $$

by (27).

Now set

$$ B(z)=z. $$

Scanning the binary digits of $n$, replace $B$ by

$$ B \leftarrow A_j(B) $$

whenever $\varepsilon_j=1$.

After all digits have been processed,

$$ B =U^{[\varepsilon_0]} \circ U^{[\varepsilon_1 2]} \circ \cdots \circ U^{[\varepsilon_r 2^r]} . $$

Since iterates of the same series commute,

$$ U^{[a]}\circ U^{[b]} =U^{[a+b]} $$

by repeated application of (27). Hence

$$ B =U^{\left[\sum_{j=0}^{r}\varepsilon_j2^j\right]} =U^{[n]}. $$

Therefore $B(z)$ is the desired iterate.

To obtain only the first $N$ coefficients, every composition is performed modulo $z^N$. Exercise 11 provides an algorithm for composing two power series to precision $N$, so each step above can be carried out with truncated series.

The complete procedure is therefore:

  1. Compute

$$ A_0=U,\qquad A_{j+1}=A_j\circ A_j \pmod{z^N} $$

for $0\le j<r$. 2. Initialize $B=z$. 3. For each $j$ with $\varepsilon_j=1$, replace

$$ B\leftarrow A_j(B)\pmod{z^N}. $$ 4. Output $B$.

The result is $U^{[n]}(z)$ through terms of degree $N-1$.

Verification

The construction of the series $A_j$ is correct because induction gives

$$ A_j=U^{[2^j]}. $$

The update step is correct because after processing digits

$\varepsilon_0,\ldots,\varepsilon_s$,

$$ B = U^{\left[\sum_{j=0}^{s}\varepsilon_j2^j\right]}. $$

This statement is proved by induction on $s$. When

$\varepsilon_s=0$, $B$ is unchanged. When $\varepsilon_s=1$,

$$ B \leftarrow U^{[2^s]}\circ U^{\left[\sum_{j=0}^{s-1}\varepsilon_j2^j\right]} = U^{\left[\sum_{j=0}^{s}\varepsilon_j2^j\right]} $$

by (27).

After the final digit is processed,

$$ B=U^{[n]}. $$

This completes the proof.

Notes

The method is the exact analogue of binary exponentiation. Equation (28) supplies the repeated-squaring step, while equation (27) supplies the multiplication step. The only additional ingredient is a routine for composing power series modulo $z^N$, such as the algorithm of Exercise 11.