General post

Composition of formal power series

Originally posted 2024-11-22

Last updated 2024-11-22

Background

Note: This is supplemental material originally meant to go with the meetup on classical methods for power series coefficient extraction. I wasn’t able to find a simple, sufficient, self-contained treatment of computational coefficient extraction for power series compositions using iterative methods, so hope this fills that void.¹ Understanding integer compositions will make it easier to follow the “composition product terms” discussed below.²

Definition of power series composition

Let $h(z) = f(g(z))$ , where $f$ , $g$ , and $h$ are all formal power series. We want to compute the coefficients of $h(z)$ , $h_n = [z^n] h(z)$ for all $n$ . Unfortunately, doing so efficiently (in space and time) is non-trivial. We’ll develop a solution by working with the original series definition and inspecting the coefficients that come out at each step. Start by writing

h(z) = \sum_{m=0}^{\infty} f_m \left(\sum_{k=0}^{\infty} g_k z^k\right)^m

and expanding terms. We find that

\begin{aligned} h_0 &= [z^0] \sum_{m=0}^{\infty} f_m \left( \sum_{k=0}^{\infty} g_k z^k \right)^m \\ &= \sum_{m=0}^{\infty} f_m [z^0] \left( \sum_{k=0}^{\infty} g_k z^k \right)^m \\ &= f_0 + \sum_{m=1}^{\infty} f_m [z^0] \left( \sum_{k=0}^{\infty} g_k z^k \right)^m \\ &= f_0 + \sum_{m=1}^{\infty} f_m g_0^m \\ &= \sum_{m=0}^{\infty} f_m g_0^m \end{aligned}

Notice that this is a series that does not converge if $\vert g_0 \vert > 0$ (we require all coefficients to be finite; technically, to stabilize). Therefore, we conclude that we must have $g_0 = 0$ for the composition $h$ to converge. This immediately tells us that $h_0 = f_0$ , since the only term that contributes a constant is the one with $(g(x))^0$ .

Computing coefficients by brute force

Using similar logic for expansion of the other factors, let’s compute a few coefficients, $h_n$ . Recall that $g_0 = 0$ , so we can start the inner sum at $1$ . Now, note that $[z^k] \left( \sum_{m=1}^{\infty} g_m z^m \right)^n$ is essentially asking for the sum of all products $g_{l_1} g_{l_2} \cdots g_{l_n}$ such that $\sum_{i=0}^{n} l_i = k$ . Obviously, each $0 \lt l_i \leq n$ , but we also know that $g_0 = 0$ , so only positive indices actually contribute. Let’s put this in context of the total sum for the coefficient of interest:

\begin{aligned} h_n &= \sum_{m = 0}^{\infty} f_m [z^n] \left( \sum_{k=1}^{\infty} g_k z^k \right)^m \end{aligned}

Now it becomes clearer that each $f_m$ contributes one term for each (strong) $m$ -composition of $n$ .³ Each such integer composition corresponds to a certain product of coefficients of $g(z)$ . For example, $g_1 g_3 g_1$ is a $3$ -composition of $5$ . It has $3$ parts and the product of the implied powers of $z$ is $z^5$ , since $1 + 3 + 1 = 5$ . Consequently this contributes to the coefficient $h_5$ after being multiplied by the outer coefficient $f_3$ .

Here are the first few coefficients fully expanded:

\begin{aligned} h_0 &= f_0 \\ h_1 &= f_1g_1 \\ h_2 &= f_1g_2 + f_2(g_1g_1) \\ h_3 &= f_1g_3 + f_2(g_1g_2 + g_2g_1) + f_3(g_1g_1g_1) \\ h_4 &= f_1g_4 + f_2(g_1g_3 + g_2g_2 + g_3g_1) + f3(g_1g_1g_2 + g_1g_2g_1 + g_2g_1g_1) + f4(g_1g_1g_1g_1) \\ h_5 &= f_1g_5 + f_2(g_1g_4 + g_2g_3 + g_3g_2 + g_4g_1) + f_3(g_1g_1g_3 + g_1g_2g_2 + g_1g_3g_1 + g_2g_1g_2 + g_2g_2g_1 + g_3g_1g_1) + f_4(g_1g_1g_1g_2 + g_1g_1g_2g_1 + g_1g_2g_1g_1 + g_2g_1g_1g_1) + f_5(g_1g_1g_1g_1g_1) \\ \dots \end{aligned}

It may jump out at you that visiting each partition (which grows exponentially in $\sqrt{n}$ ) and multiplying by its corresponding multinomial coefficient would be more efficient than visiting each composition (which grows exponentially in $n$ ). However, this obscures a natural recurrence relation which helps us actually compute the $h_n$ “efficiently”.⁴

Deriving a recurrence

Suppose we’re computing $h_n$ where $n \gt 0$ . This will involve a sum of the form

\begin{aligned} h_n &= \sum_{m=0}^{n} f_m c_{m, n} \end{aligned}

where $c_{m, n}$ is precisely the sum-of-products of coefficients from $g(z)$ that we’ve been discussing. In particular, it counts the coefficient products corresponding to each distinct $m$ -partition of $n$ . We want an explicit recurrence for $c_{m, n}$ .

Note that $c_{1, n} = g_n$ (there’s only one $1$ -composition of $n$ ) and $c_{n, n} = g_1^n$ (there’s only one $n$ -composition of $n$ ).

Now, consider the “composition product sum”, $c_{m, n}$ for $1 \lt m \lt n$ . If you look at the row computing $h_n$ for a fixed $n$ , you can see that each composition product sum of $m$ parts builds off of the product sums of the rows of $m-1$ parts from some row before it. Specifically, consider the operation of prepending a given factor $g_i$ to any $m-1$ composition of $n-i$ for all valid $g_i$ . This generates all valid $m$ -composition products of $n$ . Moreover, since we’re prepending a factor to the each term in a given set of composition products (rather than, for example, replacing one of the previous factors with $g_i$ ), this has the effect of scaling the corresponding sub-composition by $g_i$ .

A concrete example will probably help this stick. Consider the $3$ -compositions of $4$ , i.e. the composition terms of the $f_3$ term that go into $h_4$ . If we choose $g_1$ as the first factor, then this leaves us with $2$ -compositions of $3$ . We read these off of the $f_2$ term of $h_3$ : $g_1 g_2 + g_2 g_1$ . Similarly, we could have a first factor of $g_2$ by combining this with the $2$ -compositions of $2$ . As before, we read these off of the $f_2$ term of $h_2$ : $g_1 g_1$ . Putting these all together, we have the terms of the $f_3$ term of $h_4$ : $g_1 g_1 g_2 + g_1 g_2 g_1 + g_2 g_1 g_1$ . As desired, this matches exactly what we found by brute force under $h_4$ .

We now have our recurrence relation:

c_{m, n} = \sum_{i=1}^{n-m+1} g_i c_{m - 1, n - i}

Note that the sum is over $1 \le i \le 0$ because we can only have valid $m$ -compositions of $n$ when $m \le n$ . Thus we can only have $m-1$ -compositions of $n-i$ if $m-1 \le n-i \Rightarrow i \le n - m + 1$ .

Space and runtime

The recurrence shown above suggests a direct algorithm for computing each $h_n$ :

Maintain a ragged array of $c_{m, n}$ coefficients computed so far.
Bootstrap with $c_{0, 0} = 1$ .
Compute $h_n$ as: $\sum_{m=0}^{n} f_m \sum_{i=1}^{n-m+1} g_i c_{m-1, n-i}$

Note: Each implicit value of $c_{m, n}$ computed above should be stored in the array to avoid exponentially duplicated work.

If we’re computing coefficient $h_n$ , then we effectively have to compute all coefficients $h_i$ for $0 \le i \le n$ , so we’re using $O(n^2)$ storage space for the composition product sum array.⁵

At each step $i$ , the double sum costs us $O(n^2)$ runtime, for a total of $O(n^3)$ to get $h_n$ . While this could potentially work fine up to a few thousand (without many deeper levels of composition), it will quickly become unwieldy for larger $n$ . Unfortunately, I’m not aware of any good solutions in arbitrary precision. You seem to be stuck with lossy DFT-based solutions or asymptotic analysis.

I did, however see many references to advanced or impractical techniques such as those by Brent and Kung, as well as Bell polynomials and Faà di Bruno’s formula. ↩
Any time I use the phrase “integer composition”, this should be understood as a “strong integer composition”, since we do not allow zeros in this context. ↩
To see why, note that we’re taking an infinite multinomial raised to the power $m$ and, consequently, need to select $m$ terms in total whose powers of $z$ sum to $n$ , i.e., the power of interest. ↩
We will not be using FFT-based methods here because these increase the complexity and prerequisites dramatically, but also because the goal is to compute exact numbers in arbitrary precision. In particular, in the combinatorial classes we will be dealing with, the output coefficients are always intergral (except for some intermediate fractional coefficients as observed in function composition). ↩
We’re also using $O(n)$ storage for the original $f_i$ s and $g_i$ s since these are reused while computing each sum. While it doesn’t affect the asymptotics, it does mean that you don’t get much benefit in a setting where the coefficients of $F(z)$ and $G(z)$ are streamed. ↩