Binomial Convolutions for Rational Power Series

Ira M. Gessel^{1} Department of Mathematics Brandeis University Waltham, MA 02453 USA

                  gessel@brandeis.edu

Ishan Kar College of Letters and Science University of California Berkeley, CA 94720 USA

                 ishankar@berkeley.edu

Abstract

The binomial convolution of two sequences (a_{n}) and (b_{n}) is the sequence whose nth term is ∑_{k=0}^{n} \binom{n}{k}a_{k}b_{n−k}. If (a_{n}) and (b_{n}) have rational generating functions, then so does their binomial convolution. We discuss an efficient method, using resultants, for computing this rational generating function and give several examples involving Fibonacci and tribonacci numbers and related sequences. We then describe a similar method for computing Hadamard products of rational generating functions. Finally, we describe two additional methods for computing binomial convolutions and Hadamard products of rational power series, one using symmetric functions and one using partial fractions.

1Supported by a grant from the Simons Foundation (#427060).

1 Introduction

The ordinary convolution of two sequences (α_{n}) and (β_{n}) is the sequence (γ_{n}) defined by γ_{n} = ∑_{k=0}^{n} α_{k}β_{n−k}. Ordinary convolutions are closely related to ordinary generating functions: if (γ_{n}) is the ordinary convolution of (α_{n}) and (β_{n}) then

(∑_{n=0}^{∞} α_{n}x^{n})(∑_{n=0}^{∞} β_{n}x^{n}) = ∑_{n=0}^{∞} γ_{n}x^{n}.

The binomial convolution of two sequences (a_{n}) and (b_{n}) is the sequence (c_{n}) defined by c_{n} = ∑_{k=0}^{n} \binom{n}{k}a_{k}b_{n−k}. The exponential generating function ∑_{n=0}^{∞} c_{n}x^{n}/n! for (c_{n}) is the product of the exponential generating functions for (a_{n}) and (b_{n}). So binomial convolutions of sequences with simple exponential generating functions are easily dealt with.

We are concerned here with binomial convolutions of sequences with rational ordinary generating functions. Let us define an operation ⊙ on formal power series in x with coefficients in a field of characteristic 0 by

∑_{n=0}^{∞} a_{n}x^{n} ⊙ ∑_{n=0}^{∞} b_{n}x^{n} = ∑_{n=0}^{∞} c_{n}x^{n},  (1)

where c_{n} = ∑_{k=0}^{n} \binom{n}{k}a_{k}b_{n−k}. In other words, (1) is equivalent to

(∑_{n=0}^{∞} a_{n} \frac{x^{n}}{n!} )(∑_{n=0}^{∞} b_{n} \frac{x^{n}}{n!} )= ∑_{n=0}^{∞} c_{n} \frac{x^{n}}{n!} .  (2)

We call A(x)⊙B(x) the binomial product of A(x) and B(x). (This operation is sometimes called the Hurwitz product [1, 9] or shuffle product [3].) It follows from the equivalence of (1) and (2) that ⊙ is associative, with identity element 1, and that the inverse with respect to ⊙ of 1/(1 − αx) is 1/(1 + αx). As we will see in Theorem 6, if A(x) and B(x) are rational, then so is A(x) ⊙ B(x). We will discuss several methods for computing binomial products of rational power series and give examples.

In Section 4 we describe a method for computing binomial products of rational power series efficiently, using resultants. In Section 5 we give several examples of binomial convolution formulas involving Fibonacci, tribonacci, and Perrin numbers.

We discuss in Section 6 a similar method, using resultants, for computing Hadamard products of rational power series. Then in Section 7 we discuss a different method, using symmetric functions, for computing binomial products and Hadamard products of rational power series, and in Section 8 we discuss another method for binomial and Hadamard products, using partial fractions to compute constant terms of Laurent series.

1.1 Computing binomial convolutions

For rational power series with denominators of degree 2, it is not difficult to compute binomial convolutions directly. For example, for the Fibonacci numbers F_{n}, with

∑_{n=0}^{∞} F_{n}x^{n} = \frac{x}{1 − x − x^{2}} ,

we have the explicit formula F_{n} = (α^{n} −β^{n})/\sqrt{5}, where α = (1+\sqrt{5})/2 and β = (1−\sqrt{5})/2. Thus since α + 1 = α^{2} and β + 1 = β^{2}, we have

∑_{k=0}^{n} \binom{n}{k}F_{k} = \frac{1}{\sqrt{5}} ∑_{k=0}^{n} \binom{n}{k}(α^{k} − β^{k})
= ((α + 1)^{n} − (β + 1)^{n})/\sqrt{5}
= (α^{2n} − β^{2n})/\sqrt{5} = F_{2n},

as is well known. (Note that our computation does not use the explicit values of α and β, but only the quadratic equation that they satisfy.) In the same way, we can prove the binomial convolution

∑_{k=0}^{n} \binom{n}{k}F_{k}F_{n−k} = \frac{1}{5}(2^{n}L_{n} − 2),  (3)

where L_{n} is the Lucas number F_{n−1} + F_{n+1}, as shown by Church and Bicknell [6].

More challenging are sequences whose generating functions have a denominator of degree 3. One well-known example is the tribonacci sequence, defined by T_{−1} = T_{0} = 0, T_{1} = 1, and T_{n} = T_{n−1} + T_{n−2} + T_{n−3} for n ≥ 2, with generating function

∑_{n=0}^{∞} T_{n}x^{n} = \frac{x}{1 − x − x^{2} − x^{3}} .

(There are several different conventions for the index of the first nonzero tribonacci number. We start with T_{1}, following Komatsu [20], but some authors start with T_{0} or T_{2}.) The tribonacci numbers are sequence \underline{A000073} in the On-Line Encyclopedia of Integer Sequences (OEIS) [15].

Komatsu [21] found a formula, given in (22), for the numbers

∑_{k=0}^{n} \binom{n}{k} T_{k}T_{n−k}  (4)

by using the exponential generating function for T_{n}, which involves the zeros of a cubic polynomial. We will give a simpler proof of Komatsu’s formula, and a generalization, in Section 5.3.

Prodinger [26] proved Komatsu’s formula by first computing the binomial product

\frac{x}{1 − x − x^{2} − x^{3}} ⊙ \frac{x}{1 − x − x^{2} − x^{3}} = \frac{1}{11}(\frac{1 + x + 10x^{2}}{1 − 2x − 4x^{2} − 8x^{3}} − \frac{1 + x − 8x^{2}}{1 − 2x + 2x^{3}} ).

He used the following approach to computing binomial products of rational power series. Suppose we know the generating functions ∑_{n=0}^{∞} a_{n}x^{n} = A(x) and ∑_{n=0}^{∞} b_{n}x^{n} = B(x). Then a straightforward computation gives

∑_{n=0}^{∞} x^{n} ∑_{k=0}^{n} \binom{n}{k}a_{k} ∑_{m=k}^{∞} y^{m}b_{m−k} = \frac{1}{1 − x}A (\frac{xy}{1 − x}) B(y).  (5)

The binomial convolution ∑_{k=0}^{n} \binom{n}{k}a_{k}b_{n−k} is the coefficient of x^{n}y^{n} in (5), so the binomial product A(x) ⊙ B(x) may be obtained by extracting the diagonal from (5). For A(x) and B(x) rational, Prodinger did this using Hautus and Klarner’s residue method [13]. We will discuss a more efficient method for computing the diagonal of (5) in Section 8.

Ekhad and Zeilberger [8] used another efficient method for computing binomial products of rational power series by solving a system of linear equations, using the fact that a proper rational power series in x with a denominator of degree d is determined by the coefficients of x^{i} for 0 ≤ i ≤ 2d. (Recall that a proper rational function is one in which the degree of the numerator is less than the degree of the denominator.)

Cerlienco, Mignotte, and Piras [5, Section A IV 1] sketched another approach, using determinants, for computing binomial and Hadamard products.

Our main method for computing binomial products is based on resultants of polynomials. If the denominator of A(x) is ∏_{i}(1 − α_{i}x) and the denominator of B(x) is ∏_{j}(1 − β_{j}x) then A(x) ⊙ B(x) is a rational function with denominator (not necessarily in lowest terms) ∏_{i,j}(1 − α_{i}β_{j}x), and this product can be computed from the coefficients of the denominators of A(x) and B(x) in terms of a resultant, with no need to factor these denominators into linear factors. We use a similar resultant method to compute the denominator of the Hadamard product of two rational generating functions in Section 6.

Alecci, Barbero, and Murru [1, Theorem 5] have also observed that resultants can be used to compute binomial products.

2 On rational power series over a field

We first prove a useful fact that clarifies the status of rational power series over a field. Another proof was given by Klazar and Horský [19, Theorem 3]. The result can also be derived using properties of determinants; for example, Equation (10) of Cerlienco, Mignotte, and Piras [5] yields a determinantal formula for the denominator of a rational power series in terms of its coefficients.

Theorem 1. Suppose that

A(x) = ∑_{n=0}^{∞} a_{n}x^{n} = \frac{p(x)}{q(x)},

where each a_{n} is in a field F, and p(x) and q(x) are polynomials with coefficients in an extension field G of F. Then there exist polynomials P (x) and Q(x) with coefficients in F such that

A(x) = \frac{P (x)}{Q(x)}.

In order to prove Theorem 1, we first prove two lemmas from linear algebra.

Lemma 2. Let F ⊆ G be fields. If a subset W of F^{n} is linearly independent over F then it is linearly independent over G.

Proof. Let W ⊂ F^{n} be linearly independent over F. Since W is finite, we may let the elements of W be (w_{i,1}, w_{i,2}, …, w_{i,n}), for i = 1, …, |W |, and consider the matrix M = (w_{i,j})_{i=1,...,|W |;j=1,...,n}. We can convert M to a matrix M′ in row-echelon form by operations that preserve the span of its rows over both F and G, and the dimension of the row space of M′ over both F and G is the number of nonzero rows.

Lemma 3. Let F ⊆ G be fields. Suppose that S is a set of vectors in F^{n} and that there is a nonzero vector t in G^{n} such that s · t = 0 for all s ∈ S, where · is the usual dot product defined by (s_{1}, …, s_{n}) · (t_{1}, …, t_{n}) = s_{1}t_{1} + ⋯ + s_{n}t_{n}. Then there is a nonzero vector u ∈ F^{n} such that s · u = 0 for all s ∈ S.

Proof. The hypothesis implies that the dimension of the span of S over G is less than n. Since the dimension of the span of S (over either F or G) is the size of a largest linearly independent subset, Lemma 2 implies that the dimension of S over F is also less than n.

The lemma then follows from a basic fact of linear algebra: If K is a field and U is a set of vectors in K^{n} then the dimension of the span of U over K is less than n if and only if there is a nonzero vector t ∈ K^{n} such that u · t = 0 for every u ∈ U.

Proof of Theorem 1. Let q(x) = q_{0} +q_{1}x+⋯+q_{m}x^{m}. Since q(x) ∑_{n=0}^{∞} a_{n}x^{n} is a polynomial in x, it follows that for some N,

a_{i}q_{m} + a_{i+1}q_{m−1} + ⋯ + a_{i+m}q_{0} = 0

for all i ≥ N. Then by Lemma 3, there exist Q_{0}, …, Q_{m} in F, not all 0, such that

a_{i}Q_{m} + a_{i+1}Q_{m−1} + ⋯ + a_{i+m}Q_{0} = 0

for all i ≥ N. Thus

(Q_{0} + Q_{1}x + ⋯ + Q_{m}x^{m})A(x)

is a polynomial in x with coefficients in F, so A(x) is a quotient of polynomials with coefficients in F.

In what follows, we will be working with power series with coefficients in a field F, and we will derive formulas for these power series that express them as quotients of polynomials with coefficients in an extension field of F. Theorem 1 guarantees that we can express these power series as quotients of polynomials with coefficients in F. Moreover, Theorem 1 allows us to use the term “rational power series over the field F ” to mean both a rational power series with coefficients in F and a power series which is a quotient of polynomials with coefficients in F.

3 Rationality of binomial products

We now prove that the binomial product of two rational power series is rational. There are more direct proofs but our approach proves some useful formulas along the way. For other proofs see Fliess [9, Proposition 3], Bacher [3, Proposition 3.2], and Alecci, Barbero, and Murru [1, Theorem 5]. Another proof, pointed out by Alin Bostan, that the binomial product of rational power series is rational follows from the fact that ∑_{n=0}^{∞} a_{n}x^{n} is rational if and only if ∑_{n=0}^{∞} a_{n}x^{n}/n! is a linear combination of series of the form x^{i}e^{αx}, as these are clearly closed under multiplication.

Lemma 4. For all nonnegative integers j and k and all α, β ∈ F we have

\frac{x^{j}}{(1 − αx)^{j+1}} ⊙ \frac{x^{k}}{(1 − βx)^{k+1}} = \binom{j + k}{j}\frac{x^{j+k}}{(1 − (α + β)x)^{j+k+1}} .  (6)

Proof. For any α ≠ 0 we define f_{n}(α, k) for n, k ∈ ℕ (where ℕ is the set of nonnegative integers) by f_{n}(α, k) = α^{n−k}\binom{n}{k}. We define f_{n}(0, k) to be the limit lim_{α→0} f_{n}(α, k), so

f_{n}(0, k) = \begin{cases}1, & if n = k; \\ 0, & otherwise. \end{cases}

(More precisely, for fixed n and k, f_{n}(α, k) is a polynomial in α. If n < k then this polynomial is identically 0 and if n ≥ k this polynomial is a constant times a nonnegative power of α; in all cases f_{n}(0, k) is obtained by setting α = 0 in this polynomial.) It is easy to check that

∑_{n=0}^{∞} f_{n}(α, k)x^{n} = \frac{x^{k}}{(1 − αx)^{k+1}}

and

∑_{n=0}^{∞} f_{n}(α, k)\frac{x^{n}}{n!} = \frac{x^{k}}{k!} e^{αx}

for all α, including α = 0. Thus

(∑_{n=0}^{∞} f_{n}(α, j)\frac{x^{n}}{n!} )(∑_{n=0}^{∞} f_{n}(β, k)\frac{x^{n}}{n!} ) = \binom{j + k}{j}\frac{x^{j+k}}{(j + k)!}e^{(α+β)x},

and (6) follows from the equivalence of (1) and (2).

We note a useful consequence of (6).

Corollary 5. For any power series A(x) and any β ∈ F, we have

A(x) ⊙ \frac{1}{1 − βx} = \frac{1}{1 − βx}A (\frac{x}{1 − βx}) .

Proof. By linearity, it is sufficient to prove the formula for A(x) = x^{j}, which is the case α = 0, k = 0 of (6).

Theorem 6. Let A(x) and B(x) be rational power series over F. We may express A(x) as

A(x) = \frac{R(x)}{∏_{i=1}^{m}(1 − α_{i}x)},

where R(x) is a polynomial of degree less than m and the α_{i} lie in some extension field of F, but they need not be distinct nor nonzero. We may express B(x) similarly as

B(x) = \frac{S(x)}{∏_{j=1}^{n}(1 − β_{j}x)}.

Then

A(x) ⊙ B(x) = \frac{T (x)}{∏_{i=1}^{m} ∏_{j=1}^{n}(1 − (α_{i} + β_{j})x)}.  (7)

for some polynomial T (x) of degree at most mn − 1.

Proof. We first find the usual partial fraction decomposition of A(1/x)/x,

A(1/x)/x = U(x) + ∑_{i=1}^{M} ∑_{k=0}^{e_{i}} \frac{A_{ik}}{(x − a_{i})^{k+1}}

where U(x) is a polynomial, the a_{i} and A_{ik} lie in some extension field of F, the a_{i} are distinct (one of them may be 0), and the e_{i} are nonnegative integers. Replacing x with 1/x and dividing by x gives us

A(x) = \frac{U(1/x)}{x} + ∑_{i=1}^{M} ∑_{k=0}^{e_{i}} \frac{A_{ik}x^{k}}{(1 − a_{i}x)^{k+1}} .  (8)

Since A(x) has a power series expansion, U(x) must be 0. (We note that A(x) is proper if and only if a_{i} ≠ 0 for all i.)

Similarly we have the partial fraction expansion

B(x) = ∑_{j=1}^{N} ∑_{l=0}^{f_{j}} \frac{B_{jl}x^{l}}{(1 − b_{j}x)^{l+1}} .

Then

A(x) ⊙ B(x) = ∑_{i=1}^{M} ∑_{j=1}^{N} ∑_{k=0}^{e_{i}} ∑_{l=0}^{f_{j}} A_{ik}B_{jl}\binom{k + l}{k}\frac{x^{k+l}}{(1 − (a_{i} + b_{j})x)^{k+l+1}} by (6)

= ∑_{i=1}^{M} ∑_{j=1}^{N} \frac{P_{ij}(x)}{(1 − (a_{i} + b_{j})x)^{e_{i}+f_{j}+1}}  (9)

where P_{ij}(x) is a polynomial of degree at most e_{i} + f_{j}.

It follows from (9) that

A(x) ⊙ B(x) = \frac{Q(x)}{∏_{i=1}^{M} ∏_{j=1}^{N}(1 − (a_{i} + b_{j})x)^{e_{i}+f_{j}+1}}  (10)

for some polynomial Q(x) of degree less than ∑_{i=1}^{M} ∑_{j=1}^{N}(e_{i} + f_{j} + 1). In particular, A(x) ⊙ B(x) is rational.

Now choose a sequence α_{1}, …, α_{m}, where m = ∑_{i=1}^{M}(e_{i} + 1), consisting of e_{i} occurrences of a_{i} for each i, so that ∏_{i=1}^{m}(1 − α_{i}x) = ∏_{i=1}^{M}(1 − a_{i}x)^{e_{i}+1} and similarly choose β_{1}, …, β_{n}, where n = ∑_{j=1}^{N}(f_{j} + 1), so that ∏_{j=1}^{n}(1 − β_{j}x) = ∏_{j=1}^{N}(1 − b_{j}x)^{f_{j}+1}. Then

∏_{i=1}^{m} ∏_{j=1}^{n}(1 − (α_{i} + β_{j})x) = ∏_{i=1}^{M} ∏_{j=1}^{N}(1 − (a_{i} + b_{j})x)^{(e_{i}+1)(f_{j}+1)}.

Since (e_{i} + 1)(f_{j} + 1) − (e_{i} + f_{j} + 1) = e_{i}f_{j} ≥ 0, multiplying the numerator and denominator of the right side of (10) by ∏_{i=1}^{M} ∏_{j=1}^{N}(1 − (a_{i} + b_{j})x)^{e_{i}f_{j}} gives the desired representation (7).

We may restate Theorem 6 in a less elegant but more computationally practical way by separating out the α_{i} and β_{j} that are equal to 0.

Corollary 7. Let

A(x) = \frac{R(x)}{∏_{i=1}^{m}(1 − α_{i}x)},

and

B(x) = \frac{S(x)}{∏_{j=1}^{n}(1 − β_{j}x)},

where R(x) and S(x) are polynomials and the α_{i} and β_{j} are all nonzero. Let

u = max(deg R(x) + 1 − m, 0) and v = max(deg S(x) + 1 − n, 0).

Then

A(x) ⊙ B(x) = \frac{T (x)}{(∏_{i=1}^{m}(1 − α_{i}x))^{v}(∏_{j=1}^{n}(1 − β_{j}x))^{u} ∏_{i=1}^{m} ∏_{j=1}^{n}(1 − (α_{i} + β_{j})x)}. (11)

for some polynomial T (x) of degree less than (u + m)(v + n).

Note that A(x) ⊙ B(x) need not be proper even if A(x) and B(x) are; for example by Lemma 6 we have

\frac{x^{j}}{(1 − x)^{j+1}} ⊙ \frac{x^{k}}{(1 + x)^{k+1}} = \binom{j + k}{j}x^{j+k}.

As is clear from the proof of Theorem 6, the denominator in lowest terms of A(x) ⊙ B(x) will in many cases be a proper divisor of that given by (7) or (11), but these formulas, together with Theorem 10, which allows us to compute ∏_{i=1}^{m} ∏_{j=1}^{n}(1 − (α_{i} + β_{j})x) without factoring the denominators of A(x) and B(x), are efficient enough in all of our applications.

One important special case is when A(x) and B(x) have the same denominator. If we assume, for simplicity, that A(x) and B(x) are proper and have no repeated zeros, then the partial fraction expansion (9) shows that A(x) ⊙ B(x) can be written with denominator ∏_{i=1}^{n}(1 − 2α_{i}x) ∏_{1≤i<j≤n}(1 − (α_{i} + α_{j})x).

As a first example, with proper rational functions and no repeated factors, we have

\frac{x}{(1 − x)(1 − 2x)} ⊙ \frac{x}{(1 − 3x)(1 − 5x)} = \frac{x^{2}(2 − 11x)}{(1 − 4x)(1 − 5x)(1 − 6x)(1 − 7x)}.

An example with an improper rational function is

\frac{x^{3}}{1 − x} ⊙ \frac{1}{1 − 2x} = \frac{x^{3}}{(1 − 2x)^{3}(1 − 3x)}.

With the notation of Corollary 7, here we have m = 1, n = 1, u = 3, and v = 0. Finally, a more complicated example is

\frac{x^{2}}{(1 − x)^{2}} ⊙ \frac{x^{2}}{(1 − 2x)^{2}} = \frac{x^{4}(6 − 30x + 49x^{2} − 27x^{3})}{(1 − x)^{2}(1 − 2x)^{2}(1 − 3x)^{3}}.

Here the denominator given by Corollary 7 is (1 − x)^{2}(1 − 2x)^{2}(1 − 3x)^{4} but because of the multiple factors in the denominators (1 − x)^{2} and (1 − 2x)^{2}, the actually denominator is a proper divisor of this product.

An alternative approach to binomial products of improper rational power series is to decompose a rational function as a polynomial plus a proper rational function. For proper rational power series we can take all the α_{i} and β_{j} in Theorem 6 to be nonzero; equivalently, in Corollary 7 we have u = v = 0. We then need to compute binomial products of polynomials with rational power series. By linearity, it is sufficient to compute x^{m} ⊙ A(x) in terms of A(x). The next lemma tells us how to do this.

Lemma 8. Let A(x) be a power series in x. Then for any nonnegative integer m we have

x^{m} ⊙ A(x) = \frac{x^{m}}{m!} \frac{d^{m}}{dx^{m}} (x^{m}A(x)).  (12)

Proof. By linearity, it is sufficient to prove (12) in the case A(x) = x^{n}. In this case we have

\frac{x^{m}}{m!} \frac{d^{m}}{dx^{m}} (x^{m}A(x)) = \frac{x^{m}}{m!} \frac{d^{m}}{dx^{m}} x^{m+n} = \binom{m + n}{m}x^{m+n} = x^{m} ⊙ A(x).

4 Resultants

4.1 The resultant of two polynomials

Given two polynomials

A(x) = ∑_{i=0}^{m} a_{i}x^{i} = a_{m} ∏_{i=1}^{m}(x − α_{i})

and

B(x) = ∑_{j=0}^{n} b_{j}x^{k} = b_{n} ∏_{j=1}^{n}(x − β_{j}),

their resultant with respect to the variable x may be defined by

Res(A(x), B(x), x) = a_{m}^{n}b_{n}^{m} ∏_{i,j}(α_{i} − β_{j}).  (13)

It is well known [12] that Res(A(x), B(x), x) can be computed as a determinant. Lemma 9. The resultant can be computed by

Res(A(x), B(x), x) = det Syl(A(x), B(x)),  (14)

where Syl(A(x), B(x)) is the Sylvester matrix of size (n + m) × (n + m),

Syl(A(x), B(x)) = [\underbrace{\begin{aligned}a_{m} & ⋯ & a_{0} \\ ... & & ... \\ & a_{m} & ⋯ & a_{0} \\ b_{n} & ⋯ & b_{0} \\ ... & & ... \\ & b_{n} & ⋯ & b_{0} \end{aligned}}_{m+n}]\begin{matrix}\} n \\ \} m \end{matrix}

Note that unlike formula (13), the determinant formula (14) expresses the resultant in terms of the coefficients of A and B.

4.2 Computing the binomial product with resultants

We now describe how to compute the product ∏_{i=1}^{m} ∏_{j=1}^{n}(1 − (α_{i} + β_{j})x) in Corollary 7 without factoring the denominators. Theorem 10. Let U(x) = ∏_{i=1}^{m}(1 − α_{i}x) and V (x) = ∏_{j=1}^{n}(1 − β_{j}x). Then

∏_{i=1}^{m} ∏_{j=1}^{n}(1 − (α_{i} + β_{j})x) = (−1)^{mn} Res ((1 − y)^{m}U(\frac{x}{1 − y} ), y^{n}V (\frac{x}{y}), y) .

Proof. We have

(1 − y)^{m}U(\frac{x}{1 − y} ) = ∏_{i=1}^{m}(1 − y − α_{i}x) = (−1)^{m} ∏_{i=1}^{m}(y − (1 − α_{i}x))

and

y^{n}V (x/y) = ∏_{j=1}(y − β_{j}x),

so the result follows from (13).

Corollary 7 and Theorem 10 together give us a procedure for computing binomial products of rational power series that can be efficiently implemented on a computer algebra system. By Corollary 7 we can express the binomial convolution A(x)⊙B(x) as a quotient of polynomials T (x)/D(x), where D(x) may be computed explicitly with the help of Theorem 10, and we have a bound on the degree of T (x), deg T (x) ≤ d. Then T (x) can be computed from

T (x) = (A(x) ⊙ B(x))D(x)

by computing the coefficients up to x^{d} on the right side.

As a simple example, let us compute in detail the binomial product of the generating functions for the Fibonacci and Pell sequences. The generating function for the Fibonacci sequence is A(x) = x/(1 − x − x^{2}) and that for the Pell sequence is B(x) = x/(1 − 2x − x^{2}). Since both are proper, by Corollary 7 we know that the binomial product A(x) ⊙ B(x) may be expressed as T (x)/D(x) where T (x) has degree at most 3 and D(x) is a polynomial of degree 4 computed by Theorem 10.

With the notation of Theorem 10, we have m = n = 2, U(x) = 1 − x − x^{2}, and V (x) = 1 − 2x − x^{2}, so

(1 − y)^{2}U (\frac{x}{1 − y} ) = y^{2} + (x − 2)y − x^{2} − x − 1

and

y^{2}V (x/y) = y^{2} − 2xy − x^{2}.

The Sylvester matrix is

[\begin{matrix}1 & x − 2 & −x^{2} − x + 1 & 0 \\ 0 & 1 & x − 2 & −x^{2} − x + 1 \\ 1 & −2x & −x^{2} & 0 \\ 0 & 1 & −2x & −x^{2} \end{matrix}]

with determinant D(x) = 1 − 6x + 7x^{2} + 6x^{3} − 9x^{4}. We compute directly that

A(x) ⊙ B(x) = 2x^{2} + 9x^{3} + 40x^{4} + ⋯

So the numerator is

T (x) = (2x^{2} + 9x^{3} + 40x^{4} + ⋯)(1 − 6x + 7x^{2} + 6x^{3} − 9x^{4}) = 2x^{2} − 3x^{3}.

Thus

\frac{x}{1 − x − x^{2}} ⊙ \frac{x}{1 − 2x − x^{2}} = \frac{2x^{2} − 3x^{3}}{1 − 6x + 7x^{2} + 6x^{3} − 9x^{4}} .

5 Examples

5.1 Fibonacci numbers

We first discuss some binomial convolution formulas for Fibonacci numbers. We omit the details of the computations of binomial products, which are done by the computer using the procedure described in Section 4.2, together with partial fraction expansion.

Recall that the Fibonacci numbers have the generating function

∑_{n=0}^{∞} F_{n}x^{n} = \frac{x}{1 − x − x^{2}}  (15)

and the Lucas numbers L_{n} = F_{n−1} + F_{n+1} have the generating function

∑_{n=0}^{∞} L_{n}x^{n} = \frac{2 − x}{1 − x − x^{2}}

Let us first prove the identity of Church and Bicknell [6] mentioned in Section 1.1. (See also Prodinger [26].) We find that

\frac{x}{1 − x − x^{2}} ⊙ \frac{x}{1 − x − x^{2}} = \frac{2x^{2}}{1 − 3x − 2x^{2} + 4x^{3}} = \frac{1}{5} (\frac{2 − 2x}{1 − 2x − 4x^{2}} − \frac{2}{1 − x}) .

Since

\frac{2 − 2x}{1 − 2x − 4x^{2}} = \frac{2 − (2x)}{1 − (2x) − (2x)^{2}} = ∑_{n=0}^{∞} 2^{n}L_{n}x^{n},

we have

∑_{k=0}^{n} \binom{n}{k}F_{k}F_{n−k} = \frac{1}{5}(2^{n}L_{n} − 2).  (16)

This identity may be generalized in many ways. We describe one generalization, based on the multisected Fibonacci generating functions [14]

∑_{n=0}^{∞} F_{pn+q}x^{n} = \frac{F_{q} + (−1)^{q}F_{p−q}x}{1 − L_{p}x + (−1)^{p}x^{2}} .  (17)

and

∑_{n=0}^{∞} L_{pn+q}x^{n} = \frac{L_{q} − (−1)^{q}L_{p−q}x}{1 − L_{p}x + (−1)^{p}x^{2}} .  (18)

Here p and q may be arbitrary integers, and the Fibonacci and Lucas numbers are extended to negative integer indices by F_{−n} = (−1)^{n−1}F_{n} and L_{−n} = (−1)^{n}L_{n}. Letting f(x) = x/(1 − x − x^{2}), we find that

5f(ax) ⊙ f(bx) = \frac{2 − (a + b) x}{1 − (a + b) x − (a + b)^{2} x^{2}} − \frac{2 − (a + b) x}{1 − (a + b) x − (a^{2} − 3ab + b^{2}) x^{2}} . (19)

The first term on the right in (19) is

∑_{n=0}^{∞}(a + b)^{n}L_{n}x^{n}.

We would like to choose a and b so that the denominator in the second term on the right in (19) is 1 − L_{p}(cx) + (−1)^{p}(cx)^{2} for some p and c, so that we can apply (17) or (18). A calculation that we omit suggests that we should take a = F_{p−1} and b = F_{p+1} for some p. Once we have these values for a and b, it is easy to verify the result: We have a + b = F_{p−1} + F_{p+1} = L_{p} and

a^{2} − 3ab + b^{2} = (b − a)^{2} − ab
= (F_{p+1} − F_{p−1})^{2} − F_{p+1}F_{p−1} = F_{p}^{2} − F_{p+1}F_{p−1}
= (−1)^{p−1}

by Cassini’s identity for Fibonacci numbers. Thus

5f(F_{p−1}x) ⊙ f(F_{p+1}x) = \frac{2 − L_{p}x}{1 − L_{p}x − L_{p}^{2}x^{2}} − \frac{2 − L_{p}x}{1 − L_{p}x + (−1)^{p}x^{2}}

which by (15) and (18) is equal to

∑_{n=0}^{∞} L_{p}^{n}F_{n}x^{n} − ∑_{n=0}^{∞} L_{pn}x^{n}.

Thus

∑_{k=0}^{n} \binom{n}{k}F_{p−1}^{k}F_{p+1}^{n−k}F_{k}F_{n−k} = \frac{1}{5}(L_{p}^{n}L_{n} − L_{pn}).  (20)

Church and Bicknell’s identity (16) is the case p = 0 of (20).

As another example, we have

\frac{x}{1 − x − x^{2}} ⊙ \frac{1}{1 − x^{2}} = \frac{x − x^{2}}{1 − 2x − 3x^{2} + 4x^{3} − x^{4}}

= \frac{1}{2} (\frac{x}{1 + x − x^{2}} + \frac{x}{1 − 3x + x^{2}} )

= \frac{1}{2} ∑_{n=0}^{∞}(−1)^{n−1}F_{n}x^{n} + \frac{1}{2} ∑_{n=0}^{∞} F_{2n}x^{n},

where we have used the case p = 2, q = 0 of (17), so

∑_{k=0}^{⌊n/2⌋} \binom{n}{2k}F_{n−2k} = \frac{1}{2}((−1)^{n−1}F_{n} + F_{2n}).

The squares of the Fibonacci numbers have the well-known generating function

∑_{n=0}^{∞} F_{n}^{2}x^{n} = \frac{x − x^{2}}{1 − 2x − 2x^{2} + x^{3}}  (21)

which we will rederive as (29). We have

\frac{x − x^{2}}{1 − 2x − 2x^{2} + x^{3}} ⊙ \frac{10}{1 − 5x^{2}} = \frac{x − 3x^{2} − 2x^{4} − 6x^{5}}{1 − 4x − 15x^{2} + 50x^{3} + 35x^{4} − 114x^{5} + 36x^{6}}

= \frac{2 − 3x}{1 − 3x + x^{2}} − \frac{4 + 4x}{1 + 2x − 4x^{2}} + \frac{2 − 3x}{1 − 3x − 9x^{2}}

= ∑_{n=0}^{∞} L_{2n}x^{n} + ∑_{n=0}^{∞}(−2)^{n+1}L_{n}x^{n} + ∑_{n=0}^{∞} 3^{n}L_{n}

and we find the binomial convolution

10 ∑_{k=0}^{⌊n/2⌋} \binom{n}{2k}5^{k}F_{n−2k}^{2} = L_{2n} + (3^{n} + (−2)^{n+1})L_{n}.

5.2 Second order recurrent sequences

Let g(x) = (c + dx)/(1 − ax − bx^{2}). We find that

g(x) ⊙ g(x) = \frac{p}{1 − ax} + \frac{q + rx}{1 − 2ax − 4bx^{2}} ,

where p, q, and r are certain rational functions of a, b, c, and d, too complicated to be worth writing out here. If we put some restrictions on the parameters, we get a much nicer formula.

Theorem 11. Let

g(x) = ∑_{n=0}^{∞} G_{n}x^{n} = \frac{2 − ax}{1 − ax − bx^{2}} .

Then

g(x) ⊙ g(x) = \frac{2}{1 − ax} + g(2x),

so

∑_{k=0}^{n} \binom{n}{k}G_{k}G_{n−k} = 2a^{n} + 2^{n}G_{n}.

Special cases of Theorem 11 are a = b = 1 (Lucas numbers; this identity was given by Church and Bicknell [6]), a = 1, b = 2 (Jacobsthal-Lucas numbers, \underline{A014551}), and a = 2, b = 1 (companion Pell numbers, \underline{A002203}).

5.3 Tribonacci numbers

As noted in the introduction, the tribonacci numbers T_{n} may be defined by the generating function

∑_{n=0}^{∞} T_{n}x^{n} = \frac{x}{1 − x − x^{2} − x^{3}} .

Komatsu [20] (see also Komatsu and Li [22] and Komatsu [21]) gave the binomial convolution formula

∑_{k=0}^{n} \binom{n}{k}T_{k}T_{n−k} = \frac{1}{22}(2^{n}T_{n}^{(2,3,10)} + 2 ∑_{k=0}^{n} \binom{n}{k}(−1)^{k}T_{k}^{(−1,2,7)}).  (22)

Here the numbers T_{n}^{(s_{0},s_{1},s_{2})} satisfy the tribonacci recurrence T_{n} = T_{n−1} + T_{n−2} + T_{n−3} for n ≥ 3 with initial values T_{0}^{(s_{0},s_{1},s_{2})} = s_{0}, T_{1}^{(s_{0},s_{1},s_{2})} = s_{1}, and T_{2}^{(s_{0},s_{1},s_{2})} = s_{2}, so the ordinary tribonacci numbers are T_{n} = T_{n}^{(0,1,1)}. It is easy to check that the generating function for T_{n}^{(s_{0},s_{1},s_{2})} is

∑_{n=0}^{∞} T_{n}^{(s_{0},s_{1},s_{2})}x^{n} = \frac{a + bx + cx^{2}}{1 − x − x^{2} − x^{3}} ,

where a = s_{0}, b = s_{1} −s_{0}, and c = s_{2} −s_{1} −s_{0}. (Equivalently, T_{n}^{(s_{0},s_{1},s_{2})} = aT_{n+1} +bT_{n} +cT_{n−1} with these values of a, b, and c.)

We now derive Komatsu’s formula with binomial products. Let t(x) = x/(1−x−x^{2} −x^{3}) be the tribonacci generating function. Then we compute

t(x) ⊙ t(x) = \frac{2x^{2}(1 − x − x^{2} − 2x^{3})}{1 − 4x + 2x^{3} + 12x^{4} − 8x^{5} − 16x^{6}}
= \frac{1 + x + 10x^{2}}{11(1 − 2x − 4x^{2} − 8x^{3})} − \frac{1 + x − 8x^{2}}{11(1 − 2x + 2x^{3})},  (23)

as found by another method by Prodinger [26]. The first term in (23) is

\frac{1}{22} ∑_{n=0}^{∞} 2^{n}T_{n}^{(2,3,10)}x^{n}.

For the second term in (23) we can check easily that

\frac{1 + x − 8x^{2}}{1 − 2x + 2x^{3}} = \frac{1}{1 − x} ⊙ \frac{1 + 3x − 6x^{2}}{1 + x − x^{2} + x^{3}}

and that

\frac{1 + 3x − 6x^{2}}{1 + x − x^{2} + x^{3}} = ∑_{n=0}^{∞} T_{n}^{(1,−2,−7)}(−x)^{n} = − ∑_{n=0}^{∞} T_{n}^{(−1,2,7)}(−x)^{n},

proving Komatsu’s formula.

There is nothing really special about tribonacci numbers here; something similar will hold in general for generating functions with a cubic denominator:

Theorem 12. Let r(x) and s(x) be proper rational functions with the same denominator D(x) = 1+Ax+Bx^{2} +Cx^{3}, where C ≠ 0. Suppose that D(x) = (1−α_{1}x)(1−α_{2}x)(1−α_{3}x) where 2α_{i} ≠ α_{j} + α_{k} for all i, j, k ∈ {1, 2, 3}. Then for some polynomials u(x) and v(x) of degree at most 2,

r(x) ⊙ s(x) = \frac{u(x)}{D(2x)} + \frac{1}{1 − Ax} ⊙ \frac{v(x)}{D(−x)}.  (24)

Proof. The denominator of r(x) ⊙ s(x) may be factored as D_{1}(x)D_{2}(x), where D_{1}(x) = (1−2α_{1}x)(1−2α_{2}x)(1−2α_{3}x) and D_{2}(x) = (1−(α_{1} +α_{2})x)(1−(α_{1} +α_{3})x)(1−(α_{2} +α_{3})x). The condition on the α_{i} implies that D_{1}(x) and D_{2}(x) are relatively prime, so we have the partial fraction expansion,

r(x) ⊙ s(x) = \frac{u(x)}{D_{1}(x)} + \frac{w(x)}{D_{2}(x)}

for some polynomials u(x) and w(x) of degree at most 2. Clearly D_{1}(x) = D(2x). Since α_{1} + α_{2} + α_{3} = −A, we have

D_{2}(x) = (1 + (α_{1} + A)x)(1 + (α_{2} + A)x)(1 + (α_{3} + A)x).

Therefore 1/(1 + Ax) ⊙ w(x)/D_{2}(x) has denominator D(−x) and so may be expressed as v(x)/D(−x) for some v(x) of degree at most 2; i.e.,

\frac{1}{1 + Ax} ⊙ \frac{w(x)}{D_{2}(x)} = \frac{v(x)}{D(−x)}

Since 1/(1 − Ax) and 1/(1 + Ax) are inverses with respect to ⊙, this implies that

\frac{w(x)}{D_{2}(x)} = \frac{1}{1 − Ax} ⊙ \frac{v(x)}{D(−x)}.

Note that the condition on the α_{i} is equivalent to D(x) having no repeated zeros and D_{1}(x) and D_{2}(x) being relatively prime, and this will be easy to check in any particular example. If D_{1}(x) are D_{2}(x) are not relatively prime we will have a similar but slightly different partial fraction expansion.

5.4 Perrin numbers

Another interesting third-order recurrent sequence is the Perrin sequence (\underline{A001608}) with generating function

P (x) = ∑_{n=0}^{∞} P_{n}x^{n} = \frac{3 − x^{2}}{1 − x^{2} − x^{3}} .

Here the parameter A of Theorem 12 is 0, so (24) will simplify. However, the result turns out to be even nicer than we have any reason to expect.

We have

P (x) ⊙ P (x) = 3\frac{3 − 11x^{2} − 15x^{3} + 4x^{4} + 4x^{5}}{1 − 5x^{2} − 7x^{3} + 4x^{4} + 4x^{5} − 8x^{6}}

= \frac{3 − 4x^{2}}{1 − 4x^{2} − 8x^{3}} + \frac{6 − 2x^{2}}{1 − x^{2} + x^{3}}

= P (2x) + 2P (−x).

Thus we have the unexpectedly simple formula

∑_{k=0}^{n} \binom{n}{k}P_{k}P_{n−k} = (2^{n} + 2(−1)^{n})P_{n}.

More generally, if we set Q(x) = (3 − x^{2})/(1 − x^{2} − ax^{3}), where a is arbitrary, then Q(x) ⊙ Q(x) = Q(2x) + 2Q(−x).

Surprisingly, the Jacobsthal number generating function J(x) = x/(1 − x − 2x^{2}) satisfies a very similar formula, 3J(x) ⊙ J(x) = J(2x) + 2J(−x).

We give one more curious identity, involving a fourth-order recurrent sequence. Let

R(x) = \frac{1 − 2x^{3}}{1 − 8x^{3} + 4x^{4}} .

Then

R(x) ⊙ R(x) = \frac{1}{4}(R(2x) + P (4x^{2})),

where P (x) is the Perrin sequence generating function.

6 The Hadamard product of rational power series

The Hadamard product [16] of the power series A(x) = ∑_{n=0}^{∞} a_{n}x^{n} and B(x) = ∑_{n=0}^{∞} b_{n}x^{n} is

defined by

A(x) ∗ B(x) = ∑_{n=0}^{∞} a_{n}b_{n}x^{n}.  (25)

As we shall see, Hadamard products have some similarity to binomial products. We give here a brief account of the computation of Hadamard products of rational power series using resultants. Another approach to Hadamard products of rational power series using determinants has been given by Potekhina and Tolovikov [25, 24].

It is well known and easy to prove (see, e.g., Stanley [28, Theorem 4.1.1]) that f(x) is a proper rational power series if and only if the coefficient of x^{n} in f(x) is of the form ∑_{i} P_{i}(n)α_{i}^{n} where each P_{i}(n) is a polynomial in n and each α_{i} is nonzero. Clearly the product of two functions of this type is also of this type, so the Hadamard product of proper rational power series is a proper rational power series, and it follows easily that the Hadamard product of any two rational power series is rational.

There is an interesting explicit formula, somewhat analogous to Lemma 4, for a special case of Hadamard products. We have

\frac{x^{i}}{(1 − ax)^{m+1}} ∗ \frac{x^{j}}{(1 − bx)^{n+1}} = ∑_{k=0}^{∞} \binom{m + k − i}{k − i}\binom{n + k − j}{k − j}a^{k−i}b^{k−j}x^{k}

= \frac{∑_{k=0}^{∞} \binom{m+j−i}{k−i}\binom{n+i−j}{k−j}a^{k−i}b^{k−j}x^{k}}{(1 − abx)^{m+n+1}} .  (26)

The numerator of (26) may be expressed as a Hadamard product; it is x^{i}(1 + ax)^{m+j−i} ∗ x^{j}(1 + bx)^{n+i−j}.

Equation (26) is equivalent to a classical formula of Euler for hypergeometric series [2, p. 2, Equation (2)]. See also Kar [16, Theorem 4.1] for another proof of the case i = j = 0. This formula actually holds for all m and n, where the binomial coefficient \binom{a}{b} is defined to be a(a − 1) ⋯ (a − b + 1)/b! if b is a nonnegative integer, and is 0 otherwise.

If m and n are nonnegative integers then the sum in the numerator of (26) has only finitely many nonzero terms. In particular, if m and n are nonnegative integers, i ≤ m + j, and j ≤ n + i, which is always the case when (26) represents a Hadamard product of proper rational functions, then the nonzero terms in the numerator range from k = max(i, j) to k = min(n + i, m + j).

It follows from (26) (or otherwise) that if A(x) is a proper rational function with denominator (1 − ax)^{m+1} and B(x) is a proper rational function with denominator (1 − bx)^{n+1} then A(x) ∗ B(x) is a proper rational function that may be written with denominator (1 − abx)^{(m+1)(n+1)}. Thus we have the following analogue of Theorem 6.

Theorem 13. Let A(x) be a rational power series with denominator ∏_{i=1}^{m}(1 − α_{i}x) and let B(x) be a rational power series with denominator ∏_{j=1}^{n}(1 − β_{j}x). Then A(x) ∗ B(x) is a rational power series that may be written with denominator ∏_{i=1}^{m} ∏_{j=1}^{n}(1−α_{i}β_{j}x). Moreover, if A(x) and B(x) are proper then so is A(x) ∗ B(x).

We now give the analogue for Hadamard products of Theorem 10.

Theorem 14. Let U(x) = ∏_{i=1}^{m}(1 − α_{i}x) and V (x) = ∏_{j=1}^{n}(1 − β_{j}x). Then

∏_{i=1}^{m} ∏_{j=1}^{n}(1 − α_{i}β_{j}x) = (−1)^{mn} Res(U(y), y^{n}V (x/y), y).

Proof. We have

U(y) = (−1)^{m} ∏_{i=1}^{m} α_{i} ∏_{i=1}^{m}(y − α_{i}^{−1})

and

y^{n}V (x/y) = ∏_{j=1}^{m}(y − β_{j}x).

Thus

Res(U(y), y^{n}V (x/y), y) = (−1)^{mn}(∏_{i=1}^{m} α_{i})^{n} ∏_{i=1}^{m} ∏_{j=1}^{n}(α_{i}^{−1} − β_{j}x)

= (−1)^{mn} ∏_{i=1}^{m} ∏_{j=1}^{n}(1 − α_{i}β_{j}x).

As with binomial products, the numerator of A(x) ∗ B(x) can be found by multiplying the denominator with the first mn terms of the infinite series expansion of A(x) ∗ B(x).

As an example, we use this method to find the Hadamard product of

A(x) = \frac{x}{1 − ax − bx^{2}}

and

B(x) = \frac{x}{1 − cx − dx^{2}} .

Examples of well-known second-order recurrent sequences with generating functions of this type are Fibonacci (\underline{A000045}), Pell (\underline{A000129}), and Jacobsthal (\underline{A001045}).

To compute A(x) ∗ B(x) we take U(x) = 1 − ax − bx^{2} = (1 − α_{1}x)(1 − α_{2}x) and V (x) = 1−cx−dx^{2} = (1−β_{1}x)(1−β_{2}x) in Theorem 14. For the denominator of A(x)∗B(x) we have

∏_{i=1}^{2} ∏_{j=1}^{2}(1 − α_{i}β_{j}x) = Res(1 − ay − by^{2}, y^{2} − cxy − dx^{2}, y)
= 1 − acx − (a^{2}d + bc^{2} + 2bd)x^{2} − abcdx^{3} + b^{2}d^{2}x^{4}.  (27)

Multiplying the denominator (27) by the initial terms of

A(x) ∗ B(x) = x + acx^{2} + (a^{2} + b)(c^{2} + d)x^{3} + (a^{3} + 2ab)(c^{3} + 2cd)x^{4} + ⋯

we find the numerator x − bdx^{3}. So

A(x) ∗ B(x) = \frac{x − bdx^{3}}{1 − acx − (a^{2}d + bc^{2} + 2bd)x^{2} − abcdx^{3} + b^{2}d^{2}x^{4}} .  (28)

A combinatorial proof of (28) was given by Shapiro [27] and further combinatorial proofs of similar Hadamard product identities were given by Kim [17, 18].

As a special case of (28) we have the well-known generating function for the squares of the Fibonacci numbers used in Section 5,

∑_{n=0}^{∞} F_{n}^{2}x^{n} = \frac{x − x^{3}}{1 − x − 4x^{2} − x^{3} + x^{4}} = \frac{x − x^{2}}{1 − 2x − 2x^{2} + x^{3}} .  (29)

Frontczak, Goy, and Shattuck [10, 11] give several formulas for Hadamard products of two or three second-order recurrent sequences that can be easily be derived by our methods.

7 Symmetric functions

In this section we briefly describe an alternative approach to computing denominators for binomial products and Hadamard products of rational power series. This method has been used by Dvornicich and Traverso [7] to compute the resultants corresponding to these denominators. Bostan, Flajolet, Salvy, and Schost [4] reformulated this approach with a view to computational efficiency, obtaining a “nearly optimal” algorithm for computing these resultants.

Let U(x) = ∏_{i=1}^{m}(1 − α_{i}x) and V (x) = ∏_{j=1}^{n}(1 − β_{j}x). We want to express ∏_{i,j}(1 − (α_{i} + β_{j})x) and ∏_{i,j}(1 − α_{i}β_{j}x) in terms of the coefficients of U(x) and V (x).

Let e_{k}(u_{1}, …, u_{n}) denote the ith elementary symmetric polynomial in u_{1}, …, u_{n},

e_{k}(u_{1}, …, u_{n}) = ∑_{1≤i_{1}<i_{2}<···<i_{k}≤n} u_{i_{1}}u_{i_{2}} … u_{i_{n}},

with e_{0}(u_{1}, …, u_{n}) = 1. Then the coefficient of x^{k} in ∏_{i=1}^{n}(1 − u_{i}x) is equal to

(−1)^{k}e_{k}(u_{1}, …, u_{n}),

and e_{k}(u_{1}, …, u_{n}) = 0 for k > n. Thus to find the denominator polynomials for the binomial and Hadamard products, it is sufficient to express e_{k}(α_{1} + β_{1}, …, α_{i} + β_{j}, …, α_{m} + β_{n}) and e_{k}(α_{1}β_{1}, …, α_{i}β_{j}, …, α_{m}β_{n}) in terms of the e_{k}(α) and e_{k}(β), where for any symmetric polynomial f, f(α) means f(α_{1}, …, α_{m}) and similarly for f(β).

To do this we use the power sum symmetric polynomials p_{j}(u_{1}, …, u_{n}) = u_{1}^{j} + ⋯ + u_{n}^{j}. They are related to the elementary symmetric polynomials by the generating function formula

∑_{k=0}^{∞} e_{k}z^{k} = exp(∑_{j=1}^{∞}(−1)^{j−1} \frac{p_{j}}{j} z^{j}),  (30)

where we are omitting the arguments to e_{k} and p_{j}, which upon differentiating with respect to z yields Newton’s recurrence

ke_{k} = ∑_{i=1}^{k}(−1)^{i−1}e_{k−i}p_{i},  (31)

from which the elementary and power sum symmetric polynomials can easily be determined from each other. (With a computer algebra system it may be easier, though less efficient, to use (30) directly.) We will also need the value p_{0}(u_{1}, …, u_{n}) = n.

Let us write p_{k}(α ⊙ β) for

p_{k}(α_{1} + β_{1}, …, α_{i} + β_{j}, …, α_{m} + β_{n})

and p_{k}(α ∗ β) for

p_{k}(α_{1}β_{1}, …, α_{i}β_{j}, …, α_{m}β_{n}).

It is easy to see that

p_{k}(α ∗ β) = p_{k}(α)p_{k}(β),  (32)

and for p_{k}(α ⊙ β) we have

p_{k}(α ⊙ β) = ∑_{i,j}(α_{i} + β_{j})^{k} = ∑_{i,j} ∑_{l=0}^{k} \binom{k}{l}α_{i}^{l}β_{j}^{k−l}

= ∑_{l=0}^{k} \binom{k}{l}p_{l}(α)p_{k−l}(β).  (33)

We can then compute p_{k}(α) and p_{k}(β) for k from 1 to mn from the coefficients of U(x) and V (x) using (31), then compute p_{k}(α ⊙ β) and p_{k}(α ∗ β) from (32) and (33), and finally compute e_{k}(α ⊙ β) and e_{k}(α ∗ β) from (31).

As an example, we consider again the Fibonacci and Pell sequences; A(x) = x/(1−x−x^{2}) and B(x) = x/(1−2x−x^{2}). Here we have e_{1}(α) = 1, e_{2}(α) = −1, e_{1}(β) = 2, and e_{2}(β) = −1. Using (30) or (31) we compute p_{n}(α), which are the Lucas numbers, and p_{n}(β), which are the companion Pell numbers (\underline{A002203}), for n from 0 to 4. We then compute p_{n}(α ∗ β) and p_{n}(α ⊙ β) from (32) and (33), and finally we compute e_{n}(α ∗ β) and e_{n}(α ⊙ β) from (30) and (31). The values that occur in this computation are listed in Table 1.

Thus the denominator of A(x) ∗ B(x) is 1 − 2x − 7x^{2} − 2x^{3} + x^{4} and the denominator of A(x) ⊙ B(x) is 1 − 6x + 7x^{2} + 6x^{3} − 9x^{4}.

We note that the generating function A(x) ∗ B(x) = (x − x^{3})/(1 − 2x − 7x^{2} − 2x^{3} + x^{4}) can be found in the OEIS [15] at \underline{A001582}. See also Mező [23, Remark 9].

\begin{matrix}n & 0 & 1 & 2 & 3 & 4 \\ e_{n}(α) & 1 & 1 & −1 & 0 & 0 \\ e_{n}(β) & 1 & 2 & −1 & 0 & 0 \\ p_{n}(α) & 2 & 1 & 3 & 4 & 7 \\ p_{n}(β) & 2 & 2 & 6 & 14 & 34 \\ p_{n}(α ∗ β) & 4 & 2 & 18 & 56 & 238 \\ p_{n}(α ⊙ β) & 4 & 6 & 22 & 72 & 278 \\ e_{n}(α ∗ β) & 1 & 2 & −7 & 2 & 1 \\ e_{n}(α ⊙ β) & 1 & 6 & 7 & −6 & −9 \end{matrix}

Table 1: Fibonacci and Pell sequence computation

8 Partial fractions

Here we describe briefly another method, using partial fractions, for computing Hadamard and binomial products of rational power series. This method was used by Xin [29, Section 1-2] for Hadamard products.

We start with a simple example, the Hadamard product of the generating functions for Fibonacci and Pell numbers. Let A(x) = x/(1 − x − x^{2}) = ∑_{n=0}^{∞} F_{n}x^{n} and B(x) = x/(1 − 2x − x^{2}) = ∑_{n=0}^{∞} P_{n}x^{n}. Then

f = A(t)B(x/t) = ∑_{m,n=0}^{∞} F_{m}P_{n}x^{n}t^{m−n}.  (34)

With series like (34) involving infinitely many negative powers of variables, we must be careful about the power series ring in which we are working. Here we are in the ring ℚ((t))[[x]] of formal power series in x with coefficients that are Laurent series in t. So although elements of this ring may have arbitrary negative powers of t, the coefficient of any power of x involves only finitely many negative powers of t.

The constant term in t in (34) is

∑_{n=0}^{∞} F_{n}P_{n}x^{n} = A(x) ∗ B(x).

To compute it, we start with the partial fraction expansion of A(t)B(x/t) as a rational function of t,

\frac{1}{1 − t − t^{2}} \frac{x − x^{3} + (x^{3} + 2x^{2}) t}{1 − 2x − 7x^{2} − 2x^{3} + x^{4}} + \frac{1}{t^{2} − 2xt − x^{2}} \frac{x^{3} − x^{5} + (x^{3} + 2x^{2}) t}{1 − 2x − 7x^{2} − 2x^{3} + x^{4}}

= \frac{1}{1 − 2x − 7x^{2} − 2x^{3} + x^{4}} (R + S), (35)

where

R = \frac{x − x^{3} + (x^{3} + 2x^{2}) t}{1 − t − t^{2}}

and

S = \frac{(x^{3} − x^{5})t^{−2} + (2x^{2} + x^{3})t^{−1}}{1 − 2xt^{−1} − x^{2}t^{−2}} .

We want the coefficient of t^{0} in (35). Note that as elements of ℚ((t))[[x]], R contains only nonnegative powers of t and S contains only negative powers of t. Thus the constant term in t of A(t)B(x/t) may be obtained by deleting S from (35) and then setting t = 0, yielding A(x) ∗ B(x) = (x − x^{3})/(1 − 2x − 7x^{2} − 2x^{3} + x^{4}).

We can look at the partial fraction approach in a slightly different way, which shows how it is related to resultants. With A(x) and B(x) as above, to find the partial fraction expansion in t of A(t)B(x/t), we need to find two polynomials L(t) and M(t) in t, each of degree at most 1, with coefficients that are rational functions of x, such that

\frac{L(t)}{1 − t − t^{2}} + \frac{M(t)}{t^{2} − 2xt − x^{2}} = A(t)B(x/t) = \frac{xt^{2}}{(1 − t − t^{2})(t^{2} − 2xt − x^{2})}.

Equivalently,

L(t)(−x^{2} − 2xt + t^{2}) + M(t)(1 − t − t^{2}) = xt^{2}.  (36)

If we set L(t) = a + bt and M(t) = c + dt, then equating coefficients of powers of t in (36) gives the equivalent system

[\begin{matrix}−x^{2} & −2x & 1 & 0 \\ 0 & −x^{2} & −2x & 1 \\ 1 & −1 & −1 & 0 \\ 0 & 1 & −1 & −1 \end{matrix}] [\begin{matrix}a \\ b \\ c \\ d \end{matrix}] = [\begin{matrix}0 \\ 0 \\ x \\ 0 \end{matrix}]  (37)

Note that this matrix is essentially the Sylvester matrix for −x^{2} − 2xt + t^{2} and 1 − t − t^{2}, as polynomials in t. The constant term in t of A(t)B(x/t) is obtained by setting t = 0 in L(t) = a + bt so we may obtain this coefficient by solving for a in (37). By Cramer’s rule, the denominator of a, as a rational function of x, will be the determinant of the matrix in (37), so we see that this approach to Hadamard products is related to that of Section 6.

Our next theorem shows that the partial fraction method works in general for Hadamard products. For simplicity we assume that our rational power series are proper.

Theorem 15. Let A(x) and B(x) be proper rational power series with coefficients in a field F. Then A(t)B(x/t), as a rational function of t, has a partial fraction expansion of the form U(x, t) + V (x, t) in which, as elements of the Laurent series field F ((t))((x)), U(x, t) has only nonnegative powers of t and V (x, t) has only negative powers of t. Thus A(x) ∗ B(x) is obtained by setting t = 0 in U(x, t).

Proof. Let A(x) = N_{A}(x)/D_{A}(x) and B(x) = N_{B}(x)/D_{B}(x), where N_{A}(x), D_{A}(x), N_{B}(x), and D_{B}(x) are polynomials and D_{B}(x) has degree n. Then

A(t)B(x/t) = \frac{N_{A}(t) · t^{n}N_{B}(x/t)}{D_{A}(t) · t^{n}D_{B}(x/t)}  (38)

We first note that D_{A}(t) and t^{n}D_{B}(x/t) are relatively prime as polynomials in t. One way to see this is to factor D_{A}(t) and t^{n}D_{B}(x/t) into linear factors over some extension field G of F. Then the linear factors of D_{A}(t) are all of the form t − α for α ∈ G while the linear factors of t^{n}D_{B}(x/t) are all of the form t − βx for β ∈ G.

Since A(x) and B(x) are proper, it follows easily from (38) that A(t)B(x/t), as a rational function of t, is proper. Thus A(t)B(x/t) has a partial fraction expansion of the form

\frac{R_{A}(t, x)}{D_{A}(t)} + \frac{R_{B}(t, x)}{t^{n}D_{B}(x/t)},

where R_{A}(t, x) and R_{B}(t, x) are polynomials in t with coefficients that are rational functions of x, and the degree of R_{B}(t, x), as a polynomial in t, is less than n. Then t^{−n}R_{B}(t, x) is a polynomial in t^{−1} with no constant term, with coefficients that are rational functions of x, and 1/D_{B}(x/t), as an element of the Laurent series ring F ((t))((x)), has no positive powers of t. Thus t^{−n}R_{B}(t, x)/D_{B}(x/t) has only negative powers of t.

We know that R_{A}(t, x) is a polynomial in t, so R_{A}(t, x)/D_{A}(t) is a power series in t. Thus R_{A}(t, x)/D_{A}(t) is the sum of all terms in A(t)B(x/t) with nonnegative powers of t, and so the constant term in t in A(t)/B(x/t) may be obtained by setting t = 0 in R_{A}(t, x)/D_{A}(t).

While in principle the partial fraction method may be no more efficient than other methods, it may be easier in practice when using a computer algebra system in which partial fraction expansion is already built in. If the denominator polynomials D_{A}(x) and D_{B}(x) are not irreducible, the built-in partial fraction function may give us a more refined decomposition than we need, but we can still apply the same method to this decomposition, or alternatively, we may use the extended Euclidean algorithm, or solve a system of equations by whatever method is convenient, to find the two-term partial fraction expansion.

We can compute binomial products in a similar way, using Prodinger’s observation, discussed in Section 1.1, that the binomial product of A(x) and B(x) is the diagonal of (5). It follows that A(x) ⊙ B(x) is the constant term in t in

\frac{1}{1 − t}A (\frac{x}{1 − t}) B (\frac{x}{t})  (39)

which can also be computed by partial fraction expansion. (There is an analogue of Theorem 15 for binomial products.)

For example, if we take A(x) to be the Fibonacci generating function x/(1 − x − x^{2}) and B(x) to be the Pell generating function x/(1 − x − 2x^{2}) then expanding (39) by partial fractions in t gives

\frac{1}{1 − t}A (\frac{x}{1 − t}) B (\frac{x}{t}) = \frac{x^{2}}{1 − 6x + 7x^{2} + 6x^{3} − 9x^{4}} (R(x, t) + S(x, t)),

where

R(x, t) = \frac{2 − 5x + x^{2} + 3x^{3} + (x − 1)t}{1 − x − x^{2} + (x − 2)t + t^{2}}

and

S(x, t) = \frac{(1 − x)t + 2x^{2} − 3x^{3}}{t^{2} − 2xt − x^{2}} = \frac{(1 − x)t^{−1} + (2x^{2} − 3x^{3})t^{−2}}{1 − 2xt^{−1} − x^{2}t^{−2}} .

Thus

A(x) ⊙ B(x) = \frac{x^{2}}{1 − 6x + 7x^{2} + 6x^{3} − 9x^{4}} R(x, 0) = \frac{2x^{2} − 3x^{3}}{1 − 6x + 7x^{2} + 6x^{3} − 9x^{4}} ,

as we found by a different method in Section 4.2.

9 Acknowledgment

We would like to thank Alin Bostan for helpful suggestions and for bringing references [1, 4, 5] to our attention.

References

[1] Gessica Alecci, Stefano Barbero, and Nadir Murru, Some notes on the algebraic structure of linear recurrent sequences. arXiv:2302.13867 [math.NT], 2023. To appear, Ricerche di Matematica.

[2] W. N. Bailey, Generalized Hypergeometric Series, Cambridge Tracts in Mathematics and Mathematical Physics, No. 32, Stechert-Hafner, 1964.

[3] Roland Bacher, On exponentials of exponential generating series, Algebra Number Theory 4 (2010), 919–942.

[4] Alin Bostan, Philippe Flajolet, Bruno Salvy, and Éric Schost, Fast computation of special resultants, J. Symbolic Comput. 41 (2006), 1–29.

[5] L. Cerlienco, M. Mignotte, and F. Piras, Suites récurrentes linéaires: propriétés algébriques et arithmétiques, Enseign. Math. 33 (1987), 67–108.

[6] C. A. Church and Marjorie Bicknell, Exponential generating functions for Fibonacci identities, Fibonacci Quart. 11 (1973), 275–281.

[7] Roberto Dvornicich and Carlo Traverso, Newton symmetric functions and the arithmetic of algebraically closed fields. In Applied Algebra, Algebraic Algorithms and Error-Correcting Codes (Menorca, 1987 ), Vol. 356 of Lecture Notes in Comput. Sci., pp. 216–224. Springer, 1989.

[8] Shalosh B. Ekhad and Doron Zeilberger, Automatic generation of convolution identities for C-finite sequences. Arxiv preprint arXiv:2108.02918 [math.CO], 2021.

[9] Michel Fliess, Sur divers produits de séries formelles, Bull. Soc. Math. France 102 (1974), 181–191.

[10] R. Frontczak, T. Goy, and M. Shattuck, Identities relating six members of the Fibonacci family of sequences, Carpathian Math. Publ. 14 (2022), 6–19.

[11] Robert Frontczak, Taras Goy, and Mark Shattuck, Fibonacci-Lucas-Pell-Jacobsthal relations, Ann. Math. Inform. 55 (2022), 28–48.

[12] I. M. Gel’fand, M. M. Kapranov, and A. V. Zelevinsky, Discriminants, Resultants, and Multidimensional Determinants, Birkhäuser, 1994.

[13] M. L. J. Hautus and D. A. Klarner, The diagonal of a double power series, Duke Math. J. 38 (1971), 229–235.

[14] Verner E. Hoggatt, Jr. and Janet Crump Anaya, A primer for the Fibonacci numbers. XI. Multisection generating functions for the columns of Pascal’s triangle, Fibonacci Quart. 11 (1973), 85–90, 104.

[15] OEIS Foundation Inc., The On-Line Encyclopedia of Integer Sequences, Published electronically at https://oeis.org, 2023.

[16] Ishan Kar, A new method to compute the Hadamard product of two rational functions, Rose-Hulman Undergraduate Mathematics Journal 23 (2022), 12 pp.

[17] Jong Hyun Kim, Hadamard products and tilings, J. Integer Sequences 12 (2009), Article 09.7.4.

[18] Jong Hyun Kim, Hadamard Products, Lattice Paths, and Skew Tableaux, Ph. D. thesis, Brandeis University, 2011.

[19] Martin Klazar and Richard Horský, Are the Catalan numbers a linear recurrence sequence?, Amer. Math. Monthly 129 (2022), 166–171.

[20] Takao Komatsu, Convolution identities for tribonacci numbers, Ars Combin. 136 (2018), 199–210.

[21] Takao Komatsu, Convolution identities for Tribonacci-type numbers with arbitrary initial values, Palest. J. Math. 8 (2019), 413–417.

[22] Takao Komatsu and Rusen Li, Convolution identities for tribonacci numbers with symmetric formulae, Math. Rep. (Bucur.) 21 (71) (2019), 27–47.

[23] István Mező, Several generating functions for second-order recurrence sequences, J. Integer Sequences 12 (2009), Article 09.3.7.

[24] Elena A. Potekhina, Application of Hadamard product to some combinatorial and probabilistic problems, Discrete Math. Appl. 27 (2017), 177–186.

[25] Elena A. Potekhina and Mikhail I. Tolovikov, Oscillating random walk and the Hadamard product of rational functions, Discrete Math. Appl. 24 (2014), 29–44.

[26] Helmut Prodinger, Convolution identities for tribonacci numbers via the diagonal of a bivariate generating function, Palest. J. Math. 10 (2021), 440–442.

[27] Louis W. Shapiro, A combinatorial proof of a Chebyshev polynomial identity, Discrete Math. 34 (1981), 203–206.

[28] Richard P. Stanley, Enumerative Combinatorics. Vol. I, Wadsworth & Brooks/Cole, 1986.

[29] Guoce Xin, The Ring of Malcev–Neumann Series and the Residue Theorem, Ph. D. Thesis, Brandeis University, 2004.

2020 Mathematics Subject Classification: Primary 05A19, Secondary 05A15. Keywords: binomial convolution, Hurwitz product, Hadamard product, resultant, rational power series.

(Concerned with sequences \underline{A000045}, \underline{A000073}, \underline{A000129}, \underline{A001045}, \underline{A001582}, \underline{A001608}, \underline{A002203}, \underline{A002203}, and \underline{A014551}.)

Received April 24 2023; revised version received December 17 2023. Published in Journal of Integer Sequences, January 4 2024.

Return to Journal of Integer Sequences home page.