多项式求导
原理
如果有
那么由于求导的加法原则,每一项可以单独考虑,则有
实现
void Der(int *f, int *g, int lf) {
For (i, 1, lf) g[i - 1] = 1ll * i * f[i] % Mod; g[lf] = 0;
}
多项式积分
原理
如果有
同样由于积分的加法原则,每一项可以单独考虑,则有(常数项可以忽略)
实现
void Int(int *f, int *g, int lf) {
g[0] = 0;
For (i, 1, lf + 1)
g[i] = 1ll * f[i - 1] * inv[i] % Mod;
}
多项式乘法
原理
定义两个多项式 (A(x) = sum_{i = 0}^{n} a_i x^i, B(x) = sum_{i = 0}^m b_i x^i) 。那么乘积为 (C(x) = sum_{i = 0}^{n} sum_{j = 0}^m a_ib_jx^{i + j}) 。
我们考虑利用快速傅里叶变换解决,我们求出 (omega_n^k) 代入的所有的点值,然后用逆变换就可以求出所有系数了。
对于模意义下的我们利用 (g^{frac{p - 1}n}) 代替 (omega) (其中 (g) 为原根)就可以了。
实现
我们接下来都默认用费马质数来作为模数进行 ( ext{NTT}) 。
int rev[Maxn], W[Maxn], len;
void NTT(int *P, int opt) {
Rep (i, len) if (i < rev[i]) swap(P[i], P[rev[i]]);
for (int i = 2, p; p = i >> 1, i <= len; i <<= 1) {
W[0] = 1; W[1] = fpm(~opt ? g : invg, (Mod - 1) / i);
For (k, 2, p - 1) W[k] = 1ll * W[k - 1] * W[1] % Mod;
for (int j = 0; j < len; j += i) Rep (k, p) {
int u = P[j + k], v = 1ll * P[j + k + p] * W[k] % Mod;
P[j + k] = (u + v) % Mod; P[j + k + p] = (u - v + Mod) % Mod;
}
}
if (!~opt) {
int invn = fpm(len, Mod - 2);
Rep (i, len) P[i] = 1ll * P[i] * invn % Mod;
}
}
void Prepare(int lc) {
int cnt = -1; for (len = 1; len <= lc; len <<= 1) ++ cnt;
Rep (i, len) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << cnt);
}
int A[Maxn], B[Maxn], C[Maxn];
void Mult(int *a, int *b, int *c, int la, int lb) {
int lc = la + lb; Prepare(lc);
Rep (i, len) A[i] = i <= la ? a[i] : 0; NTT(A, 1);
Rep (i, len) B[i] = i <= lb ? b[i] : 0; NTT(B, 1);
Rep (i, len) C[i] = 1ll * A[i] * B[i] % Mod; NTT(C, -1);
For (i, 0, lc) c[i] = C[i];
}
多项式牛顿迭代
原理
给出函数 (f) ,找到一个多项式 (A),使得 (f(A(x)) equiv 0 pmod {x^n}) 。
在 (pmod x) 时,可以求出 (A(x)) 的常数项。
已经求出 (f(A(x)) equiv 0 mod {x^n}) ,现在求 (f(B(x)) equiv 0 mod {x^{2n}}) ,(B(x)) 与 (A(x)) 前 (n) 项系数相同。
在 (mod {x^{2n}}) 下,在 (A(x)) 处展开 (f(B(x))) 。
不难发现在 (i ge 2) 的时候,最低项已经超过了 (x^{2n}) ,所以只有 (i = 0, 1) 要考虑,也就是
大多数情况下,多项式复合是一个比较难求解的一个东西,但对于一些特殊的 (f) 可以方便求出。
多项式求逆
原理
好像这个不好套牛迭,我们现推算了。
设 (F(x)) 是 (P(x)) 在 (pmod {x^n}) 下的逆元,那么表示出来就是
当 (n = 1) 求逆元即可,假设我们求出在 (mod x^{frac n 2}) 意义下的逆元 (G(x)) ,那么有
我们有
那么递归求解即可,复杂度 (displaystyle mathcal T(n) = mathcal T(frac n2) + mathcal O(n log n) = mathcal O(n log n)) 。
实现
由于这个是大部分多项式操作的一个基础操作(加减乘除)之一,所以要尽量加快速度。
可以先 ( ext{NTT}) 然后做点值乘法,能少一半常数。(注意一开始需要是 (2^k ge len) 形式)
void Inv(int *f, int *g, int lf) {
if (lf == 1) return void(g[0] = fpm(f[0], Mod - 2));
Inv(f, g, lf >> 1); Prepare(lf << 1);
Rep (i, len) A[i] = i < lf ? f[i] : 0; NTT(A, 1);
Rep (i, len) B[i] = i < lf ? g[i] : 0; NTT(B, 1);
Rep (i, len) C[i] = 1ll * A[i] * B[i] % Mod * B[i] % Mod; NTT(C, -1);
Rep (i, lf) g[i] = (g[i] * 2ll + Mod - C[i]) % Mod;
}
多项式开根
原理
求 (F^2(x) equiv P(x) pmod {x^n}) ,此时的 (f(A(x)) = A^2 - P) 我们要求它的零点。
那么我们套用牛迭的式子就得到了
复杂度 (displaystyle mathcal T(n) = mathcal T(frac n2) + mathcal O(n log n) = mathcal O(n log n)) 。 。。主定理是真的强大
注意常数项非 (1) 的时候,通常都会保证存在模意义下的二次剩余,然后用 ( ext{BSGS}) 求出它是 (g^k) ,然后开根就是 (g^{k/2}) 。
实现
void Sqrt(int *f, int *g, int lf) {
if (lf == 1) return void(g[0] = getsqrt(f[0]));
Sqrt(f, g, lf >> 1); static int tmp[Maxn], tinv[Maxn];
Rep (i, lf) tmp[i] = 2 * g[i] % Mod; Inv(tmp, tinv, lf);
Mult(f, tinv, tmp, lf, lf);
const int inv2 = fpm(2, Mod - 2);
Rep (i, lf) g[i] = (1ll * g[i] * inv2 % Mod + tmp[i]) % Mod;
}
多项式ln
原理
求 (G(x) = ln F(x)) ,我们有
那么我们求导然后再求逆即可,复杂度 (mathcal O(n log n)) 。
常数项需要是 (1) ,其他形如 (ln x) 的无理数在模意义下无法表示。
实现
void Ln(int *f, int *g, int lf) {
static int der[Maxn], tmp[Maxn];
Der(f, der, lf); Inv(f, tmp, lf);
Mult(der, tmp, tmp, lf, lf); Int(tmp, g, lf);
}
多项式exp
原理
可以套用牛顿迭代和多项式 (ln) 解决。
分治的做法虽然多了个 (mathcal O(log)) 但好写许多(而且通常由于常数优势会比较快,除非你会论文哥那个牛迭优化),原理如下
我们考虑把 (A(x)) 先求导成 (A'(x)) ,令 (mid = lfloor frac{l + r}2 floor) 然后考虑分治求出 ([l, mid]) 的 (B(x)) ,然后考虑对于 ((mid, r]) 的贡献,也就是 (A_{0 sim r - l}) 和 (B_{l sim mid}) 卷然后贡献上去。
注意 (A(0) = 0) 才在模意义下有含义,此时 (B(0) = 1) ,以及积分需要平移一位并乘上 (frac{1}{n}) 。复杂度 (mathcal O(n log^2 n)) 。
实现
void Exp(int *a, int *b, int l, int r) {
if (l == r) return void(b[l] = (l ? 1ll * b[l] * inv[l] % Mod : 1));
int mid = (l + r) >> 1;
Exp(a, b, l, mid); static int tmp[Maxn];
Mult(a, b + l, tmp, r - l, mid - l);
For (i, mid + 1, r) b[i] = (b[i] + tmp[i - l - 1]) % Mod;
Exp(a, b, mid + 1, r);
}
void Exp(int *f, int *g, int lf) {
static int der[Maxn];
Der(f, der, lf); der[lf] = 0;
memset(g, 0, sizeof(int) * (lf + 1));
Exp(der, g, 0, lf);
}
多项式快速幂
原理
求 (G(x) = F^k(x)) ,我们取对数有 (ln G(x) = k ln F(x)) 也就是 (G(x) = exp(k ln F(x))) 。
但是 (F(0)) 可能不为 (1) ,那么我们把 (F(x)) 除掉 (F(0)) 即可。那么还可能为 (0) 我们平移即可,也就是除掉 (x^a) 。
然后 (exp) 结束记得乘回来它的 (k) 次幂即可。
代码
不想判 (0) 了。。
void Pow(int *f, int *g, int lf, int power) {
static int tf[Maxn], tmp[Maxn], invf; assert(f[0]);
invf = fpm(f[0], Mod - 2);
Rep (i, lf) tf[i] = 1ll * f[i] * invf % Mod;
Ln(tf, tmp, lf); Rep (i, lf) tmp[i] = 1ll * tmp[i] * power % Mod; Exp(tmp, g, lf);
invf = fpm(fpm(invf, Mod - 2), power);
Rep (i, lf) g[i] = 1ll * g[i] * invf % Mod;
}
多项式复合逆
原理
设 (G(F(x)) = x) 求 (F(x)) 的一项。此时 (G(x)) 可能是和 (F) 同阶的一个多项式,套用牛顿迭代项数会特别多,根本无法求解。
但要求单项即 ([x^n]F(x)) 有一个可以快速计算的式子
至于证明可以参看这篇博客 。
实现
int Lagrange(int *f, int lf, int x) {
static int tmp[Maxn], res[Maxn];
Rep (i, lf) res[i] = f[i + 1]; Inv(res, tmp, lf); Power(tmp, res, n);
return 1ll * res[x - 1] * fpm(x, Mod - 2) % Mod;
}