形式幂级数与生成函数
形式幂级数
设 $z$ 是未定元，我们定义一个形式幂级数是一个无穷级数：
$$ A(z)=a_0+a_1z^1+a_2z^2+…+a_nz^n+…=\sum\limits_{i\ge 0} a_iz^i $$
其中 $\langle a_0,a_1,a_2,…\rangle$ 是一个无穷的实系数序列。
在这里，我们并不关心 $z$ 是否代入了某个值，它如同一个“占位符”，并不参与我们的运算。
生成函数
对于任意一个序列 $\langle g_0,g_1,g_2,…,g_n\rangle$，我们可以认为其不存在的项全部等于 0，然后将其看作一个无穷的数列。这样，我们可以定义这个数列的普通生成函数（Ordinary Generating Function, OGF）是一个形式幂级数：
$$G(z)=\sum\limits_{i\ge 0} g_iz^i$$
可以定义这个数列的指数生成函数（Exponential Generating Funcion, EGF）是一个形式幂级数：
$$G(z)=\sum\limits_{i\ge 0} \dfrac{g_iz^i}{i!}$$
这样我们就通过统一化的语言来表达了一个数列，而这就是生成函数的意义所在：解决数列问题的工具。
从另一种角度来讲，生成函数也可以看作一个次数是正无穷的多项式，因此生成函数的乘法即是多项式相乘。
一般情况下，我们只关心某个生成函数的前 $n$ 项系数，记作 $A(z)\text{ }\text{mod}\text{ } z^n$。这个形式对应了一个有限项的系数序列。
卷积
乘法卷积
对于数列 $f,g$，我们定义它们的乘法卷积（或者直接称作卷积）是一个数列 $c$，满足：
$$ c_i=\sum\limits_{j=0}^i f_jg_{i-j}$$
同时记作 $c=fg$。容易发现，两个数列乘法卷积的 OGF 是它们各自的 OGF 相乘得到的结果。
二项卷积
对于数列 $f,g$，我们定义它们的二项卷积是一个数列 $c$，满足：
$$c_i=\sum\limits_{j=0}^i \dbinom{i}{j}f_jg_{i-j}$$
同时记作 $c=fg$，这里乘法的定义与上面的不同，但是惯用记号一样，需要根据语境来区分。容易发现，两个数列二项卷积的 EGF 是它们的 EGF 相乘得到的结果。
狄利克雷卷积
对于数列 $f,g$，我们定义它们的狄利克雷卷积是一个数列 $c$，满足：
$$c_i=\sum\limits_{k|i} f_{k}g_{i/k}$$
记号同上。狄利克雷卷积在数论中很常见。
类似的，我们可以定义异或卷积 $(fg)_i=\sum\limits_{j\oplus k=i} f_jg_k$，或卷积 $(fg)_i=\sum\limits_{j\cup k=i}f_jg_k$，与卷积 $(fg)_i=\sum\limits_{j\cap k=i} f_jg_k$，其中我们用 $\oplus$ 代表按位异或，用 $\cup$ 代表按位或，用 $\cap$ 代表按位与。
除去狄利克雷卷积可以直接 $O(n\ln n)$ 计算之外，其它卷积计算的朴素实现均是 $O(n^2)$ 的。因此我们需要通过下面的这些手段去计算卷积。
离散傅里叶变换
离散傅里叶变换（Discrete Fourier Transform, DFT）是对多项式进行从系数表示到点值表示的变换。设 $f(x)$ 是一个 $n-1$ 次多项式，那么它的 DFT 对应一个长度为 $n$ 的序列 $\mathcal F(f)$，满足：
$$ \mathcal F(f)_i=f(x_i) $$
其中 $x_i(i=0,1,2,…,n-1)$ 是代入的点值。
快速傅里叶变换（Fast Fourier Transform, FFT）通过 $n$ 次单位根 $\omega_n$ 来进行 DFT，其中 $\omega_n$ 满足 $\omega_n^n=1$；在复数域下，这样的数字一共有 $n$ 个。
如果令 $\omega_n^1=\cos \frac{2\pi}{n}+i\sin\frac{2\pi}{n}$，那么它们分别是 $\omega_n^0,\omega_n^1,\omega_n^2,…,\omega_n^{n-1}$。
由于 $n$ 次单位根的一些特殊性质，所以当 $n=2^k$ 时，我们可以通过分治来求解点值，以及通过点值来求解系数，后者称作逆变换（Inverse Discrete Fourier Transfrom, IDFT）。
这样就得到了一个在 $O(n\log n)$ 的时间内进行 DFT 的算法（FFT）。需要注意的一点是，如果用 FFT 计算序列卷积，得到的是长度为 $n$ 的循环卷积 $(fg)_i=\sum\limits_{j+k\text{ }\text{mod }n=i}f_jg_k$，用生成函数记作：
$$F(z)G(z)\text{ }\text{mod }z^n-1$$
FFT 的代码实现可以在模板梳理中找到。
快速数论变换
- 原根：在 $\mathbb F_R$ 下，如果一个数 $g^i$ 一共有 $\varphi(R)$ 种不同的取值，或者说 $g$ 模 $R$ 的阶是 $\varphi(R)$，那么称 $g$ 是 $\mathbb F_R$ 下的原根。
可以证明，原根存在的充要条件是 $R=2,4,p^a,2p^a$，其中 $p$ 是一个素数。可以发现，素数的原根是总存在的。
考虑做长度为 $n=2^k$ 的 DFT，有没有什么方式可以避免复数运算？
对于形如 $P=a2^r+1(r\ge k)$ 的素数，设 $g$ 是 $\mathbb F_P$ 下的原根，那么 $g_n=g^{\varphi(P)/n}$ 有着与 $n$ 次单位根 $\omega_n$ 相同的性质，即 $g_n^n=1$。根据原根的性质，我们可以发现 $g_n^0,g_n^1,…,g_n^{n-1}$ 是 $n$ 个两两不同的数字，且都满足 $(g_n^l)^n\equiv g^{l\varphi(P)}\equiv 1\text{ }(\text{mod }P)$。于是，我们可以用 $g_n$ 来代替 $\omega_n$ 进行多项式变换，这样就得到了快速数论变换，本质上是 FFT 在模域下的变种。
一个质数的原根一定存在，但并不是只要存在原根就可以 NTT，另一个条件是 $n|\varphi(P)$，形象的理解即可以等分圆周。
代码同样可以在模板梳理中找到。
值得一提的是，这种计算点值的方法与单位根反演有着异曲同工之妙。单位根反演即是这个柿子：
$$[k|n]=\dfrac{1}{k}\sum\limits_{i=0}^{k-1}\omega_k^{in}$$
其中 $\omega_k$ 是 $k$ 次单位根，模义下可用 $g_k$ 来代替；这可以看作是计算了一个点值，在一些题目中会具有优良的性质。
离散沃尔什变换
显然，位运算卷积不满足卷积结果的点值是原点值相乘。那么我们可以通过一些手段构造某种变换，来让我们变换出的序列满足点值相乘的规律。这些变换分别是：
- 按位或
$$\mathcal{F}(f)_i=\sum\limits_{j\cup i=i}a_j$$
- 按位与
$$\mathcal{F}(f)_i=\sum\limits_{j\cap i=i}a_j$$
- 按位异或
$$\mathcal{F}(f)_i=\sum\limits_{j=0}^{n-1} (-1)^{|i\cap j|} f_j$$
其中 $|x|$ 代表 $x$ 在二进制下为 1 的位数。可以证明，这些变换满足与 DFT 同样的性质，即点值的乘积是卷积结果的点值。
对于异或的变换，还有一个重要的柿子是 $\mathcal{F}(\mathcal{F}(f))_i=n·f_i$，可以通过代入化简来证明。
这些变换都可以通过类似 FFT 的方法在 $O(n\log n)$ 的时间内得到。代码可以参考模板梳理。
微积分
形式幂级数的微分定义为逐项微分，即：
$$\dfrac{\text d}{\text dz} \left(\sum\limits_{i\ge 0} a_iz^i \right)=\sum\limits_{i\ge 0} (i+1)a_{i+1}z^i$$
同样的，形式幂级数的积分定义为逐项积分。我们可以定义一个形式幂级数的定积分为：
$$\int \left(\sum\limits_{i\ge 0}a_iz^i \right)\text dz=\sum\limits_{i>0}\dfrac{a_{i-1}}{i}z^i$$
我们可以用微积分来分析一个幂级数的闭合形式，这在涉及指标变换的时候较为方便。
乘法逆
对于幂级数 $A(z)$，我们定义它的乘法逆是一个幂级数 $B(z)$ 满足 $A(z)B(z)\equiv 1 \text{ }(\text{mod } z^n)$，记作 $B(z)=A^{-1}(z)$。
当 $A(z)$ 的常数项不为零时，它的乘法逆是总存在的。证明如下：
设 $A(z)=\sum\limits_{i=0}^{n-1} a_ix^i, B(z)=\sum\limits_{i=0}^{n-1} b_ix^i$，通过乘法逆的定义可得：
$$\begin{align} &b_0=\dfrac{1}{a_0}\newline &\sum\limits_{j=1}^{n-1} a_jb_{i-j}=0\text{ }(i>0)\end{align}$$
对 $i$ 归纳可以得到：
$$ b_i=-b_0\sum\limits_{j=0}^{i-1}a_{i-j}b_j\text{ }(i>0) $$
从而只要 $b_0$ 存在就可以构造出 $B(z)$。同时这也给出了朴素 $O(n^2)$ 求逆的方法，下一步考虑怎么优化。
设 $B_0(z)=A^{-1}(z)\text{ }\text{mod } z^n$ 是一个只有前 $n$ 项系数有效的幂级数，可以如下倍增：
$$\begin{align} B_0(z)A(z)&\equiv 1\text{ }(\text{mod } z^n)\newline (B_0(z)A(z)-1)^2&\equiv 0\text{ }(\text{mod }z^{2n})\newline B_0^2(z)A^2(z)-2B_0(z)A(z)+1&\equiv 0\text{ }(\text{mod }z^{2n}) \newline B_0(z)(2-B_0(z)A(z))&\equiv B(z)\text{ }(\text{mod }z^{2n}) \end{align}$$
其中 $B(z)=A^{-1}(z)\text{ }\text{mod } z^{2n}$ 是一个只有前 $2n$ 项系数有效的幂级数。
依此倍增，复杂度 $T(n)=T(n/2)+O(n\log n)=O(n\log n)$，实现可以在代码梳理中找到。
指对运算
泰勒展开
对一个高阶可导的函数 $f(x)$，它在 $x=x_0$ 处的泰勒展开式为：
$$f(x)=f(x_0)+f’(x_0)(x-x_0)+\dfrac{f’’(x_0)(x-x_0)^2}{2}+…+\dfrac{f^{(n)}(x-x_0)^n}{n!}+…$$
这是一个无穷级数。也可以只保留前 $n+1$ 项写作：
$$f(x)=f(x_0)+f’(x_0)(x-x_0)+\dfrac{f’’(x_0)(x-x_0)}{2}+…+\dfrac{f^{(n)}(x_0)}{n!}+r_n(x)$$
其中 $r_n(x)$ 指代一个余项。
麦克劳林公式
在泰勒展开的公式中，取 $x_0=0$，就得到了麦克劳林公式：
$$f(x)=f(0)+f’(0)x+\dfrac{f’’(0)x^2}{2!}+…+\dfrac{f^{(n)}(0)x^n}{n!}+…$$
亦即 $f(x)$ 在 $x=0$ 处的泰勒展开。
形式幂级数的乘方是有定义的，那么我们可以通过泰勒展开来定义如何对形式幂级数进行各种诡异的函数变换，诸如指对函数，正余弦函数，等等。
指数函数和对数函数
由于 $(e^x)’=e^x$，根据麦克劳林公式，可以得到：
$$\exp x=e^x=1+x+\dfrac{x^2}{2!}+\dfrac{x^3}{3!}…$$
这是一个无穷级数。于是可以定义一个形式幂级数的指数函数为：
$$\exp A(z)=\sum\limits_{i\ge 0}\dfrac{A^i(z)}{i!}$$
同理，因为：
$$\ln’ x=x^{-1},\ln ''x=-x^{-2},…,,\ln^{(n)} x=(-1)^{n-1}(n-1)!(1+x)^{-n}$$
可以定义一个形式幂级数的对数函数是：
$$\ln (1+A(z))=-\sum\limits_{i>0}\dfrac{(-1)^iA^i(z)}{i}$$
注意这里的常数项是未定义的，一般认为其等于 0。
考虑如何求对数函数。设 $B(z)=\ln A(z)$，两边求导得到：
$$B’(z)=\ln’ A(z)=\dfrac{A’(z)}{A(z)}$$
即是链式法则。于是我们只需要求逆即可，复杂度 $O(n\log n)$。
多项式指数的计算需要用到牛顿迭代。
牛顿迭代法
考虑计算一个函数 $f(x)$ 的零点。设 $x_0$ 是一个已知的近似解，在 $x_0$ 处对 $f(x)$ 泰勒展开，可以得到 $f(x)$ 的一个近似表达：
$$ f(x)\sim f(x_0)+f’(x_0)(x-x_0) $$
令 $f(x_0)+f’(x_0)(x-x_0)=0$，解得 $x=x_0-\dfrac{f(x_0)}{f’(x_0)}$。取 $x$ 作为下一个近似解 $x_2$，又可以进行同样的操作。可以证明，在某些情况下，这样就可以逼近 $f(x)$ 精确的零点。
设 $B(z)=\exp A(z)$，那么有 $\ln B(z)-A(z)=0$。设 $G(B(z))=\ln B(z)-A(z)$ 是一个以 $B(z)$ 为变元的形式幂级数，$B_0(z)=\exp A(z)\text{ }\text{mod }z^n$ 是已知的近似解，那么可以如下倍增（迭代）：
$$\begin{align} B(z)&\equiv B_0(z)-\dfrac{\ln B_0(z)-A(z)}{1/B_0(z)}\text{ }(\text{mod } z^{2n})\newline &\equiv B_0(z)(1-\ln B_0(z)+A(z))\text{ }(\text{mod } z^{2n})\end{align}$$
依此迭代，复杂度 $T(n)=T(n/2)+O(n\log n)=O(n\log n)$，实现可以在代码梳理中找到。
应用
一道栗题
请构造一个序列 $a_1,a_2,…,a_k$，满足 $\forall i\in[1, k] a_i\in[1, n]$，使得对这个序列做无限背包计数后得到的序列的前 $n$ 项在模 $P$ 意义下是 $f_1,f_2,…,f_n$，其中 $P$ 是一个给定的质数。序列的长度没有限制，可以证明解总存在且唯一。
$n\le 2^{18},10^6\le p\le 2^{30}$
设 $b_i$ 表示序列中是否包含 $i$ 这个数字，$F(z)=\sum\limits_{i=1}^n f_iz^i$，那么有：
$$\begin{align} F(z)&=\prod\limits_{i=1}^n \dfrac{1}{(1-z^i)^{b_i}}\newline \ln F(z)&=-\sum\limits_{i=1}^nb_i\ln(1-z^i)\newline \ln F(z)&=\sum\limits_{i=1}^nb_i\sum\limits_{j>0}\dfrac{z^{ij}}{j}\newline &=\sum\limits_{i=1}^n z^i\sum\limits_{j|i}\dfrac{jb_j}{i} \end{align}$$
比较系数得 $nf_n=\sum\limits_{i|n}ib_i$，亦即 $f=b*1$，莫比乌斯反演得 $b_n=\dfrac{1}{n}\sum\limits_{j|n}jf_j\varphi(n/j)$。
复杂度 $O(n(\log n+\ln n))$，需要任意模数卷积，这里用拆系数 FFT 实现。
代码：
#include<bits/stdc++.h>
using namespace std;
#define LL long long
#define LDB long double
const int CN = (1 << 18) + 3;
const LDB PI = acos(-1);
int P;
int read(){
int s = 0, ne = 1; char c = getchar();
for(;c < '0' || c > '9';c = getchar()) if(c == '-') ne = -1;
for(;c >= '0' && c <= '9';c = getchar()) s = (s << 1) + (s << 3) + c - '0';
return s * ne;
}
int add(int x, int y) {return x + y >= P ? x + y - P : x + y;}
int qp(int a, int b){
int r = 1;
while(b){
if(b & 1) r = 1ll * r * a % P;
a = 1ll * a * a % P, b >>= 1;
}
return r;
}
int invx(int x) {return qp(x, P - 2);}
namespace POLY{
class COMP{
public: LDB x, y;
COMP operator + (const COMP &o) const{
COMP r = *this; r.x += o.x, r.y += o.y;
return r;
}
COMP operator - (const COMP &o) const{
COMP r = *this; r.x -= o.x, r.y -= o.y;
return r;
}
COMP operator * (const COMP &o) const{
COMP r; r.x = x * o.x - y * o.y, r.y = x * o.y + y * o.x;
return r;
}
} ;
COMP mk(LDB a, LDB b) {COMP o; o.x = a, o.y = b; return o;}
COMP conj(COMP o) {o.y = -o.y; return o;}
int rev[CN << 2];
void cg(COMP t[], int n){
for(int i = 0; i < n; i++) rev[i] = (rev[i >> 1] >> 1) + (i & 1) * (n >> 1);
for(int i = 0; i < n; i++) if(i < rev[i]) swap(t[i], t[rev[i]]);
}
void fft(COMP t[], int n, int tp){
cg(t, n);
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1; COMP gn = mk(cos(2 * PI / (LDB)w), sin(2 * tp * PI / (LDB)w));
for(int i = 0; i < n; i += w){
COMP g = mk(1, 0);
for(int j = i; j < i + l; j++){
COMP u = t[j], v = t[j + l] * g;
t[j] = u + v, t[j + l] = u - v, g = g * gn;
}
}
}
if(tp ^ 1){
for(int i = 0; i < n; i++) t[i].x /= (LDB)n, t[i].y /= (LDB)n;
}
}
COMP p[CN << 2], q[CN << 2], x[CN << 2], y[CN << 2], z[CN << 2], w[CN << 2];
void conv(int a[], int b[], int n){ // a = a * b
int B = (1 << 15) - 1, N = 1; while(N < (n << 1)) N <<= 1;
for(int i = 0; i < N; i++) p[i] = q[i] = mk(0, 0);
for(int i = 0; i < n; i++) p[i] = mk(a[i] >> 15, a[i] & B); // k1 r1
for(int i = 0; i < n; i++) q[i] = mk(b[i] >> 15, b[i] & B); // k2 r2
fft(p, N, 1), fft(q, N, 1);
for(int i = 0; i < N; i++){
int j = (N - 1) & (N - i);
COMP k1, r1, k2, r2;
k1 = (p[i] + conj(p[j])) * mk(0.5, 0);
r1 = (p[i] - conj(p[j])) * mk(0, -0.5);
k2 = (q[i] + conj(q[j])) * mk(0.5, 0);
r2 = (q[i] - conj(q[j])) * mk(0, -0.5);
x[i] = k1 * k2, y[i] = r1 * r2, z[i] = k1 * r2, w[i] = k2 * r1;
}
for(int i = 0; i < N; i++) p[i] = x[i] + y[i] * mk(0, 1);
for(int i = 0; i < N; i++) q[i] = z[i] + w[i] * mk(0, 1);
fft(p, N, -1), fft(q, N, -1);
for(int i = 0; i < n; i++){
LL X = (LL)(0.5 + p[i].x), Y = (LL)(0.5 + p[i].y), Z = (LL)(0.5 + q[i].x), W = (LL)(0.5 + q[i].y);
X = (X % P + P) % P, Y = (Y % P + P) % P, Z = (Z % P + P) % P, W = add((W % P + P) % P, Z);
a[i] = add((X << 30) % P, add((W << 15) % P, Y));
}
for(int i = n; i < N; i++) a[i] = 0;
}
int c[CN << 2];
void inv(int a[], int b[], int n){
for(int i = 0; i < (n << 1); i++) b[i] = c[i] = 0;
b[0] = invx(a[0]);
for(int w = 2; w < (n << 1); w <<= 1){
for(int i = 0; i < w; i++) c[i] = a[i];
conv(c, b, w);
for(int i = 0; i < w; i++) c[i] = add(0, P - c[i]); c[0] = add(2, c[0]);
conv(b, c, w);
}
}
int ia[CN << 2];
void ln(int a[], int b[], int n){
for(int i = 0; i < (n << 1); i++) ia[i] = b[i] = 0;
inv(a, ia, n);
for(int i = 0; i < (n << 1); i++) c[i] = 0;
for(int i = 0; i < n - 1; i++) c[i] = 1ll * (i + 1) * a[i + 1] % P;
conv(c, ia, n);
for(int i = 1; i < n; i++) b[i] = 1ll * c[i - 1] * invx(i) % P;
}
}
int p[CN], mu[CN]; bool np[CN];
void sieve(int n){
np[1] = 1, mu[1] = 1;
for(int i = 2; i <= n; i++){
if(!np[i]) p[++p[0]] = i, mu[i] = P - 1;
for(int j = 1; j <= p[0] && i * p[j] <= n; j++){
int x = i * p[j]; np[x] = 1;
if(i % p[j]) mu[x] = P - mu[i];
else break;
}
}
}
int n, a[CN << 2], b[CN << 2];
int main()
{
// freopen("_in.in", "r", stdin);
n = read() + 1, P = read();
a[0] = 1; for(int i = 1; i < n; i++) a[i] = read();
POLY :: ln(a, b, n);
sieve(n), memset(a, 0, sizeof(a));
for(int i = 1; i < n; i++) b[i] = 1ll * b[i] * i % P;
for(int i = 1; i < n; i++) for(int j = 1; i * j < n; j++)
a[i * j] = add(a[i * j], 1ll * mu[i] * b[j] % P);
int cnt = 0;
for(int i = 1; i < n; i++) cnt += (!!a[i]);
printf("%d\n", cnt);
for(int i = 1; i < n; i++) if(a[i]) printf("%d ", i);
return 0;
}又一道栗题
有一个整数 $x$，初值是 $[0,n]$ 内的一个随机值，其中为 $i$ 的概率为 $p_i$。
要进行 $m$ 轮操作，每次把 $x$ 等概率随机地替换成 $[0,x]$ 中的一个数字，求最后对 $\forall i\in[0, n]$，$x$ 的取值为 $i$ 的概率。对 $998244353$ 取模。
$n\le 10^5,m\le 10^{18}$
设 $f[k,i]$ 表示经过 $k$ 轮，得到 $i$ 的概率，有递推关系 $f[k,i]=\sum\limits_{j\ge i}^n f[k-1,j]/j+1$。
设 $F(k,z)=\sum\limits_{i=0}^n f[k,i]z^i$ 是 $f$ 的 OGF，根据上面的递推关系，有：
$$\begin{aligned} F(k,z)&=\sum\limits_{i=0}^n z^i\sum\limits_{j=i}^n \dfrac{f[k-1,j]}{j+1}\newline &=\sum\limits_{j=0}^n \dfrac{f[k-1,j]}{j+1}\sum\limits_{i=0}^j z^i \newline &=\dfrac{1}{z-1} \sum\limits_{i=0}^n f[k-1,i]\dfrac{z^{i+1}-1}{i+1} \end{aligned}$$
考虑 $f[k-1,i]$ 转移到了哪里，容易发现它在 $F(k,·)$ 中是 $i+1$ 项系数，这个变换不容易快速实现，那么继续推导。
根据牛顿-莱布尼茨公式，可以推得：
$\int_1^z t^i\text{ }\text dt=\dfrac{z^{i+1}}{z-1}-\dfrac{1}{z-1}$，那么原式化为：
$$\begin{aligned} F(k,z)&=\dfrac{1}{z-1}\sum\limits_{i=0}^n f[k-1,i]\int_1^z t^i\text{ }\text dt\newline &=\dfrac{1}{z-1}\int_1^z\left( \sum\limits_{i=0}^nf[k-1,i]t^i\right)\text{ }\text dt\newline &= \dfrac{1}{z-1}\int_1^z F(k-1,t)\text{ }\text dt \end{aligned}$$
考虑原式中 $z^{i+1}-1$ 很是鬼畜，考虑把定积分的下指标改为从 0 开始，那么设 $G(k,z)=F(k,z+1)$，有：
$$\begin{aligned} G(k,z)&=\dfrac{1}{z}\int_1^{z+1} F(k-1,t)\text{ }\text dt \newline &=\dfrac{1}{z}\int_0^z F(k-1,t+1)\text{ }\text d(t+1)\newline &=\dfrac{1}{z}\int_0^z G(k-1,t)\text{ }\text dt \newline &= \dfrac{1}{z}\sum\limits_{i=0}^n \dfrac{g[k-1,i]}{i+1}z^{i+1}\newline &=\sum\limits_{i=0}^n \dfrac{g[k-1,i]}{i+1}z^i \end{aligned}$$
其中 $g[k,i]$ 满足 $G(k,z)=\sum\limits_{i=0}^n g[k,i]z^i$，即是对应的序列。
考虑系数的变化：$k\to k+1:g[k,i]\to g[k+1,i]/(i+1)$，那么 $m$ 次变化之后就是 $g[m,i]=\dfrac{g[0,i]}{(i+1)^m}$ 这个就可以近似线性直接求出来。
那么只需要类似于 DFT 和 IDFT 那样子，把 $F$ 变成 $G$，系数修改完之后再变回来即可。
考虑怎么变，根据定义，得到：
$$\begin{aligned} G(k,z)&=F(k,z+1)\newline \Leftrightarrow \sum\limits_{i=0}^n g[k,i]z^i&=\sum\limits_{i=0}^n f[k,i](z+1)^i\newline &=\sum\limits_{i=0}^n\sum\limits_{j=0}^i\dbinom{i}{j}z^jf[k,i] \newline &= \sum\limits_{i=0}^n z^i \sum\limits_{j=i}^n \dbinom{j}{i}f[k,j]\end{aligned}$$
考虑 $k=0$，得：
$$\begin{aligned} g[0,i]&=\sum\limits_{j=i}^n \dbinom{j}{i}p_j\newline &= \dfrac{1}{i!}\sum\limits_{t=0}^{n-i} \dfrac{(i+t)!p_{i+t}}{t!}\end{aligned}$$
这里是差为定值的卷积，依然可以化成序列卷积，那么构造 $A(z)=\sum\limits_{i=0}^n \dfrac{1}{i!},B(z)=\sum\limits_{i=0}^n b_ix^i$，其中有 $p_ii!=b_{n-i}$，即人为令下标的和为定值，那么容易发现 $g[0,i]=\dfrac{1}{i!}[z^{n-i}]A(z)B(z)$。NTT 即可，复杂度 $O(n\log n)$。
对于 $g\to f$ 的转化，根据二项式反演，有：
$$\begin{aligned} g[k,i]&=\sum\limits_{j=i}^n\dbinom{j}{i}f[k,j] \newline \Leftrightarrow f[k,i]&=\sum\limits_{j=i}^n (-1)^{j-i}\dbinom{j}{i}g[k,j] \end{aligned}$$
依然是差为定值的卷积，同理可仿照上面求出，总复杂度 $O(n\log n)$。
代码：
#include<bits/stdc++.h>
using namespace std;
#define LL long long
const int P = 998244353;
const int CN = 1e5 + 10;
LL read(){
LL s = 0, ne = 1; char c = getchar();
for(;c < '0' || c > '9';c = getchar()) if(c == '-') ne = -1;
for(;c >= '0' && c <= '9';c = getchar()) s = (s << 1) + (s << 3) + c - '0';
return s * ne;
}
int add(int x, int y) {return x + y >= P ? x + y - P : x + y;}
int qp(int a, int b){
int r = 1;
while(b){
if(b & 1) r = 1ll * r * a % P;
a = 1ll * a * a % P, b >>= 1;
}
return r;
}
int inv(int x) {return qp(x, P - 2);}
int rev[CN << 2];
void cg(int t[], int n){
for(int i = 0; i < n; i++) rev[i] = (rev[i >> 1] >> 1) + (i & 1) * (n >> 1);
for(int i = 0; i < n; i++) if(i < rev[i]) swap(t[i], t[rev[i]]);
}
void ntt(int t[], int n, int tp){
cg(t, n);
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1, gn = qp(tp ? 3 : inv(3), (P - 1) / w);
for(int i = 0; i < n; i += w){
int g = 1;
for(int j = i; j < i + l; j++){
int u = t[j], v = 1ll * g * t[j + l] % P;
t[j] = add(u, v), t[j + l] = add(u, P - v), g = 1ll * g * gn % P;
}
}
}
if(!tp){
int in = inv(n);
for(int i = 0; i < n; i++) t[i] = 1ll * t[i] * in % P;
}
}
void conv(int A[], int B[], int n){ // A = A * B
int N = 1; while(N < (n << 1)) N <<= 1;
ntt(A, N, 1), ntt(B, N, 1);
for(int i = 0; i < N; i++) A[i] = 1ll * A[i] * B[i] % P;
ntt(A, N, 0);
for(int i = n; i < N; i++) A[i] = B[i] = 0;
}
int A[CN << 2], B[CN << 2];
int n, m, p[CN], fac[CN], ifac[CN], g[CN];
int main()
{
// freopen("_in.in", "r", stdin);
n = read() + 1, m = read() % (P - 1);
fac[0] = 1; for(int i = 1; i < n; i++) fac[i] = 1ll * fac[i - 1] * i % P;
ifac[n - 1] = inv(fac[n - 1]); for(int i = n - 2; i + 1; i--) ifac[i] = 1ll * ifac[i + 1] * (i + 1) % P;
for(int i = 0; i < n; i++) p[i] = read();
for(int i = 0; i < n; i++) A[n - i - 1] = 1ll * fac[i] * p[i] % P;
for(int i = 0; i < n; i++) B[i] = ifac[i];
conv(A, B, n);
for(int i = 0; i < n; i++) g[i] = 1ll * ifac[i] * A[n - i - 1] % P;
for(int i = 0; i < n; i++) g[i] = 1ll * g[i] * inv(qp(i + 1, m)) % P;
for(int i = 0; i < n; i++) A[n - i - 1] = 1ll * fac[i] * g[i] % P;
for(int i = 0; i < n; i++) B[i] = i & 1 ? P - ifac[i] : ifac[i];
conv(A, B, n);
for(int i = 0; i < n; i++) printf("%d ", int(1ll * ifac[i] * A[n - i - 1] % P)); puts("");
return 0;
}杂项
参考
- 《生成函数的运算与组合计数问题》金策，IOI中国国家候选队论文2015
- 《再探快速傅里叶变换》毛啸，IOI中国国家候选队论文2016
- 《浅谈 OI 中常用的一些生成函数运算的合法与正确性》，rpy’s Blog