指数生成函数

引入

指数生成函数（exponential generating function，EGF）是解决排列类计数问题的有力工具。序列 \(a\) 的指数生成函数定义为形式幂级数： \(\hat{F}(x) = \sum_{n=0}^{\infty} a_n \frac{x^n}{n!}\)。

前置知识

多重集排列数

假设共有 \(n\) 个球，其中有 \(k\) 种不同种类的球，第 \(i\) 类球有 \(n_i\) 个，则这个 \(n\) 个球的全排列为： \(\frac{n!}{n_1!n_2!\cdots n_k!}\)。

考虑如何证明：

若 \(n\) 个球都不相同，那么得到一种排列之后，任意交换两个球的位置，都会得到一个新的排列。
但是在上述背景中，如果交换到了相同类别的球，将得到原来的排列。每种类别内部的排列可以知道一种方案被重复了 \(n_1!n_2!\cdots n_k!\) 次，故最终拿 \(n!\) 除以 \(n_1!n_2!\cdots n_k!\) 即是答案。

普通生成函数

指数型生成函数的意义在于解决排列类计数问题，而普通型生成函数旨在解决组合类计数问题。

考虑普通生成函数，在组合类计数问题中，在多个集合中分别选取若干个物品，每个集合都有一个选取方案数 \(\sum a_i x^i\)，其中 \(a_i\) 表示在这个集合中选取 \(i\) 个物品的方案。然后多个GF相乘起来得到了最终的GF，假设需要求解的物品总数为 \(m\) 的方案，那么 \([x^m]F(x)\) 也就是答案。

而排列类计数，难免需要涉及到多重集排列数，也就是我们需要考虑到去重这个问题。若从集合 \(m\) 中选出了 \(b_i\) 个物品，如果物品有不同，多重集排列数公式中的分母是不变化的，我们根本无法区分每个集合选取的数量，也就不知道对于每一类集合要乘的数是多少。

既然最后计算没有办法，那么我们换「提前计算」，这种思路在编程题计算中非常常见，我们习惯称为「费用提前打」。我们选取物品在排列，就算其他集合选取的数量多少，并没有关系。

故我们在原本的OGF（普通生成函数）的基础上，将单项式 \(a_i x^i\) 改为 \(a_i \frac{x^i}{i!}\)，也就得到了EGF（指数生成函数）。

值得注意的是，最终我们得到的总的GF中 \(x^n\) 前面的系数并不是直接的答案，我们计算出来的只有分母 \(b_1!b_2!\cdots b_k!\)，最终分子的 \(n!\) 还需要我们自己乘上去。

写成 \(F(x)\) 是将其看成OGF，此时 \([x^n]F(x) = \frac{a_n}{n!}\)，也就是我们求解出来的值。
\(\hat{F}(x)\) 是将其看成EGF，\([x^n]\hat{F}(x) = a_n\)，也就是将 \(\frac{x^n}{n!}\) 看成一个整体了。故最终的答案为 \([x^m]\hat{F}(x) = m! \times [x^m]F(x)\)，如果我们代码中写为多项式展开，求解出来的 \([x^m]F(x)\) 此时需要记得乘上 \(m!\)。

取球

假设有 \(3\) 个箱子，每个箱子里的球种类不同，箱子内的球种类相同，三个箱子分别含有 \(4\)、 \(3\)、 \(6\) 个球。求解拿出 \(1 \sim 8\) 个球出来排列分别有多少种方案。

箱子 \(1\) 的GF为： \((1 + \frac{x}{1!} + \frac{x^2}{2!} + \frac{x^3}{3!} + \frac{x^4}{4!})\)。
箱子 \(2\) 的GF为： \((1 + \frac{x}{1!} + \frac{x^2}{2!} + \frac{x^3}{3!})\)。
箱子 \(3\) 的GF为： \((1 + \frac{x}{1!} + \frac{x^2}{2!} + \frac{x^3}{3!} + \frac{x^4}{4!} + \frac{x^5}{5!} + \frac{x^6}{6!})\)。

最后相乘可以得到总的GF：

\[ \begin{aligned} \hat{F}(x) &= \left(1 + \frac{x}{1!} + \frac{x^2}{2!} + \frac{x^3}{3!} + \frac{x^4}{4!}\right)\left(1 + \frac{x}{1!} + \frac{x^2}{2!} + \frac{x^3}{3!}\right)\left(1 + \frac{x}{1!} + \frac{x^2}{2!} + \frac{x^3}{3!} + \frac{x^4}{4!} + \frac{x^5}{5!} + \frac{x^6}{6!}\right)\\ &= 1 + \frac{3}{1}x + \frac{9}{2}x^2 + \frac{9}{2}x^3 + \frac{10}{3}x^4 + \frac{77}{40}x^5 + \frac{643}{720}x^6 + \frac{1}{180}x^7 + \frac{61}{5760}x^8 + \frac{5}{21168}x^9 + \frac{17}{8640}x^{10} + \frac{13}{17280}x^{11} + \frac{1}{103680}x^{12} + \frac{1}{103680}x^{13} \end{aligned} \]

故，可以从最终的展开式中，看出：

拿 \(2\) 个球出来排列的方案数为 \(\frac{9}{2} \times 2! = 9\) 种。
拿 \(3\) 个球出来排列的方案数为 \(\frac{9}{2} \times 3! = 27\) 种。
拿 \(4\) 个球出来排列的方案数为 \(\frac{10}{3} \times 4! = 80\) 种。
拿 \(7\) 个球出来排列的方案数为 \(\frac{1}{180} \times 7! = 196\) 种。
拿 \(11\) 个球出来排列的方案数为 \(\frac{13}{17280} \times 11! = 39270\) 种。

但是将 \(F(x)\) 写成类似OGF的表达形式，不能直观的看出结果，所以我们一般将 \(\frac{a_i}{i!}\) 分开，将其中的 \(\frac{1}{i!}\) 提取出来，将其转化为EGF的形式：

\[ \begin{aligned} \hat{F}(x) &= 1 + \frac{3}{1!}x + \frac{9}{2!}x^2 + \frac{27}{3!}x^3 + \frac{80}{4!}x^4 + \frac{231}{5!}x^5 + \frac{643}{6!}x^6 + \frac{708}{7!}x^7 + \frac{4270}{8!}x^8 + \frac{3912}{9!}x^9 + \frac{21000}{10!}x^{10} + \frac{39270}{11!}x^{11} + \frac{6006}{12!}x^{12} + \frac{7007}{13!}x^{13} \end{aligned} \]

这样从分子上，便可以更加清晰看到方案。

基本运算

对于两个EGF的加法、减法运算，因为其本身还是一个多项式，直接相同系数相加即可。对于 \(F(x),G(x)\)，所以有：

\[ \begin{aligned} \hat{F}(x)\hat{G}(x) &= \sum_{i=0}^{\infty} a_i \frac{x^i}{i!} \sum_{k=0}^{\infty} b_k \frac{x^k}{k!}\\ &= \sum_{i=0}^{\infty} x^i \sum_{j=0}^{i} a_j b_{i - j} \frac{1}{j!(i - j)!}\\ &= \sum_{i=0}^{\infty} x^i \sum_{j=0}^{i} a_j b_{i - j} \frac{1}{j!(i - j)!} \cdot \frac{i!}{i!}\\ &= \sum_{i=0}^{\infty} \frac{x^i}{i!} \sum_{j=0}^{i} \dbinom{i}{j} a_j b_{i - j} \end{aligned} \]

故可以得到 \(\hat{F}(x)\hat{G}(x)\) 是序列 \((\sum_{j=0}^{i} \dbinom{i}{j} a_j b_{i - j})\) 的指数生成函数。

封闭式

类似OGF，EGF也有封闭形式，我们根据指数函数 \(e^x\) 的泰勒展开式，可以得到： \(e^x = 1 + x + \frac{x^2}{2!} + \frac{x^3}{3!} + \cdots = \sum_{i=0}^{\infty} \frac{x^i}{i!}\)。所以，这也就是这个生成函数，叫做指数生成函数的原因。经过了OGF的学习，我们不难知道对于求解某一类的单项 \([x^n]F(x)\)，当然EGF也是这一个思路，然后多个GF的封闭式相乘，最后展开，求其中的单项 \([x^n]F(x)\)。

我们依旧先介绍常见的序列的封闭式，然后再进行总结。

按照序号，下文中分别称为例式1、例式2、例式3……

从考虑序列 \(a = \langle 1,1,1,1,1,1,\cdots \rangle\)，其EGF为 \(\hat{F}(x) = \sum_{i=0}^{\infty} \frac{x^i}{i!}\)。

考虑上面 \(e^x\) 的泰勒展开式，不难得出其封闭式： \(\hat{F}(x) = e^x\)。
考虑序列 \(a = \langle 1,p,p^2,p^3,p^4,\cdots \rangle\)，即等比数列，其EGF为 \(\hat{F}(x) = \sum_{i=0}^{\infty} p^i \frac{x^i}{i!}\)。

易得其封闭式为 \(\hat{F}(x) = e^{px}\)。

考虑序列 \(a = \langle 1,-1,1,-1,1,\cdots \rangle\)，即 \(a_i = (-1)^i\)，其EGF为 \(\hat{F}(x) = \sum_{i=0}^{\infty} (-1)^i \frac{x^i}{i!} = e^{-x}\)。
考虑序列 \(a = \langle 0,1,0,1,0,1,\cdots \rangle\)，即 \(a_i = [2 \mid i]\)，其EGF为：

\(\hat{F}(x) = \sum_{i=0}^{\infty} [2 \mid i] \frac{x^i}{i!} = \sum_{i=0}^{\infty} \frac{x^{2i}}{(2i)!}\)。

可以得到 \(a_i = [2 \mid i] = \frac{1 + (-1)^i}{2}\)，也就是将例式1和例式3相加，然后除以2，便是答案。

当然，如果直接变换过于晦涩，可以考虑直接看序列本身，相比较较容易理解一些。

所以有 \(\hat{F}(x) = \frac{e^x + e^{-x}}{2}\)。

考虑序列 \(a = \langle 0,1,0,1,0,\cdots \rangle\)，即 \(a_i = [2 \nmid i]\)，其EGF为：

\(\hat{F}(x) = \sum_{i=0}^{\infty} [2 \nmid i] \frac{x^i}{i!} = \sum_{i=0}^{\infty} \frac{x^{2i + 1}}{(2i + 1)!}\)。

\(a_i = [2 \nmid i] = \frac{1 - (-1)^i}{2}\)，使用例式1和例式3的原始序列，将其发现其存在下标上两处需要为1，和数相反的情况，相同的为偶数下标，此处需要为0，相反的为奇数下标，此处需要为1。

所以有 \(\hat{F}(x) = \frac{e^x - e^{-x}}{2}\)。

考虑序列 \(a_i = A_n^n = n!\)，其中 \(n\) 为正整数常数，其EGF为 \(\hat{F}(x) = \sum_{i=0}^{\infty} i! \frac{x^i}{i!}\)。

易得其封闭式为 \(\hat{F}(x) = \frac{1}{1 - x}\)。

考虑 \(\ln(1 + x)\) 的泰勒展开式，\(\ln(1 + x) = x - \frac{1}{2}x^2 + \frac{1}{3}x^3 - \cdots\)。

其EGF为 \(\hat{F}(x) = \ln(1 + x) = \sum_{i=1}^{\infty} (-1)^{i - 1} \frac{x^i}{i} = \sum_{i=1}^{\infty} (-1)^{i - 1} (i - 1)! \frac{x^i}{i!}\)。

对应的序列为： \(a_i = (-1)^{i - 1} (i - 1)!, i > 0\)，特别地 \(a_0 = 0\)，即 \(\frac{a_i}{i!} = \frac{(-1)^{i - 1}}{i}\)。

考虑 \(-\ln(1 - x) = \sum_{i=1}^{\infty} \frac{x^i}{i}\)。

其EGF为 \(\hat{F}(x) = -\sum_{i=1}^{\infty} (-1)^{i - 1} \frac{(-x)^i}{i} = \sum_{i=1}^{\infty} \frac{x^i}{i}\)。

对应的序列为： \(a_i = (i - 1)!, i > 0\)，特别地 \(a_0 = 0\)，即 \(\frac{a_i}{i!} = \frac{1}{i}\)。

考虑 \(\sin x\) 的泰勒展开式，\(\sin x = x - \frac{1}{3!}x^3 + \frac{1}{5!}x^5 - \frac{1}{7!}x^7 + \cdots\)。

其EGF为： \(\hat{F}(x) = \sum_{i=0}^{\infty} (-1)^i \frac{x^{2i + 1}}{(2i + 1)!}\)。

对应的序列 \(a = \langle 0,1,0,-1,0,1,0,-1,\cdots \rangle\)，也就是以 \(0,1,0,-1\) 为循环节的序列。 - 考虑 \(\cos x\) 的泰勒展开式，\(\cos x = 1 - \frac{1}{2!}x^2 + \frac{1}{4!}x^4 - \frac{1}{6!}x^6 + \cdots\)。

其EGF为： \(\hat{F}(x) = \sum_{i=0}^{\infty} (-1)^i \frac{x^{2i}}{(2i)!}\)。

对应的序列 \(a = \langle 1,0,-1,0,1,0,-1,\cdots \rangle\)，也就是以 \(1,0,-1\) 为循环节的序列。

考虑 \(\arcsin x\) 的泰勒展开式，\(\arcsin x = x + \frac{1}{2} \cdot \frac{x^3}{3} + \frac{1 \times 3}{2 \times 4} \cdot \frac{x^5}{5} + \frac{1 \times 3 \times 5}{2 \times 4 \times 6} \cdot \frac{x^7}{7} + \cdots\)。

其EGF为：

\[ \begin{aligned} \hat{F}(x) &= x + \sum_{i=1}^{\infty} \frac{x^{2i + 1}}{(2i + 1)!} \left( \prod_{j=1}^{i} \frac{2j - 1}{2j} \right)\\ &= \sum_{i=0}^{\infty} \frac{(2i)!}{(i!)^2} \cdot \frac{1}{2i + 1} \cdot \frac{x^{2i + 1}}{(2i + 1)!}\\ &= \sum_{i=0}^{\infty} \frac{( (2i)! )^2}{4^i (i!)^2} \cdot \frac{x^{2i + 1}}{(2i + 1)!} \end{aligned} \]

对应的序列为：

\[ \begin{cases} a_i = \frac{(i - 1)!!}{2^{\frac{i - 1}{2}} \left( \frac{i - 1}{2} \right)!} & , i \mod 2 = 1 \\ a_i = 0 & , i \mod 2 = 0 \end{cases} \]

封闭式总结

此处对上一段「封闭式推导」进行了总结，最终归纳为如下一段，笔者的目的是为了方便自己为了避免冗余以及记忆负担，相同类型的式子仅选取了最终的拓展形式，和上文中的式子区分了，按照序号下称为式1、式2、式3、……、式9。

\(e^{px} = \sum_{i=0}^{\infty} p^i \frac{x^i}{i!}\)。
\(\frac{e^x + e^{-x}}{2} = \sum_{i=0}^{\infty} \frac{x^{2i}}{(2i)!}\)。
\(\frac{e^x - e^{-x}}{2} = \sum_{i=0}^{\infty} \frac{x^{2i + 1}}{(2i + 1)!}\)。
\((1 + x)^n = \sum_{i=0}^{n} \dbinom{n}{i} \frac{x^i}{i!} \cdot i! = \sum_{i=0}^{n} \dbinom{n}{i} x^i\)。
\(\ln(1 + x) = \sum_{i=1}^{\infty} (-1)^{i - 1} \frac{(i - 1)!}{i} \frac{x^i}{i!} = x - \frac{1}{2}x^2 + \frac{1}{3}x^3 - \cdots\)。
\(-\ln(1 - x) = \sum_{i=1}^{\infty} \frac{(i - 1)!}{i} \frac{x^i}{i!} = x + \frac{1}{2}x^2 + \frac{1}{3}x^3 + \cdots\)。
\(\sin x = \sum_{i=0}^{\infty} (-1)^i \frac{x^{2i + 1}}{(2i + 1)!} = x - \frac{1}{3!}x^3 + \frac{1}{5!}x^5 - \frac{1}{7!}x^7 + \cdots\)。
\(\cos x = \sum_{i=0}^{\infty} (-1)^i \frac{x^{2i}}{(2i)!} = 1 - \frac{1}{2!}x^2 + \frac{1}{4!}x^4 - \frac{1}{6!}x^6 + \cdots\)。
\(\arcsin x = \sum_{i=0}^{\infty} \frac{( (2i)! )^2}{4^i (i!)^2 (2i + 1)} \frac{x^{2i + 1}}{(2i + 1)!} = x + \frac{1}{2} \cdot \frac{x^3}{3} + \frac{1 \times 3}{2 \times 4} \cdot \frac{x^5}{5} + \frac{1 \times 3 \times 5}{2 \times 4 \times 6} \cdot \frac{x^7}{7} + \cdots\)。

指数生成函数与普通生成函数

如何理解指数生成函数？我们定义序列 \(a\) 的指数生成函数是：

\[ F(x)=\sum_{n\ge 0}a_n\frac{x^n}{n!} \]

但 \(F(x)\) 实际上也是序列 \(\left\langle \frac{a_n}{n!} \right\rangle\) 的普通生成函数。

这两种理解没有任何问题。也就是说，不同的生成函数只是对问题理解方式的转变。

EGF 中多项式 exp 的组合意义

前置知识：多项式 exp

如果您还没有学习多项式 exp 请先跳过这里，这是由 exp 理解引出的意义，从某种意义上来说可以加深对 EGF 的理解。

EGF 中 \(f^n(x)\) 的 \(f\) 默认是一个 EGF，那么我们首先考虑任意两个 EGF 的乘积

\[ \hat{H}(x) = \hat{F}(x)\hat{G}(x) = \sum_{n\geq 0} \left[ \sum_{i = 0}^n\binom {n}{i}f_ig_{n-i} \right] \frac{x^n}{n!} \]

对于两个 EGF 相乘得到的 \([x^k]\hat{H}(x)\)，实际上是一个卷积。而如果考虑多个 EGF 相乘得到的 \([x^k]\hat{H}(x)\)，实际上就是对每个 EGF 选择一项 \(x^{a_i}\) 使得 \(\sum_ia_i=k\) 时每种情况系数的和。从集合的角度来理解就是把 \(n\) 个有标号元素划分为 \(k>0\) 个有标号集合的方案数。

如果 \(k=0\) 则系数显然为原 EGF 各项常数的积，但是多项式 \(\exp\) 中某些要求导致 \(\exp\) 的 \(f(x)\) 常数项必须为 \(0\)，具体的原因在下文中会做出一些说明。

多项式系数定义（具体请参考排列组合一栏中的多项式系数组合意义）里是默认集合有序的，但是 \(\exp(f(x))\) 中 \(f^k(x)\)， \(k\) 个 \(f(x)\) 相乘得到的 EGF 相同，而集合划分显然是无序的，因此其系数应该乘上 \(\frac{1}{k!}\)。

设 \(F_k(n)\) 为 \(n\) 个有标号元素划分成 \(k\) 个非空无序集合（因为是 \(\exp\) 所以有非空要求）的情况，\(f_i\) 为 \(i\) 个元素组成一个集合时，\(i\) 个元素的有限集上特定组合结构的数量（是原有的 EGF，对这一个集合元素计数的方案，仅仅与该集合大小有关），那么 \(F_k(n)\)

\[ F_k(n)=\frac{n!}{k!}\sum_{\sum_{i}^ka_i=n}\prod_{j=1}^{k}\frac{f_{a_j}}{a_j!} \]

设 \(f_n\) 的 EGF 为 \(\hat{F}(x)\)，即：

\[ \hat{F}(x) = \sum_{n \geq 0} f_n\frac{x^n}{n!} \]

设 \(F_k(n)\) 的 EGF 为 \(G_k(x)\)，则：

\[ \begin{aligned} G_k(x)&=\sum_{n \geq 0} F_k(n)\frac{x^n}{n!}\\ &=\sum_{n \geq 0} x^n\frac{1}{k!}\sum_{\sum_i^k a_i=n}\prod_{j=1}^{k}\frac{f_{a_j}}{a_j!}\\ &=\frac{1}{k!}\sum_{n \geq 0}\sum_{\sum_i^k a_i=n}\prod_{j=1}^{k}\frac{f_{a_j}x^{a_j}}{a_j!}\\ &=\frac{1}{k!}\hat{F}^k(x) \end{aligned} \]

对于所有的 \(k \geq 0\) ：

\[ \sum_{k \geq 0}G_k(x) = \sum_{k \geq 0}\frac{\hat{F}^k(x)}{k!} = \exp{\hat{F}(x)} \]

上面是从组合角度直接列式理解，我们也可以从递推方面来证明 \(\exp(f(x))\) 和 \(f(x)\) 两者间的关系。

同样设 \(F_k(n)\) 为 \(n\) 个有标号元素划分成 \(k\) 个非空集合（无标号）的情况，\(g_i\) 为 \(i\) 个元素组成一个集合内部的方案数（意义同上文中的 \(f_i\)），并令 \(G(x)\) 为 \(\{g_i\}\) 的 EGF，\(H_k(x)\) 为 \(\{F_k(n)\}\) 的 EGF。

\(n\) 个元素中取出 \(i\) 个元素作为一个单独划分出去的集合共有 \(g_i\) 种方案，剩下的 \(n-i\) 个元素构成 \(k-1\) 个集合共 \(F_{k-1}(n-i)\) 种方案，但最后的划分方案中，一个方案里的每个集合都会被枚举为单独划分出去的集合，所以重复计算了 \(k\) 次，故还需要除以 \(k\)。

\[ \begin{aligned} H_k(x) &= \sum_{n\ge 0}\cfrac{x^n}{n!}F_k(n)\\ &=\sum_{n\ge 0}\cfrac{x^n}{n!}\sum_{i=1}^{n-k+1}\binom n {i} F_{k-1}(n-i)\times g_i\times \cfrac{1}{k}\\ &=\cfrac{1}{k}\sum_{n\ge 0}\cfrac{x^n}{n!}\sum_{i=0}^{n}\binom n {i}F_{k-1}(n-i)\times g_i\\ &=\cfrac{1}{k}\cdot H_{k-1}(x)G(x) \end{aligned} \]

上界是由非空集合划分推出的 \(n-(k-1)\geq i\)（前 \(k-1\) 个集合每个集合最少有一个元素），但是如果超过枚举上界涉及的 \(F_{k-1}(n-i)\) 设为 \(0\)，那么就没有影响。

得到递推式之后可递归展开，边界为 \(k=1\) 时 \(H_1(x)=G(x)\)。

\[ \begin{aligned} H_k(x) &= \cfrac{1}{k}\cdot H_{k-1}(x)G(x)\\ &= \cfrac{1}{k}\cdot\cfrac{1}{k-1}\cdot H_{k-2}(x)G^2(x)\\ &=\cdots \\ &= \cfrac{1}{k}\cdot\cfrac{1}{k-1} \cdots\cfrac{1}{2}\cdot H_{1}(x)G^{k-1}(x)\\ &= \cfrac{1}{k!}G^{k}(x)\\ \end{aligned} \]

同样的有：

\[ \begin{aligned} \sum_{k\ge 0}H_k(x)=\sum_{k\ge 0}\cfrac{G^k(x)}{k!}=\exp G(x) \end{aligned} \]

显然 定义成划分为非空集合（\(g_0=0\)）是符合本身的意义的，如果 包含空集（\(g_0=1\)），那么对应 \([x^n]G^k\) 中就会有 \([x^n]G^y,y>k\) 的贡献（在至少一个 \(G\) 中选择常数项），有计重，得不到所求量。

从递推式的角度讲，多个 EGF 的乘积也可以看作一个类似于背包的组合（合并两组计数对象的过程）。

总结多项式 \(\exp\) 的意义就是：有标号元素构成的集合的生成集族有多少种情况，或划分为任意个非空子集的总方案数。

例题

P5748 集合划分计数

一个有 \(n\) 个元素的集合，将其分为任意个非空子集，求方案数。

注意划分出的集合间是无序的，即 \(\{\{1,2\},\{3\}\}\) 和 \(\{\{3\},\{2,1\}\}\) 算作一种方案。

由于答案可能会很大，所以要对 \(998244353\) 取模。

【数据范围】 \(T = 1000\)， \(1\le n \le 10^5\)。

这个问题的答案被称为贝尔数 \(B_n\)。

首先，我们根据上面的例题知非空集合的EGF

\[\hat{F}(x) = e^x - 1\]

那么根据集合与划分的关系，就可以得到贝尔数的EGF

\[\hat{G}(x) = \exp(\hat{F}(x))\]

多项式EXP就行了。

实现

#include <bits/stdc++.h>
using namespace std;
#define endl '\n'
typedef long long ll;
const int MAXN = 4e6 + 10;
const ll p = 998244353, g = 3;
ll qpow(ll a, ll b)
{
    static ll res;
    for (res = 1; b; b >>= 1, a = a * a % p)
        if (b & 1) res = res * a % p;
    return res;
}
ll getinv(ll x) { return qpow(x, p - 2); }
const ll gi = getinv(g);
int rev[MAXN], n;
ll inv[MAXN], F[MAXN], G[MAXN], f[MAXN];
void initinv(int n)
{
    inv[1] = 1;
    for (int i = 2; i <= n; i++) inv[i] = (p - p / i) * inv[p % i] % p;
}
void Dert(int n, ll *F, ll *G)
{
    for (int i = 1; i < n; i++) G[i - 1] = F[i] * i % p;
}
void Int(int n, ll *F, ll *G)
{
    for (int i = 0; i < n; i++) G[i + 1] = F[i] * inv[i + 1] % p;
}
void NTT(ll *A, int n, int op)
{
    for (int i = 0; i < n; i++)
        if (i < rev[i]) swap(A[i], A[rev[i]]);
    for (int i = 2; i <= n; i <<= 1)
    {
        ll g1 = qpow(op == 1 ? g : gi, (p - 1) / i);
        for (int j = 0; j < n; j += i)
        {
            ll gk = 1;
            for (int k = j; k < j + i / 2; k++)
            {
                ll x = A[k], y = A[k + i / 2] * gk % p;
                A[k] = (x + y) % p;
                A[k + i / 2] = (x - y + p) % p;
                gk = gk * g1 % p;
            }
        }
    }
    if (op == 1) return;
    ll ni = getinv(n);
    for (int i = 0; i < n; i++) A[i] = A[i] * ni % p;
}
void Inv(int n, ll *F, ll *G)
{
    static ll A[MAXN], B[MAXN];
    G[0] = getinv(F[0]);
    int len, lim;
    for (len = 1; len < (n << 1); len <<= 1)
    {
        lim = len << 1;
        copy(F, F + len, A), copy(G, G + len, B);
        for (int i = 0; i < lim; i++)
            rev[i] = (rev[i >> 1] >> 1) | ((i & 1) * len);
        NTT(A, lim, 1), NTT(B, lim, 1);
        for (int i = 0; i < lim; i++)
            G[i] = ((2 - A[i] * B[i] % p) * B[i] % p + p) % p;
        NTT(G, lim, -1);
        fill(G + len, G + lim, 0);
    }
    fill(A, A + len, 0), fill(B, B + len, 0), fill(G + n, G + len, 0);
}
void Ln(int n, ll *F, ll *G)
{
    static ll A[MAXN], B[MAXN];
    Inv(n, F, A), Dert(n, F, B);
    int lim;
    for (lim = 1; lim < (n << 1); lim <<= 1);
    for (int i = 0; i < lim; i++)
        rev[i] = (rev[i >> 1] >> 1) | ((i & 1) * (lim >> 1));
    NTT(A, lim, 1), NTT(B, lim, 1);
    for (int i = 0; i < lim; i++) A[i] = A[i] * B[i] % p;
    NTT(A, lim, -1), fill(A + n, A + lim, 0), Int(n - 1, A, G),
        fill(A, A + lim, 0), fill(B, B + lim, 0);
}
void Exp(int n, ll *F, ll *G)
{
    static ll A[MAXN], B[MAXN];
    int lim, len;
    G[0] = 1;
    for (len = 1; len < (n << 1); len <<= 1)
    {
        lim = len << 1;
        copy(F, F + len, A), Ln(len, G, B);
        for (int i = 0; i < lim; i++)
            rev[i] = (rev[i >> 1] >> 1) | ((i & 1) * len);
        NTT(A, lim, 1), NTT(B, lim, 1), NTT(G, lim, 1);
        for (int i = 0; i < lim; i++)
            G[i] = (1 - B[i] + A[i] + p) % p * G[i] % p;
        NTT(G, lim, -1), fill(A, A + lim, 0), fill(B, B + lim, 0);
    }
    fill(G + n, G + len, 0);
}
void init()
{
    initinv(4e6);
    F[1] = 1;
    ll n = 1;
    while (n <= 1e5) n <<= 1;
    Exp(n, F, G);
    G[0] = (G[0] + p - 1) % p;
    memset(F, 0, sizeof(F));
    Exp(n, G, F);
    for (ll i = 1, frac = 1; i <= n; i++, frac = frac * i % p)
        f[i] = F[i] * frac % p;
}
int main()
{
    ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
    init();
    int T, x;
    cin >> T;
    while (T--) cin >> x, cout << f[x] << endl;
    return 0;
}

P4841 [集训队作业2013] 城市规划

刚刚解决完电力网络的问题，阿狸又被领导的任务给难住了。

刚才说过，阿狸的国家有 \(n\) 个城市，现在国家需要在某些城市对之间建立一些贸易路线，使得整个国家的任意两个城市都直接或间接的连通。

为了省钱, 每两个城市之间最多只能有一条直接的贸易路径。对于两个建立路线的方案，如果存在一个城市对，在两个方案中是否建立路线不一样，那么这两个方案就是不同的，否则就是相同的。现在你需要求出一共有多少不同的方案。

好了，这就是困扰阿狸的问题。换句话说，你需要求出 \(n\) 个点的简单 (无重边无自环) 有标号无向连通图数目。

由于这个数字可能非常大, 你只需要输出方案数对 \(1004535809\) ( \(479 \times 2 ^{21} + 1\) ) 取模即可。

对于 \(100\%\) 的数据，\(n \le 130000\)

我们设无向图的EGF为 \(\hat{F}(x)\)，由于任意两点间可以选择连边和不连边，因此有

\[\hat{F}(x) = \sum_{i \geq 0} 2^{\binom{i}{2}}\frac{x^i}{i!}\]

另一方面，设无向连通图的 EGF 为 \(\hat{G}(x)\)，我们可以认为，一个任意无向图可以由多个无向连通图组成，换句话说，一个任意无向图可以划分成多个无向连通图，因此有

\[\hat{F}(x) = \exp(\hat{G}(x))\]

就有 \(\hat{G}(x) = \ln\hat{F}(x)\)。

多项式ln就行了。

实现

#include <bits/stdc++.h>
using namespace std;
#define endl '\n'
typedef long long ll;
const int MAXN = 4e5 + 10;
const ll p = 1004535809, g = 3;
ll qpow(ll a, ll b)
{
    ll res = 1;
    for (; b; b >>= 1, a = a * a % p)
        if (b & 1) res = res * a % p;
    return res;
}
ll getinv(ll x) { return qpow(x, p - 2); }
const ll gi = getinv(g);
int rev[MAXN], n;
ll inv[MAXN], F[MAXN], G[MAXN], frac[MAXN], ifrac[MAXN];
void initinv(int n)
{
    inv[1] = frac[1] = frac[0] = 1;
    for (int i = 2; i <= n; i++)
        inv[i] = (p - p / i) * inv[p % i] % p, frac[i] = frac[i - 1] * i % p;
    ifrac[n] = getinv(frac[n]);
    for (int i = n - 1; i + 1; i--) ifrac[i] = (i + 1) * ifrac[i + 1] % p;
}
void Dert(int n, ll *F, ll *G)
{
    for (int i = 1; i < n; i++) G[i - 1] = F[i] * i % p;
}
void Int(int n, ll *F, ll *G)
{
    for (int i = 0; i < n; i++) G[i + 1] = F[i] * inv[i + 1] % p;
}
void NTT(ll *A, int n, int op)
{
    for (int i = 0; i < n; i++)
        if (i < rev[i]) swap(A[i], A[rev[i]]);
    for (int i = 2; i <= n; i <<= 1)
    {
        ll g1 = qpow(op == 1 ? g : gi, (p - 1) / i);
        for (int j = 0; j < n; j += i)
        {
            ll gk = 1;
            for (int k = j; k < j + i / 2; k++)
            {
                ll x = A[k], y = A[k + i / 2] * gk % p;
                A[k] = (x + y) % p;
                A[k + i / 2] = (x - y + p) % p;
                gk = gk * g1 % p;
            }
        }
    }
    if (op == 1) return;
    ll ni = getinv(n);
    for (int i = 0; i < n; i++) A[i] = A[i] * ni % p;
}
void Inv(int n, ll *F, ll *G)
{
    static ll A[MAXN], B[MAXN];
    G[0] = getinv(F[0]);
    int len, lim;
    for (len = 1; len < (n << 1); len <<= 1)
    {
        lim = len << 1;
        copy(F, F + len, A), copy(G, G + len, B);
        for (int i = 0; i < lim; i++)
            rev[i] = (rev[i >> 1] >> 1) | ((i & 1) * len);
        NTT(A, lim, 1), NTT(B, lim, 1);
        for (int i = 0; i < lim; i++)
            G[i] = ((2 - A[i] * B[i] % p) * B[i] % p + p) % p;
        NTT(G, lim, -1);
        fill(G + len, G + lim, 0);
    }
    fill(A, A + len, 0), fill(B, B + len, 0), fill(G + n, G + len, 0);
}
void Ln(int n, ll *F, ll *G)
{
    static ll A[MAXN], B[MAXN];
    Inv(n, F, A), Dert(n, F, B);
    int lim;
    for (lim = 1; lim < (n << 1); lim <<= 1);
    for (int i = 0; i < lim; i++)
        rev[i] = (rev[i >> 1] >> 1) | ((i & 1) * (lim >> 1));
    NTT(A, lim, 1), NTT(B, lim, 1);
    for (int i = 0; i < lim; i++) A[i] = A[i] * B[i] % p;
    NTT(A, lim, -1), fill(A + n, A + lim, 0), Int(n - 1, A, G),
        fill(A, A + lim, 0), fill(B, B + lim, 0);
}
int main()
{
    ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
    initinv(4e5);
    cin >> n;
    F[0] = 1;
    for (int i = 1; i <= n; i++)
        F[i] = qpow(2, 1ll * i * (i - 1) / 2) * ifrac[i] % p;
    Ln(n + 1, F, G);
    cout << G[n] * frac[n] % p;
    return 0;
}