多项式合集
前言
本来不想多项式的，但奈何做 P6667 的时候发现很多前置知识一直不会，遂开文，正好补了往年的遗憾，也把代码翻新一下。
这个标题是每到冬天我才想做多项式题。
多项式乘法逆
P4238：给定一个多项式 $F(x)$ ，请求出一个多项式 $G(x)$ ， 满足 $F(x) * G(x) \equiv 1\pmod{x^n}$ 。系数对 $998244353$ 取模。
推导
以下除以二是向上取整。
设已经求出 $G'(x)\equiv F(x)^{-1}\pmod{x^{\frac{n}{2}}}$ ，将 $G(x)\equiv F(x)^{-1}\pmod{x^n}$ 与之相减得
$$G(x)-G'(x)\equiv 0\pmod{x^{\frac{n}{2}}} $$
式子两边同时取平方
$$G(x)^2-2\cdot G(x)\cdot G'(x)+G'(x)^2\equiv 0\pmod{x^n} $$
式子两边同时乘 $F(x)$
$$\begin{aligned}G(x)-2\cdot G'(x)+G'(x)^2\cdot F(x)&\equiv 0&\pmod{x^n}\\ G(x)&\equiv 2\cdot G'(x)-G'(x)^2\cdot F(x)&\pmod{x^n} \end{aligned}$$
代码
可以递归也可倍增，递归常数比倍增大。
递归
void inv (int *f, int *g, int m) {
if (m == 1) { g[0] = qpow(f[0], mod - 2); return; }
inv (f, g, (m + 1) >> 1);
int lim, lsiz; for (lim = 1; lim < (m << 1); lim <<= 1, lsiz++);
init (lim, lsiz);
static int tmp[pN]; cpy(tmp, f, m); clr(tmp + m, lim - m);
NTT(tmp, lim), NTT(g, lim);
for (int i = 0; i < lim; ++i) g[i] = mul(dec(2, mul(g[i], tmp[i])), g[i]);
NTT(g, lim, 0);
clr(g + m, lim);
}
倍增
void px(int *f, int *g, int lim) { for (int i = 0; i < lim; ++i) f[i] = mul(f[i], g[i]); }
void inv (int *f, int n) {
static int g[pN], r[pN], tmp[pN];
g[0] = qpow(f[0], mod - 2);
int lim = 1, lsiz = 0;
for (int len = 2; (len >> 1) <= n; len <<= 1) {
lim = len, lsiz++; init(lim, lsiz);
cpy(r, g, len >> 1); cpy(tmp, f, lim);
NTT(tmp, lim); NTT(r, lim); px(r, tmp, lim); NTT(r, lim, 0); clr(r, len >> 1);
cpy(tmp, g, len);
NTT(tmp, lim); NTT(r, lim); px(r, tmp, lim); NTT(r, lim, 0);
for (int i = (len >> 1); i < len; ++i) g[i] = dec(mul(g[i], 2), r[i]);
}
cpy(f, g, n); clr(g, lim); clr(r, lim); clr(tmp, lim);
}
多项式除法
P4512：给定一个 $n$ 次多项式 $F(x)$ 和一个 $m$ 次多项式 $G(x)$ ，请求出多项式 $Q(x)$ , $R(x)$ ，满足以下条件：
- $Q(x)$ 次数为 $n-m$ ， $R(x)$ 次数小于 $m$
- $F(x) = Q(x) \times G(x) + R(x)$
所有的运算在模 $998244353$ 意义下进行。
推导
令 $R(x)$ 次数为 $m-1$ ， $F_R(x)$ 表示 $F(x)$ 翻转系数的多项式，便有性质 $F_R(x)=x^nF(\frac{1}{x})$ ，其中 $n$ 是次数。
把 $\frac{1}{x}$ 代入 $F(x) = Q(x) \times G(x) + R(x)$
$$\begin{aligned} F\left(\frac{1}{x}\right) &= Q\left(\frac{1}{x}\right) \times G\left(\frac{1}{x}\right) + R\left(\frac{1}{x}\right)\\ x^nF\left(\frac{1}{x}\right) &= x^{n-m}Q\left(\frac{1}{x}\right) \times x^mG\left(\frac{1}{x}\right) + x^{n-m+1}\cdot x^{m-1}R\left(\frac{1}{x}\right)\\ F_R(x) &= Q_R(x)\times G_R(x)+x^{n-m+1}R_R(x)\\ F_R(x) &\equiv Q_R(x)\times G_R(x)\pmod{x^{n-m+1}}\\ Q_R(x) &\equiv \frac{F_R(x)}{G_R(x)}\pmod{x^{n-m+1}} \end{aligned}$$
知一可得二矣。
代码
void div (int *f, int *g, int n, int m) {
static int fR[pN], gR[pN];
int L = n - m + 1;
reverse(f, f + n); cpy(fR, f, L); reverse(f, f + n);
reverse(g, g + m); cpy(gR, g, L); reverse(g, g + m);
inv(gR, L); times(gR, fR, L, L, L); reverse(gR, gR + L);
times(g, gR, n, n, n);
for (int i = 0; i < m - 1; i++) g[i] = dec(f[i], g[i]);
clr (g + m - 1, L);
cpy(f, gR, L); clr(f + L, n - L);
}
模操作类似
void mof (int *f, int *g, int n, int m) {
static int fR[pN], gR[pN];
int L = n - m + 1;
reverse(f, f + n); cpy(fR, f, L); reverse(f, f + n);
reverse(g, g + m); cpy(gR, g, L); reverse(g, g + m);
inv(gR, L); times(gR, fR, L, L, L); reverse(gR, gR + L);
times(gR, g, L, m, m - 1);
for (int i = 0; i < m - 1; i++) f[i] = dec(f[i], gR[i]);
clr (f + m - 1, L); clr(gR, n); clr(fR, n);
}
多项式对数函数
P4725：给定一个多项式 $F(x)$ ，请求出一个多项式 $G(x)$ ， 满足 $G(x) \equiv \ln F(x)\pmod{x^n}$ 。系数对 $998244353$ 取模。
推导
导完积回去🥵
$$\begin{aligned} G(x) &\equiv \ln F(x)&\pmod{x^n}\\ G'(x) &\equiv \frac{F'(x)}{F(x)}&\pmod{x^n}\\ G(x) &\equiv \int\frac{F'(x)}{F(x)}&\pmod{x^n}\\ \end{aligned}$$
没什么技术含量。
代码
void der (int *f, int lim) {
for (int i = 0; i < lim; ++i) f[i - 1] = mul(f[i], i);
f[lim - 1] = 0;
}
void summa (int *f, int lim) {
for (int i = lim; i >= 1; --i) f[i] = mul(f[i - 1], inv[i]);
f[0] = 0;
}
void ln (int *f, int n) {
static int g[pN];
cpy(g, f, n); der(g, n);
inv(f, n); times(f, g, n, n, n);
summa(f, n - 1); clr(g, n);
}
多项式指数函数
P4726：给定一个多项式 $F(x)$ ，请求出一个多项式 $G(x)$ ， 满足 $G(x) \equiv \exp( F(x))\pmod{x^n}$ 。系数对 $998244353$ 取模。
泰勒展开
对于函数，泰勒展开式
$$f(x)=\sum_{i=0}\frac{f^{(i)}(x_0)}{i!}(x-x_0)^i+R_n(x) $$
其中 $R_n(x)$ 是佩亚诺余项。等式右边构造了一个多项式一样的函数与等式左边的 $f(x)$ 无限接近。详细地，以 $x=x_0$ 处为起点构造函数 $g(x)=\sum_ia_i(x-x_0)^i=f(x)$ ，等式在 $x=x_0$ 处同时求导多次，对于第 $n$ 次求导 $g^{(n)}(x_0)=n!a_n=f^{(n)}(x_0)\Rightarrow a_n=\frac{f^{(n)}(x_0)}{n!}$ 。
而对于多项式泰勒展开式相对更为简单
$$F(x)=\sum_{i=0}\frac{F^{(i)}(x_0)}{i!}(x-x_0)^i $$
牛顿迭代
对于一个函数 $f(x)$ ，想要找到它的某个零点 $x_0$ ，即满足 $f(x_0)=0$ ，可以使用牛顿迭代。
方法是，对于一个点 $(x_1,0)$ ，可以取得比它更靠近 $(x_0,0)$ 的点，该点为，经过它的 $f(x)$ 的切线 $y-f(x_1)=f'(x_1)(x-x_1)$ 与 $y=0$ 的交点，记为 $(x_2,0)$ ，有 $x_2=x_1-\frac{f(x_1)}{f'(x_1)}$ ……如此 $(x_\infty,0)$ 就是 $(x_0,0)$ 。
类似地，对于多项式 $F(x),G(x)$ 满足 $F(G(x))=0$ ，若已知 $F(H(x))\equiv 0\pmod{x^\frac{n}{2}}$ 则有 $G(x)\equiv H(x)-\frac{F(H(x))}{F'(H(x))}\pmod{x^n}$ 。
这有更为严谨的证明，把 $F(G(x))$ 在 $H(x)$ 处泰勒展开
$$F(G(x))=\sum_{i=0}\frac{F^{(i)}(H(x))}{i!}(G(x)-H(x))^i $$
$G(x)\bmod x^n$ 可作 $x^{\frac n2}$ 的解，则有
$$\begin{aligned} G(x)&\equiv H(x)&\pmod{x^{\frac n2}}\\ G(x)-H(x)&\equiv 0&\pmod{x^{\frac n2}}\\ (G(x)-H(x))^2&\equiv 0&\pmod{x^n}\\ \end{aligned}$$
所以
$$\begin{aligned} F(G(x))&\equiv F(H(x))+F'(H(x))(G(x)-H(x))&\pmod{x^n}\\ G(x)&\equiv H(x)-\frac{F(H(x))}{F'(H(x))}&\pmod{x^n} \end{aligned}$$
推导
$$\begin{aligned} G(x)&\equiv\exp(F(x))&\pmod {x^n}\\ \ln G(x)&\equiv F(x)&\pmod{x^n}\\ \end{aligned}$$
令 $P(G(x))=\ln G(x)-F(x)$ ，则
$$\begin{aligned} P(G(x))&\equiv0&\pmod {x^n}\\ G(x)&\equiv H(x)-\frac{\ln H(x)-F(x)}{\frac{1}{H(x)}}&\pmod{x^n}\\ G(x)&\equiv H(x)(1-\ln H(x)+F(x))&\pmod{x^n}\\ \end{aligned}$$
代码
void exp (int *f, int n) {
static int tmp[pN], g[pN];
g[0] = 1;
int lim = 1;
for (int len = 2; (len >> 1) <= n; len <<= 1) {
lim = len;
cpy(tmp, g, len >> 1); ln(tmp, len);
for (int i = 0; i < len; ++i) tmp[i] = dec(f[i], tmp[i]);
tmp[0] = add(tmp[0], 1);
times(g, tmp, len, len, len);
}
lim >>= 1;
cpy(f, g, n); clr(g, lim); clr(tmp, lim);
}
多项式开根
P5205&P5277：给定多项式 $F(x)$ ，求多项式 $G(x)$ 满足 $G^2(x)\equiv F(x)\pmod{x^n}$ 。系数对 $998244353$ 取模。
推导
仍是牛顿迭代，令 $P(G(x))=G^2(x)-F(x)$ ，套公式化简后得
$$G(x)\equiv\frac{H^2(x)+F(x)}{2H(x)}\pmod{x^n} $$
边界条件二次剩余即可。
代码
int N1;
struct Complex {
int a, b;
Complex (int a = 0, int b = 0): a(a), b(b) {}
};
bool operator == (Complex a, Complex b) { return a.a == b.a && a.b == b.b; }
Complex operator * (Complex a, Complex b) { return Complex(add(mul(a.a, b.a), mul(N1, mul(a.b, b.b))),add(mul(a.b, b.a), mul(a.a, b.b))); }
Complex qpowC (Complex a, ll b) {
Complex ans = Complex(1, 0);
for (; b; b >>= 1, a = a * a)
if (b & 1) ans = ans * a;
return ans;
}
mt19937 rnd(114514);
int sqrtI (int n) {
int a = rnd() % mod;
N1 = dec(mul(a, a), n);
for (; qpow(N1, (mod - 1) / 2) != mod - 1; a = rnd() % mod, N1 = dec(mul(a, a), n));
int x1 = qpowC(Complex(a, 1), (mod + 1) / 2).a, x2 = dec(mod, x1);
if (x1 > x2) swap(x1, x2);
return x1;
}
void pow2 (int *f, int lim, int lsiz) {
init(lim, lsiz);
NTT(f, lim); px(f, f, lim); NTT(f, lim, 0);
}
void sqrtp (int *f, int n) {
static int tmp[pN], g[pN];
g[0] = sqrtI(f[0]);
int lim = 1, lsiz = 0;
for (int len = 2; (len >> 1) <= n; len <<= 1) {
lim = len, lsiz++;
for (int i = 0; i < (len >> 1); ++i) tmp[i] = add(g[i], g[i]);
inv(tmp, lim); pow2(g, lim, lsiz);
for (int i = 0; i < len; ++i) g[i] = add(g[i], f[i]);
times(g, tmp, lim, lim, lim);
}
cpy(f, g, n); clr(g, lim); clr(tmp, lim);
}
多项式幂函数
P5245&P5273：给定 $F(x),k\ (k\leq 10^{10^5})$ ，求 $G(x)$ 满足： $G(x)\equiv F^k(x)\pmod{x^n}$ 。
推导
$$F^k(x)=\exp(k\ln F(x)) $$
对于常数项，若非零，整个多项式先除以常数项，做完乘回去；若为零则先左移后面移回去。
代码
void powk (int *f, int n, int k, int kphi) {
int z = 0, c;
for (; !f[z]; z++); c = qpow (f[z], mod - 2);
if (z && (1ll * z * k >= n || strlen(s + 1) > 5)) { clr(f, n); return; }
n -= z * k;
for (int i = 0; i < n; ++i) f[i] = mul(f[i + z], c);
clr(f + n, z * k); ln(f, n); for (int i = 0; i < n; ++i) f[i] = mul(f[i], k); exp(f, n);
c = qpow (c, mod - 2); c = qpow (c, kphi);
static int g[pN];
for (int i = 0; i < n; i++) g[z * k + i] = mul(f[i], c);
cpy(f, g, z * k + n); clr(g, z * k + n);
}
多项式多点求值
我就是为了这碟醋才包了这盘饺子。
P5050：给定 $F(x)$ 和 $m$ 个整数 $a_i$ ，分别求出 $F(a_i)$ 。
多项式除法取模做法
常数较大的方法。
处理分治情况，设
$$G_0(x)=\prod_{i=1}^{\lfloor\frac m2\rfloor}(x-a_i),G_1(x)=\prod_{i=\lfloor\frac m2\rfloor}^{n}(x-a_i), $$
如此代入 $i\in[1,\lfloor\frac m2\rfloor]$ 进 $G_0(x)$ 时， $G_0(i)=0$ ，右边同理。
再令 $R_0(x),R_1(x)$ 满足 $F(x)=Q_0(x)\times G_0(x)+R_0(x)=Q_1(x)\times G_1(x)+R_1(x)$ 。
分治过去只剩下 $F(x)=R_0(x),F(x)=R_1(x)$ 了。
那么先分治 NTT 求 $G(x)$ ，再多项式除法得 $R(x)$ 。
时间复杂度 $\mathcal{O}(n\log^2n)$ 。
namespace Evaluation {
using namespace poly;
int gl[17][pN];
void prework (int dep, int l, int r, int *a) {
if (l == r) {
gl[dep][l << 1] = mod - a[l];
gl[dep][l << 1 | 1] = 1;
return;
}
int mid = (l + r) >> 1;
prework(dep + 1, l, mid, a);
prework(dep + 1, mid + 1, r, a);
cpy(&gl[dep][l << 1], &gl[dep + 1][l << 1], mid - l + 2);
times(&gl[dep][l << 1], &gl[dep + 1][(mid + 1) << 1], mid - l + 2, r - mid + 1, r - l + 2);
}
int s1[pN], s2[pN];
void solve (int dep, int l, int r, int *f, int *x, int *y) {
if (r - l <= 350) {
for (int i = l; i <= r; ++i) {
int buf = 1;
for (int j = l + l; j < l + r + 1; ++j) {
y[i] = add(y[i], mul(buf, f[j]));
buf = mul(buf, x[i]);
}
}
return;
}
int mid = (l + r) >> 1;
cpy(s1, &f[l << 1], r - l + 1); mof(s1, &gl[dep + 1][l << 1], r - l + 1, mid - l + 2);
cpy(s2, &f[l << 1], r - l + 1); mof(s2, &gl[dep + 1][(mid + 1) << 1], r - l + 1, r - mid + 1);
for (int i = (l << 1); i < (r << 1 | 1); i++) f[i] = 0;
cpy(&f[l << 1], s1, mid - l + 1);
cpy(&f[(mid + 1) << 1], s2, r - mid);
clr(s1, r - l + 1); clr(s2, r - l + 1);
solve(dep + 1, l, mid, f, x, y);
solve(dep + 1, mid + 1, r, f, x, y);
}
void eval (int *f, int *x, int *y, int n, int m) {
prework(0, 0, m - 1, x);
if (n > m) mof(f, gl[0], n, m + 1);
solve(0, 0, m - 1, f, x, y);
}
}
转置算法
不会喵。