多项式全家桶
每次复习完下一次都会忘,这次下定决心一定要记下来!!!
FFT 和 NTT
做法
太简单了,看原题吧,不会退役吧。
代码
int n, m, rev[maxn];
int qpow(int x, int k, int p) {
int res = 1;
while(k) {
if(k & 1)
res = res * x % p;
k >>= 1, x = x * x % p;
}
return res;
}
void prepare(int len) {
for (int i = 0; i < len; i++) {
rev[i] = rev[i >> 1] >> 1;
if(i & 1)
rev[i] |= len >> 1;
}
}
struct Poly {
vector<int> a;
int size() {
return a.size();
}
int& operator[](int x) {
return a[x];
}
void resize(int N) {
a.resize(N);
}
void NTT(int f) {
for (int i = 1; i < a.size(); i++)
if(i < rev[i])
swap(a[rev[i]], a[i]);
for (int h = 2; h <= a.size(); h <<= 1) {
int d = qpow((f == 1 ? gb : gi), (mod - 1) / h, mod);
for (int i = 0; i < a.size(); i += h) {
int nw = 1;
for (int j = i; j < i + h / 2; j++) {
int a0 = a[j], a1 = a[j + h / 2] * nw % mod;
a[j] = (a0 + a1) % mod, a[j + h / 2] = (a0 - a1 + mod) % mod;
nw = nw * d % mod;
}
}
}
if(f == -1) {
int inv = qpow(a.size(), mod - 2, mod);
for (int i = 0; i < a.size(); i++)
a[i] = a[i] * inv % mod;
}
}
void read() {
for (int i = 0; i < a.size(); i++)
cin >> a[i];
}
void print() {
for (int i = 0; i < a.size(); i++)
cout << a[i] << " ";
cout << endl;
}
friend Poly operator*(Poly f, Poly g) {
int len = 1, t = f.size() + g.size() - 1;
while(len < t)
len <<= 1;
prepare(len);
f.resize(len), g.resize(len);
f.NTT(1), g.NTT(1);
for (int i = 0; i < len; i++)
f[i] = f[i] * g[i] % mod;
f.NTT(-1); f.resize(t);
return f;
}
}
多项式求逆
做法
假设我们已经求出来了:
那么我们要求:
我们先做差,得到:
然后平方,得到:
左右同时乘 \(A\),得到:
按这个柿子计算即可,注意没有必要把 NTT 出来的结果直接 NTT 回来,可以直接乘起来加起来。
代码
Poly get_inv(Poly f, int lim) {
if(lim == 1) {
Poly ans; ans.resize(1);
ans[0] = qpow(f[0], mod - 2, mod);
return ans;
}
Poly ans = get_inv(f, lim + 1 >> 1);
int len = 1;
while(len < lim * 2)
len <<= 1;
prepare(len);
f.resize(lim), f.resize(len), ans.resize(len);
f.NTT(1), ans.NTT(1);
for (int i = 0; i < len; i++)
ans.a[i] = ans.a[i] * (2 - f.a[i] * ans.a[i] % mod + mod) % mod;
ans.NTT(-1);
ans.resize(lim);
return ans;
}
注意,这里需要先 f.resize(lim),因为我们需要排除多余的位。
多项式开根
做法
我们假设 \(B'^2(x)=F(x)\pmod {x^{\lceil\frac{n}{2}\rceil}},B^2=F(x)\pmod{x^n}\)。
那么做差再平方得到:
直接做就可以。
代码
Poly get_sqrt(Poly f, int lim) {
if(lim == 1) {
Poly ans; ans.resize(1);
ans.a[0] = 1;
return ans;
}
Poly ans = get_sqrt(f, lim + 1 >> 1), inv = get_inv(ans, lim);
int len = 1;
while(len < lim * 2)
len <<= 1;
f.resize(lim), f.resize(len), ans.resize(len), inv.resize(len);
prepare(len);
// cout << lim << endl;
inv.NTT(1), f.NTT(1);
for (int i = 0; i < len; i++)
inv.a[i] = inv.a[i] * f.a[i] % mod;
inv.NTT(-1);
for (int i = 0; i < len; i++)
ans.a[i] = (mod + 1) / 2 * (inv.a[i] + ans.a[i]) % mod;
ans.resize(lim);
return ans;
}
一点记忆的小技巧:求逆和开根都是单式做的,所以要用倍增的想法去构造。
多项式除法
做法
我们记 \(F_R(x) = x^nF(\frac{1}{x})\),那么就有:
然后做就可以了。
代码
这里只给出核心代码,重载的运算符就没有了。
pair<Poly, Poly> get_div(Poly f, Poly g) {
int n = f.size(), m = g.size();
if(f.size() < g.size()) {
Poly q, r;
q.resize(n - m + 1); r = g;
return make_pair(q, r);
}
Poly q, r;
q = get_rev(f) * (get_inv(get_rev(g), n - m + 1));
q.resize(n - m + 1); q = get_rev(q);
r = f - g * q;
r.resize(m - 1);
return make_pair(q, r);
}
多项式 ln
做法
我们要求 \(G(x)=\ln F(x)\pmod {x^n}\)。
你可能会疑问为什么多项式会有 ln,我也不知道,但是就是对 EGF,OGF 的优化非常有用。
如何求呢?我们对两边同时求导,得到:
我们直接对 \(F\) 求逆然后乘上 \(F'\) 就可以了。
代码
int inv[maxn];
void prework_inv() {
inv[1] = 1;
for (int i = 2; i <= n; i++)
inv[i] = (mod - mod / i) * inv[mod % i] % mod;
}
Poly get_ji(Poly f) {
Poly g; g.resize(f.size() + 1);
g.a[0] = 0;
for (int i = 1; i < g.size(); i++)
g[i] = f[i - 1] * inv[i] % mod;
return g;
}
Poly get_dao(Poly f) {
Poly g; g.resize(f.size() - 1);
for (int i = 0; i < g.size(); i++)
g[i] = f[i + 1] * (i + 1) % mod;
return g;
}
Poly get_ln(Poly f) {
prework_inv();
Poly ans = get_ji(get_dao(f) * get_inv(f, f.size())); ans.resize(f.size());
return ans;
}
多项式 exp
牛顿迭代
我们先来看一个经典的问题:给定一个数 \(a\),求 \(b^2=a\),精确到小数点后若干位,并且手算。
那么你会说我会二分啊,但是比如 \(a=10^9\) 的时候二分其实挺难做的。
那么我们就可以请出牛顿迭代了。
我们可以把上面这个问题抽象为求 \(f(x) = x^2-a\) 的零点。
牛顿迭代的做法是我们先选取 \(t=a\) 作为我们的初始点,然后每次做一条 \(f(x)\) 在 \((t,f(t))\) 处的切线,找其与 \(x\) 轴的交点 \((t',0)\),然后令 \(t=t'\)。
根据很简单的知识我们知道这条直线为 \(y=f'(x)(x-t)+f(t)\)。那么零点为 \(x=t-\frac{f(t)}{f'(t)}\)。
我们发现,这样会一步步逼近我们的答案,并且次数非常小。
做法
我们对于多项式同样做牛顿迭代。
我们假设我们要求 \(ln(B(x))=A(x)\pmod {x^n}\),那么我们记 \(F(G(x)) = \ln(G(x))-A(x)\)。
那么我们如果求出了 \(ln(B_0(x)) = A(x)\pmod {x^{\lceil\frac{n}{2}\rceil}}\),那么我们就可以去做一次牛顿迭代就行。
那么此时有:
注意这里是求 \(F\) 的导,其实并不是一个复合函数,你可以把 \(G(x)\) 视为一个变量,所以不需要对里面求导。
按这个柿子做即可。
代码
Poly get_exp(Poly f, int lim) {
if(lim == 1) {
Poly ans; ans.resize(1); ans[0] = 1;
return ans;
}
Poly ans = get_exp(f, lim + 1 >> 1);
f.resize(lim), ans.resize(lim);
ans = ans * (1 - get_ln(ans) + f);
ans.resize(lim);
return ans;
}
多项式快速幂
做法
你也许会说直接类似于正常快速幂,但是复杂度 \(O(n\log^2n)\)。
其实我们可以变成 \(e^{k\ln F(x)}\) 这么求值就可以了,是单 \(\log\) 的。
代码
Poly qpow(Poly f, int k) {
return get_exp(get_ln(f) * k, f.size());
}
多项式多点求值
做法
我们考虑其实 \(f(a_i)=f(x)\pmod{x-a_i}\)。
那么我们就可以对这个东西分治,就是我们记 \(g(l,r)=\prod\limits_{i=l}^r(x-a_i)\),然后每次对左儿子和右儿子分别取模,一步步下去,复杂度 \(O(n\log^2n)\) 但是常数极大,需要一点逆天的优化(比如区间长度不大于一定值就暴力())。
据说有常数更小的厉害做法,但是我看不懂(),所以给出简单的这个做法,卡卡能过。
代码
int a[maxn];
Poly f[maxn], g;
void build(int l, int r, int t) {
if(l == r) {
f[t].resize(2), f[t][1] = 1, f[t][0] = mod - a[l];
return ;
}
int mid = l + r >> 1;
build(l, mid, t << 1), build(mid + 1, r, t << 1 | 1);
f[t] = f[t << 1] * f[t << 1 | 1];
}
void solve(int l, int r, Poly &g, int t) {
if(l + 2000 >= r) {
for (int i = l; i <= r; i++) {
int res = 0;
for (int j = g.size() - 1; j >= 0; j--) {
res = 1ll * res * a[i] % mod + g[j]; res %= mod;
}
cout << res << endl;
}
return ;
}
int mid = l + r >> 1;
Poly l1 = get_div(g, f[t << 1]).second;
solve(l, mid, l1, t << 1);
l1 = get_div(g, f[t << 1 | 1]).second;
solve(mid + 1, r, l1, t << 1 | 1);
}

浙公网安备 33010602011771号