【笔记】多项式全家桶
【笔记】多项式全家桶
Warning
- 空间记得开两倍,因为有卷积,最后结果是两多项式长度之和。
Template
p.s. 一般函数最开始是输出数组,后接输入数组,及其长度。
namespace NTT
{
const int gen = 3;
int r[N];
void ntt (int *f, int n, int sign) {
for (int i = 0; i < n; i++)
if (i < r[i]) swap(f[i], f[r[i]]);
for (int mid = 1; mid < n; mid<<=1) {
int o = qpow(gen, 1ll*sign*(mod-1)/(mid<<1)%mod+(mod-1));
for (int i = 0; i < n; i += mid<<1)
for (int j = 0, w = 1; j < mid; j++, w = 1ll*w*o%mod) {
int x = f[i+j], y = 1ll*w*f[i+mid+j]%mod;
f[i+j] = (1ll*x+y) % mod;
f[i+mid+j] = (1ll*x-y+mod) % mod;
}
}
if (sign == -1) {
int inv_n = ::inv(n);
for (int i = 0; i < n; i++) f[i] = 1ll*f[i]*inv_n%mod;
}
}
void mul (int *c, int *a, int *b, int n, int m, int (*op)(int, int)=[](int a, int b){ return int(1ll*a*b%mod); }) {
int lim = 1, len = 0;
while (lim <= n+m) lim <<= 1, len++;
for (int i = 0; i <= lim; i++)
r[i] = (r[i>>1]>>1)|((i&1)<<(len-1));
ntt(a, lim, 1), ntt(b, lim, 1);
for (int i = 0; i <= lim; i++) c[i] = op(a[i],b[i]);
ntt(c, lim, -1);
}
} // namespace NTT
namespace Poly
{
void derivative (int *g, int *f, int n) {
for (int i = 1; i <= n; i++) g[i-1] = 1ll*f[i] * i % mod;
g[n] = 0;
}
void integral (int *g, int *f, int n) {
for (int i = n; i >= 1; i--) g[i] = 1ll*f[i-1] * ::inv(i) % mod;
g[0] = 0;
}
void inv (int *g, int *f, int n) {
if (n == 0) return g[0] = ::inv(f[0]), void();
inv(g, f, n>>1);
int lim = 1; while (lim <= n<<1) lim <<= 1;
static int h[N];
for (int i = 0; i <= lim; i++) h[i] = i <= n ? f[i] : 0;
NTT::mul(g, h, g, n, n, [](int a, int b){
return int((2ll*b%mod-1ll*b*b%mod*a%mod+mod) % mod);
});
for (int i = n+1; i <= lim; i++) g[i] = 0;
}
void ln (int *g, int *f, int n) {
static int _f[N], inv_f[N];
int lim = 1; while (lim <= n<<1) lim <<= 1;
fill(_f, _f+lim+1, 0), fill(inv_f, inv_f+lim+1, 0);
derivative(_f, f, n), inv(inv_f, f, n);
NTT::mul(g, _f, inv_f, n, n);
integral(g, g, n);
}
void exp (int *g, int *f, int n) {
if (n == 0) return g[0] = 1, void(); // 题目保证 f 的常数项为 0
exp(g, f, n>>1);
int lim = 1; while (lim <= n<<1) lim <<= 1;
static int h[N]; ln(h, g, n);
for (int i = 0; i <= lim; i++) h[i] = i <= n ? ((1ll*(i==0)-h[i]+f[i]+mod)%mod) : 0;
NTT::mul(g, h, g, n, n);
for (int i = n+1; i <= lim; i++) g[i] = 0;
}
} // namespace Poly
p.s.
- 所有求
lim
的时候,理论上来说lim = n+m
的时候lim
是不需要再倍增的,但是有n+m = 1
的情况,这个时候len = 0
然后就 G 了,懒得特判所以就先打等于了。 NTT::mul
里面的*a,*b
不能相同。Poly::
里面除了derivative(),integral(),ln()
的*g,*f
可以相同,其余都不能相同。Poly::inv(),exp()
中的*g
必须是干净的。
0 前置知识
0.1 多项式求导与积分
0.1.1 求导
0.1.2 积分
单项式积分为:
\[\int x^a\mathrm dx=\dfrac 1{a+1}x^{a+1}
\]
多项式即为每个单项式积分再相加。对于每个单项式的系数,直接乘到最后的结果中。
1
1.1 FFT/NTT
1.2 多项式求逆
可以直接硬推递推式,也可以上牛顿迭代。
暂时没来得及写。
1.3 多项式对数函数 (ln)
给定 \(F(x)\),求 \(G(x)\) 满足 \(G(x)\equiv \ln F(x)\pmod{x^n}\)。
考虑到 \(\ln^\prime x=\dfrac 1x\),我们对等式两侧求导:(注意链式法则
\[G^\prime(x)=\ln^\prime F(x)=\dfrac {F^\prime(x)}{F(x)}
\]
然后再积分就得到了 \(G(x)\)。
1.4 多项式指数函数(exp)
给定 \(F(x)\),求 \(G(x)\) 满足 \(G(x)\equiv e^{F(x)}\pmod{x^n}\)。
\[G(x)\equiv e^{F(x)}\pmod{x^n}\\
\Rightarrow \ln G(x)\equiv F(x)\pmod{x^n}\\
\Rightarrow \ln (G(x))-F(x)\equiv 0\pmod{x^n}
\]
然后就可以上牛顿迭代了。
假设已经求出了 \(G_0(x)\equiv e^{F(x)}\pmod{x^{\left \lceil \frac n2 \right \rceil }}\)
由牛顿迭代:(注意现在的 \(G_0(x),F(x)\) 都是已知量,所以 \(F(x)\) 当成常数,\(G_0(x)\) 这个整体当作变量来处理
\[\begin{aligned}
G(x)&\equiv G_0(x)-\dfrac {\ln (G_0(x))-F(x)}{(\ln (G_0(x))-F(x))^\prime}\pmod {x^n}\\
&\equiv G_0(x)-\dfrac {\ln (G_0(x))-F(x)}{\frac1{G_0(x)}}\pmod {x^n}\\
&\equiv G_0(x)(1-\ln (G_0(x))+F(x))\pmod {x^n}\\
\end{aligned}
\]
1.5 多项式带余除法
给定 \(F(x),G(x)\)(项数分别为 \(n,m\)),求 \(Q(x),R(x)\)。满足:\(F(x)=G(x)Q(x)+R(x)\),且 \(Q(x)\) 的项数为 \(n-m\),\(R(x)\) 的项数小于 \(m\)。
现在问题主要在于,我们同时有两个待求的多项式。
一个很神奇的手法是,将这些多项式都翻转,严格的定义是:
\[\mathrm{rev}(F(x))=x^nF(x^{-1})
\]
于是可以推出:(过程可以参考开头的博客。