Hi everyone!

In combinatorics, you often need to compute Stirling numbers for various reason. Stirling numbers of the first kind count the number of permutations of $$$n$$$ elements with $$$k$$$ disjoint cycles, while Stirling numbers of the second kind count the number of ways to partition a set of $$$n$$$ elements into $$$k$$$ nonempty subsets. Besides that, they're often arise when you need to change between powers of $$$x$$$ and rising/falling factorials. There are three problems on Library Checker that go as follows:

**Library Checker — Stirling Number of the First Kind**. Given $$$N$$$, find $$$s(N, k)$$$ for $$$0 \leq k \leq N$$$.

**Library Checker — Stirling Number of the Second Kind**. Given $$$N$$$, find $$$S(N, k)$$$ for $$$0 \leq k \leq N$$$.

**Library Checker — Stirling Number of the First Kind (Fixed K)**. Given $$$N$$$ and $$$K$$$, find $$$s(n, K)$$$ for $$$K \leq n \leq N$$$.

Here $$$s(n, k)$$$ are Stirling numbers of the first kind, and $$$S(n, k)$$$ are Stirling numbers of the second kind. Let's solve them!

### Definitions

There are various ways to define Stirling numbers. In this blog, we will follow the algebraic approach to highlight their connection.

**Def.**: The **falling factorial** is the polynomial $$$(x)_n$$$ defined as

**Def.**: The **rising factorial** is the polynomial $$$x^{(n)}$$$ defined as

**Note**: From the definitions above, it follows that $$$x^{(n)} = (x+n-1)_n = (-1)^n (-x)_n$$$.

We will start with Stirling numbers of the second kind, as the proof is simpler for them.

**Def.**: **Stirling numbers of the second kind** are the coefficients $$$S(n, k)$$$ in the expansion of $$$x^n$$$ into falling factorials:

**Claim**: The value $$$S(n, k)$$$ equates to $$$C(n, k)$$$ the number of ways to partition $$$n$$$ *distinct* elements into $$$k$$$ non-empty sets.

**Note**: We assume that there is a unique way to partition a set of $$$0$$$ elements into $$$0$$$ non-empty sets.

**Proof**

**Def.**: **Stirling numbers of the first kind** are the coefficients $$$s(n, k)$$$ in the expansion of $$$(x)_n$$$ into powers of $$$x$$$:

**Note**: From the definition above, and the fact that $$$x^{(n)} = (-1)^n (-x)_n$$$, it follows that

**Claim**: The absolute value $$$|s(n, k)|$$$ equates to $$$c(n, k)$$$, is the number of permutations of length $$$n$$$ with $$$k$$$ disjoint cycles.

**Proof**

### Generating functions

Finding Stirling numbers with fixed $$$n$$$ or $$$k$$$ will most certainly involve working with polynomials. As it's usually most convenient to explain polynomial operations in the framework of generating functions, let's find them for Stirling numbers.

#### Stirling numbers of the second kind

As we already mentioned, $$$S(n, k)$$$ is the number of ways to partition a set of $$$n$$$ distinct objects into $$$k$$$ non-empty sets. The exponential formula tells us that if $$$F(x)$$$ is the exponential generating function for objects of some kind, then $$$\exp F(x)$$$ is the generating function for sets of such objects, also meaning that $$$\exp F(x) - 1$$$ is the generating function for non-empty sets of such objects.

As things enumerated by $$$S(n, k)$$$ can be described as "sets of $$$k$$$ non-empty sets of objects", the exponential generating function should look roughly as $$$e^{e^{x} - 1}$$$, where $$$x$$$ is a formal variable keeping track of "base" objects. However, for Stirling numbers we need to keep track of the base objects (their total amount sums to $$$n$$$), *and* of the non-empty sets, in which we split them (their amount sums to $$$k$$$).

To account for non-empty sets, we will introduce another variable $$$y$$$, and use the generating function

instead. In this way, when the external exponent is expanded into powers of $$$(e^x-1)$$$, each summand $$$\frac{(e^x-1)^k}{k!}$$$, corresponding to sets of $$$k$$$ non-empty sets will be additionally multiplied by $$$y^k$$$. Thus, we get the expansion into generating function as

##### Fixed $$$k$$$

Using the power series definition of exponent, we can then expand it as

from which we get the exponential generating function with a fixed $$$k$$$:

##### Fixed $$$n$$$

**Library Checker — Stirling Number of the Second Kind**. Given $$$N$$$, find $$$S(N, k)$$$ for $$$0 \leq k \leq N$$$.

But what if we want to expand with fixed $$$n$$$? No worries! Let $$$F(y)$$$ be the ordinary generating function of $$$S(n, k)$$$ with fixed $$$n$$$.

From the definition of Stirling numbers, considering $$$t=0, 1, \dots$$$ we get

Right-hand side can be perceived as a convolution of $$$F(y)$$$ and $$$e^y$$$, as the later expands into sum of $$$\frac{y^{t-k}}{(t-k)!}$$$. Thus,

from which follows the generating function with a fixed $$$n$$$:

Defined this way, the generating function on the right-hand side is also known as the $$$n$$$-th Touchard polynomial.

Expanded, the expression above transforms into the well-known explicit formula

which can also be obtained through the principle of inclusion and exclusion. The formula above can also be used to compute a specific $$$S(n, k)$$$ in $$$O(k \log_k n)$$$ without Fourier transform or other complicated algorithms.

#### Stirling numbers of the first kind

Generally, a permutation can be represented as a set of cycles. The exponential generating function for cycles is known to be

Now, we can mark each cycle by an additional formal variable $$$y$$$ and apply exponential formula. However, we should note that for signed Stirling numbers, we should multiply the coefficient near $$$x^n y^k$$$ with $$$(-1)^{n-k}$$$, which is done by substituting $$$x$$$ and $$$y$$$ with $$$-x$$$ and $$$-y$$$:

So, similarly to Stirling numbers of the second kind, we get the expansion

Unlike generating functions for the Stirling numbers of the second kind, this one can be checked in a very straightforward way:

##### Fixed $$$k$$$

**Library Checker — Stirling Number of the First Kind (Fixed K)**. Given $$$N$$$ and $$$K$$$, find $$$s(n, K)$$$ for $$$K \leq n \leq N$$$.

One obvious way to find the generating function with fixed $$$k$$$ is to expand

so we get the exponential generating function with a fixed $$$k$$$:

##### Fixed $$$n$$$

**Library Checker — Stirling Number of the First Kind**. Given $$$N$$$, find $$$s(N, k)$$$ for $$$0 \leq k \leq N$$$.

**Sasha and the Swaps II**. Given $$$n$$$, find the number of distinct permutations obtainable with exactly $$$k$$$ swaps for each $$$k$$$ from $$$1$$$ to $$$n-1$$$.

For a fixed $$$n$$$, we already know that the generating function would be

Computing it as is with divide and conquer approach, we would get $$$O(n \log^2 n)$$$. Can we do better? Sure thing! We can represent it as

therefore we can first compute $$$(y)_n$$$, then do a Taylor shift to find $$$(y-n)_n$$$ and multiply them in $$$O(n \log n)$$$. The overall time is

wow thats great. I had thought that we could calculate Stirling number in O(n^2).

Nice technique! I have one question: can it be generalized for the arbitrary(?) exponential families?

Let's suppose we have an exponential family $$$F$$$ and we know it's deck enumerator $$$D(x)$$$, and we are asked: what is the number of hands of weight $$$n$$$ and $$$k$$$ cards (with fixed $$$n$$$)?

For example, in this problem we have to find the number of functional graphs with $$$n$$$ vertices and $$$k$$$ components for $$$1 <= k <= n$$$.

Surely, functional graphs can be thought as an exponential family, where cards are

connectedfunctional graphs, decks are standard cards (so, basically, $$$[\frac{x^n}{n!}]D(x)$$$ is a number of connected functional graphs with $$$n$$$ vertices), and hands are all functional graphs (not necessarily connected).We know that $$$H(x) = \sum\limits_{i=1}^\infty \frac{i^i}{i!}x^i$$$, and know that $$$D(x) = \ln H(x)$$$.

So, we are able to solve for fixed $$$k$$$ in $$$O(n log n)$$$, because $$$h(n, k) = [\frac{x^n}{n!}]{\frac{{D(x)^k}}{k!}}$$$, but how can we find $$$G(y) = \sum\limits_{k=1}^\infty h(n, k){y}^{k}$$$?

Everything i tried was not successful...

Uhhh, I'm not sure I understand what you mean here. What is an exponential family? What is a deck enumerator?

My bad — taking into account an enormous amount of your blogs about combinatorics and generating functions, I had no doubt you've read Herbert S Wilf's "gfology2" (These definitions are given in section "Cards, Decks and Hands: The exponential formula").

From exponential formula we know that if we EGF of some combinatorial objects, say $$$D(x)$$$, s.t. $$$[\frac{x^n}{n!}]{D(x)} = d_{n}$$$ is a number of these objects of weight $$$n$$$, then $$$H(x) = e^{D(x)}$$$ is the EGF of sets of these objects, and, in particular, number of sets of weight $$$n$$$ and $$$k$$$ objects is $$$h(n, k) = [\frac{x^n}{n!}]({\frac{{D(x)^k}}{k!}})$$$.

So, the initial question can be formed like that: if we have EGF of some combinatorial objects $$$D(x)$$$, can we somehow find $$$G(y) = \sum\limits_{k=1}^\infty h(n, k){y}^{k}$$$ (that is, find $$$h(n, k)$$$ for fixed $$$n$$$ and $$$1 <= k <= n$$$)?

In your blog you have done the particular case with cycles and permutations, my example in comment above is about connected functional graphs and arbitrary functional graphs (it is just another example of relation "block" — "set of blocks").

Ok, so, as I understand, your question is essentially as follows:

Let $$$F(x, y) = e^{y D(x)} = H(x)^y$$$. Then, $$$[y^k] F(x, y) = \frac{D(x)^k}{k!}$$$, but what about $$$G(y) = \left[\frac{x^n}{n!}\right] F(x, y)$$$?

That's indeed an interesting generalization, and the blog doesn't offer any direct solutions. And I don't know any solution to it yet, unfortunately. Given that the problem has $$$n \leq 5000$$$ constraint, it also likely expects some quadratic, DP-based solution. But I will think of how to do it.

Yeah, you're right — that is exactly what I wanted to ask! Coming back to the CSES problem: yes, intended solution is hardly related to generating functions, I've just tried the straightforward exponential formula approach and came to the question of existance of some general methods of computing $$$h(n, k)$$$ for fixed $$$n$$$ or fixed $$$k$$$, that is why I'm here :)

I asked this on Mathoverflow (here), and was told that the problem is essentially equivalent to computing a particular polynomial in a binomial type sequence.

I also noted that we can find the value of $$$G(k)$$$ in any particular $$$k$$$ point by computing $$$[x^n] H(x)^k$$$, from which we could interpolate it in $$$O(n \log^2 n)$$$.

Well, the problem with this is that, apparently, finding $$$[x^n] H(x)^k$$$ for $$$k=1..n$$$ is a transposed problem to polynomial composition, for which best known algorithms are somewhere in $$$O(n \sqrt n \log n)$$$ range, meaning that there might be no efficient general-case algorithm.

Thanks a lot for your efforts!

This problem is a great example of case when one should give the straightforward generating functions approach up — the computational part of the solution is neither simple, nor fast enough.