1 多项式求逆
求 对于 (A(x)) 求 (B(x)) 使得 (A(x) * B(x) equiv 1 (mod n))
[AB equiv 1 (mod n)
]
[AB - 1equiv 0(mod n)
]
[(AB - 1) ^ 2 equiv 0 (mod 2n)
]
因为
[F(x) = sumlimits_{i = 0}^ {n - 1} sumlimits_{j = 0} ^ i f(j) * g(i - j)
]
(f(i) * g(i - j)) 其中必有 一项为零， 之后整个多项式次数翻倍
[A^2B^2 - 2AB + 1 equiv 0 (mod 2n)
]
[2AB - A^2B^2 equiv 1 (mod 2n)
]
[A (2B - AB^2) equiv 1(mod 2n)
]
[2B -AB^2 equiv A^{-1}(mod 2n)
]
一个多项式有逆当且仅当 (A_0) 有逆元
我们可以求出 (A^{-1} equiv A (mod x^1))
这样就可以倍增的求逆元了
(code) 代码主要慢在取模，优化还是很快的 船新版本，封装还快了3倍
#include <bits/stdc++.h>
using namespace std;
#define rg register
#define gc getchar
#define rep(i, a, b) for(int i = a; i <= b; ++i)
inline int read(){
rg char ch = gc();
rg int x = 0, f = 0;
while(!isdigit(ch)) f |= (ch == '-'), ch = gc();
while(isdigit(ch)) x = (x << 1) + (x << 3) + (ch ^ 48), ch =gc();
return f ? -x : x;
}
const int N = 4e5 + 5, mod = 998244353;
inline int ksm(int a, int b){
int ans = 1;
while(b){
if(b & 1) ans = 1ll * a * ans % mod;
b >>= 1;
a = 1ll * a * a % mod;
}
return ans;
}
const int G = 3, Gn = ksm(G, mod - 2);
int n;
struct FFT{
int to[N], a[N], b[N], A[N], B[N], n;
inline void redef(int m, int *c){
memcpy(a, c, ((n = m) + 1) * sizeof(int));
}
inline void NTT(int *a, int lim, int len, int flag){
for(int i = 0; i < lim; ++i) to[i] = (to[i >> 1] >> 1) | ((i & 1) << len - 1);
for(int i = 0; i < lim; ++i) if(to[i] > i) swap(a[i], a[to[i]]);
for(int l = 2; l <= lim; l <<= 1){
const int m = l >> 1, Gi = ksm(flag == 1 ? G : Gn, (mod - 1) / l);
for(int j = 0; j < lim; j += l){
int g = 1;
for(int i = 0; i < m; ++i, g = 1ll * g * Gi % mod){
int x = a[i + j], y = 1ll * g * a[i + j + m] % mod;
a[i + j] = x + y;
if(a[i + j] >= mod) a[i + j] -= mod;
a[i + j + m] = x + mod - y;
if(a[i + j + m] >= mod) a[i + j + m] -= mod;
}
}
}
}
inline void getni_a_to_b(){
b[0] = ksm(a[0], mod - 2);
int lim = 1, len = 0;
while(lim <= n) lim <<= 1, ++len;
int nlen = 1;
for(int xmod = 1; xmod < lim; xmod <<= 1){
++nlen;
int nlim = xmod << 2;
fill(A + (xmod << 1), A + nlim, 0); fill(B + xmod, B + nlim, 0);
memcpy(A, a, (xmod << 1) * sizeof(int)); memcpy(B, b, xmod * sizeof(int));
NTT(A, nlim, nlen, 1); NTT(B, nlim, nlen, 1);
for(int i = 0; i < nlim; ++i) A[i] = (2ll * B[i] + mod - 1ll * A[i] * B[i] % mod * B[i] % mod) % mod;
NTT(A, nlim, nlen, -1);
const int inv = ksm(nlim, mod - 2);
for(int i = 0; i < (xmod << 1); ++i) b[i] = 1ll * A[i] * inv % mod;
}
}
}T;
int g[N];
signed main(){
n = read() - 1;
rep(i, 0, n) g[i] = read();
T.redef(n, g);
T.getni_a_to_b();
rep(i, 0, n) printf("%d ", T.b[i]);
//注意求出来的逆元是在 mod n 意义下的，如果 n 变大，就不对了
gc(), gc();
return 0;
}
2 多项式开方
对于 (A(x)) 求 (B(x)) 使得 (B(x) * B(x) equiv A(x) (mod x^m))
由于
[A(x) equiv 0 (mod x^m) -> A(x)^2 equiv0 (mod x^{2m})
]
考虑倍增构造
已知
[B^2 equiv A (mod x^m)
]
[(B^2 - A)^2 equiv 0 (mod x ^ {2m})
]
[(B^2 + A)^2 equiv 4AB^2(mod x^{2m})
]
[((B^2 + A) / 2B)^2 equiv A(mod x^{2m})
]
[A^{frac 1 2} equiv (B^2 + A) / 2B (mod x^{2m})
]
对于上面的式子直接点值相除即可
复杂度 (O(n log^2n))
3 多项式牛顿迭代
(F(X)) 是一个对于多项式的函数，我们要求得多项式 (B) 使得 (F(B) equiv 0 (mod x^{2^n}))
考虑倍增解法
我们已知 (F(A) equiv 0 (mod x^{2^{n - 1}})) 怎么来解出 (B,F(B) equiv 0 (mod x^{2^n}))
对 (A) 泰勒展开
[F(B) = F(A) + frac{F'(A)} {1!} (B - A) + frac{F''(A)} {2!} (B - A)^2 +...
]
因为每个模数的答案是唯一的，且缩短最终答案一定是当前答案， (B) 的前 (2^{n - 1} - 1) 项 (= A)，这样 (B - A equiv 0 (mod x^{2^{n - 1}}), (B - A)^2 equiv 0 (mod x^{2^n}))
如果是在 (mod x^{2^n}) 意义下 后面包含 ((B - A)^2) 就都是 (0) 了
[F(B) = F(A) + F'(A)(B - A) (mod x^{2^n})
]
我们要求 (F(B) equiv 0 (mod x^{2^{n - 1}})) 解得
[B = A - frac{F(A)} {F'(A)} (mod x^{2^n})
]
4 多项式ln
求
[B = ln(A)
]
[B = int frac {dln(A)} {dA} frac {dA} {dx} dx
]
[B = int frac {A'} A dx
]
多项式求逆就行了，顺带一提 (B = int A dx) 表示 (B) 是 (A) 的不定积分，$int kx^a dx = frac k {a + 1} x^{a + 1} $
5 多项式exp
给定多项式 (A) 求
[B = e^A (mod x^n)
]
[ln(B) = A
]
设 (F(H) = ln(H) - A) 答案就是它的根
还是牛顿迭代一下
[egin{aligned}
B' &= B - frac{F(B)} {F'(B)} \
&= B - frac{ln(B) - A} {frac 1 B} \
&= B(1 - ln(B) + A)
end{aligned}]
时间复杂度 (O(n logn + T(frac n 2)) = O(n log n))
6 多项式除法/取余
给出次数分别为 (n, m) 的多项式 (A(x),B(x)) 求 (C(x))
[A(x) = B(x) * C(x) + R(x)
]
[A(frac 1 x) = B(frac 1 x) * C(frac 1 x) + R(frac 1 x)
]
同乘 (x^n)
[x^nA(frac 1 x) = x^mB(frac 1 x) * x^{n - m}C(frac 1 x) + x^nR(frac 1 x)
]
因为在多项式除法定义下，(R) 的次幂一定 (< B) 的次幂，不妨设其为 (m - 1)
设 (A^R(x) = x) 的 (A) 的次幂次方乘上 (A(frac 1 x))，可以发现这样就是把系数翻转了一下
[A^{R}(x) = B^R(x) * C^R(x) + R^R(x)n^{n - m + 1}
]
把它膜拜一下
[A^{R}(x) = B^R(x) * C^R(x) ~~ mod(x ^ {n - m + 1})
]
我们发现 (C(x)) 的最高次幂正好是 (n - m) 也就是说模了之后不影响 (C)
那多项式求逆就好了
#include <bits/stdc++.h>
using namespace std;
#define rg register
#define gc getchar
#define rep(i, a, b) for(int i = a; i <= b; ++i)
inline int read(){
rg char ch = gc();
rg int x = 0, f = 0;
while(!isdigit(ch)) f |= (ch == '-'), ch = gc();
while(isdigit(ch)) x = (x << 1) + (x << 3) + (ch ^ 48), ch = gc();
return f ? -x : x;
}
const int N = 4e5 + 5;
int A[N], B[N], a[N], aR[N], b[N], bR[N], bR_inv[N], c[N], cR[N], r[N], rR[N], to[N], n, m;
const int mod = 998244353;
inline int ksm(int a, int b){
int ans = 1;
while(b){ if(b & 1) ans = 1ll * a * ans % mod; b >>= 1; a = 1ll * a * a % mod; }
return ans;
}
const int G = 3, Gn = ksm(G, mod - 2);
inline void NTT(int *a, int lim, int flag){
rep(i, 0, lim - 1) if(to[i] > i) swap(a[i], a[to[i]]);
for(int l = 2; l <= lim; l <<= 1){
const int m = l >> 1, Gi = ksm(flag == 1 ? G : Gn, (mod - 1) / l);
for(int j = 0; j < lim; j += l){
int g = 1;
for(int i = 0; i < m; ++i, g = 1ll * g * Gi % mod){
int x = a[i + j], y = 1ll * g * a[i + j + m] % mod;
a[i + j] = x + y;
if(a[i + j] >= mod) a[i + j] -= mod;
a[i + j + m] = x + mod - y;
if(a[i + j + m] >= mod) a[i + j + m] -= mod;
}
}
}
}
inline void get_ni(int *a, int *b, int n){
int lim = 1;
while(lim <= n) lim <<= 1;
int len = 1;
b[0] = ksm(a[0], mod - 2);
for(int xmod = 1; xmod < lim; xmod <<= 1){
++len;
int nlim = xmod << 2;
rep(i, 0, nlim - 1) to[i] = (to[i >> 1] >> 1) | ((i & 1) << len - 1);
memcpy(A, a, (xmod << 1) * sizeof(int)); memcpy(B, b, xmod * sizeof(int));
fill(B + xmod, B + (xmod << 1), 0); //fill(A + (xmod << 1), A + nlim, 0);
NTT(A, nlim, 1); NTT(B, nlim, 1);
rep(i, 0, nlim - 1) A[i] = (B[i] * 2 - 1ll * A[i] * B[i] % mod * B[i] % mod + mod) % mod;
NTT(A, nlim, -1);
const int inv = ksm(nlim, mod - 2);
rep(i, 0, (xmod << 1) - 1) b[i] = 1ll * inv * A[i] % mod;
}
}
inline void mul(int *d, int *b, int n, int m, int *c){
int lim = 1, len = 0;
while(lim <= n + m) lim <<= 1, ++len;
fill(d + n + 1, d + lim, 0); fill(b + m + 1, b + lim, 0);
rep(i, 0, lim - 1) to[i] = (to[i >> 1] >> 1) | ((i & 1) << len - 1);
NTT(d, lim, 1); NTT(b, lim, 1);
rep(i, 0, lim - 1) d[i] = 1ll * d[i] * b[i] % mod;
NTT(d, lim, -1);
int inv = ksm(lim, mod - 2);
rep(i, 0, lim - 1) c[i] = 1ll * d[i] * inv % mod;
}
int main(){
n = read(), m = read();
rep(i, 0, n) a[i] = aR[n - i] = read();
rep(j, 0, m) b[j] = bR[m - j] = read();
get_ni(bR, bR_inv, n - m);
mul(aR, bR_inv, n, n - m, cR);
for(int i = n - m; ~i; --i) printf("%d ", cR[i]); puts("");
rep(i, 0, n - m) c[i] = cR[n - m - i];
mul(b, c, m, n - m, b);
for(int i = 0; i < m; ++i) printf("%d ", (a[i] - b[i] + mod) % mod);
gc(), gc();
return 0;
}
7 多项式快速幂
给定多项式 (A(X)) 求 (B(X) equiv A^k(X) mod(x ^ n))
[(A(X))^k = e^{ln(A^k(x))}
]
[=e^{kln(A(x))}
]
多项式 $ln + $ 多项式 (exp) 就行了
轮子
多项式乘法 + 多项式求逆 + 多项式开根 + 多项式(ln) + 多项式(exp) + 多项式快速幂 + FWT + 快速找原根
至于为什么没除法，因为我觉着生成函数不太可能整，而且加上不美观。。。
#include <bits/stdc++.h>
using namespace std;
#define rg register
#define gc getchar
#define rep(i, a, b) for(int i = a; i <= b; ++i)
#define per(i, a, b) for(int i = a; i >= b; --i)
#define I inline
const int N = 4e5 + 5, mod = 998244353;
I int read(){
rg char ch = gc();
rg int f = 0;
rg long long x = 0;
while(!isdigit(ch)) f |= (ch == '-'), ch = gc();
while(isdigit(ch)) x = ((x << 1) + (x << 3) + (ch ^ 48)) % mod, ch = gc();
return f ? mod - x : x;
}
I int ksm(int a, int b){
int ans = 1;
while(b){ if(b & 1) ans = 1ll * a * ans % mod; b >>= 1; a = 1ll * a * a % mod; }
return ans;
}
int G = 3, Gn = ksm(G, mod - 2);
int f[N], g[N], n, k;
I void fwt_or(int *f, int lim, int flag){
for(int l = 2; l <= lim; l <<= 1)
for(int m = l >> 1, j = 0; j < lim; j += l)
for(int i = j; i < j + m; ++i)
(f[j + m] += flag * f[j]) %= mod;
}
I void fwt_and(int *f, int lim, int flag){
for(int l = 2; l <= lim; l <<= 1)
for(int m = l >> 1, j = 0; j < lim; j += l)
for(int i = j; i < j + m; ++i)
(f[j] += flag * f[j + m]) %= mod;
}
const int inv2 = ksm(2, mod - 2);
I void fwt_xor(int *f, int lim, int flag){
for(int l = 2; l <= lim; l <<= 1)
for(int m = l >> 1, j = 0; j < lim; j += l)
for(int i = j; i < j + m; ++i){
int x = f[i], y = f[i + m];
f[i] = (x + y) % mod; f[i + m] = (x + mod - y) % mod;
if(flag == -1){
f[i] = 1ll * f[i] * inv2 % mod; f[i + m] = 1ll * f[i] * inv2 % mod;
}
}
}
I int get_phi(int x){
int len = sqrt(x);
int res = 1;
rep(i, 2, len){
if(!(x % i)){
x /= i;
res = 1ll * res * (i - 1) % mod;
while(!(x % i)) x /= i, res = 1ll * res * i % mod;
}
}
if(x != 1) res = 1ll * res * (x - 1) % mod;
return res;
}
I int find_root(int x){
int phi = get_phi(x), p = phi;
int len = sqrt(phi);
static int s[N], cnt;
cnt = 0;
rep(i, 2, len){
if(!(p % i)){
p /= i;
s[++cnt] = i;
while(!(p % i)) p /= i;
}
}
if(p != 1) s[++cnt] = p;
rep(i, 1, cnt) cout << s[i] << " "; cout << endl;
cout << phi << endl;
rep(i, 2, mod - 1){
int flag = 0;
rep(j, 1, cnt) if(ksm(i, phi / s[j]) == 1){ flag = 1; break; }
if(!flag) return i;
}
}
int fac[N], ifac[N];
I void get_fac(int n){
fac[0] = ifac[0] = 1;
rep(i, 1, n){
fac[i] = 1ll * fac[i - 1] * i % mod;
ifac[i] = 1ll * ifac[i - 1] * fac[i] % mod;
}
int inv = ksm(ifac[n], mod - 2);
per(i, n, 1){
ifac[i] = 1ll * ifac[i - 1] * inv % mod;
inv = 1ll * fac[i] * inv % mod;
}
}
struct FFT{
int A[N], B[N], c[N], b2[N], bb[N], ib2[N], sa[N], rev[N];
I void NTT(int *a, int lim, int len, int flag){
rep(i, 1, lim - 1) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (len - 1));
rep(i, 1, lim - 1) if(rev[i] > i) swap(a[i], a[rev[i]]);
for(int l = 2; l <= lim; l <<= 1){
const int m = l >> 1, Gi = ksm(flag == 1 ? G : Gn, (mod + 1) / l);
for(int j = 0; j < lim; j += l){
int g = 1;
for(int i = j; i < j + m; ++i, g = 1ll * g * Gi % mod){
int x = a[i], y = 1ll * g * a[i + m] % mod;
a[i] = (x + y) % mod;
a[i + m] = (x + mod - y) % mod;
}
}
}
}
I void mul(int *a, int *b, int na, int nb, int *c){
int lim = 1, len = 0;
while(lim <= na + nb) lim <<= 1, ++len;
memcpy(A, a, (na + 1) * sizeof(int)); memcpy(B, b, (nb + 1) * sizeof(int));
fill(A + na + 1, A + lim, 0); fill(B + nb + 1, B + lim, 0);
NTT(A, lim, len, 1); NTT(B, lim, len, 1);
rep(i, 0, lim - 1) A[i] = 1ll * A[i] * B[i] % mod;
NTT(A, lim, len, -1);
const int inv = ksm(lim, mod - 2);
rep(i, 0, na + nb) c[i] = 1ll * A[i] * inv % mod;
fill(c + na + nb + 1, c + lim, 0);
}
I void ni_ab(int *a, int *b, int n){
int lim = 1;
while(lim <= n) lim <<= 1;
b[0] = ksm(a[0], mod - 2);
for(int xmod = 1, nlen = 2; xmod < lim; xmod <<= 1, ++nlen){
int nlim = xmod << 2;
memcpy(B, b, xmod * sizeof(int)); memcpy(A, a, (xmod << 1) * sizeof(int));
fill(B + xmod, B + nlim, 0); fill(A + (xmod << 1), A + nlim, 0);
NTT(A, nlim, nlen, 1); NTT(B, nlim, nlen, 1);
rep(i, 0, nlim - 1) A[i] = ((B[i] << 1) % mod + mod - 1ll * A[i] * B[i] % mod * B[i] % mod) % mod;
NTT(A, nlim, nlen, -1);
const int inv = ksm(nlim, mod - 2);
rep(i, 0, (xmod << 1) - 1) b[i] = 1ll * A[i] * inv % mod;
}
fill(b + n + 1, b + lim, 0);
}
I void ln(int *a, int *b, int n){//bb b2
int lim = 1, len = 0;
while(lim <= n) lim <<= 1, ++len;
rep(i, 0, n - 1) bb[i] = 1ll * (i + 1) * a[i + 1] % mod;
ni_ab(a, b2, n);
mul(bb, b2, n - 1, n, b);
per(i, n, 1) b[i] = 1ll * b[i - 1] * ksm(i, mod - 2) % mod; b[0] = 0;
}
I void sqrt(int *a, int *b, int n){//bb ib2 b2
int lim = 1, len = 0;
while(lim <= n) lim <<= 1, ++len;
fill(b, b + lim, 0); fill(b2, b2 + lim, 0);
b[0] = 1;
for(int xmod = 1; xmod < lim; xmod <<= 1){
rep(i, 0, xmod - 1) b2[i] = (b[i] << 1) % mod; ni_ab(b2, ib2, (xmod << 1) - 1);
mul(b, b, xmod - 1, xmod - 1, bb);
rep(i, 0, (xmod << 1) - 1) bb[i] = (bb[i] + a[i]) % mod;
mul(bb, ib2, (xmod << 1) - 1, (xmod << 1) - 1, b);
}
}
I void exp(int *a, int *b, int n){//ib2 bb b2 c
int lim = 1; while(lim <= n) lim <<= 1;
fill(b, b + lim, 0);
b[0] = 1;
for(int xmod = 1; xmod < lim; xmod <<= 1){
ln(b, ib2, (xmod << 1) - 1); //ib2 = ln(b);
rep(i, 0, (xmod << 1) - 1) c[i] = (a[i] + mod - ib2[i]) % mod;
c[0] = (c[0] + 1) % mod;
mul(b, c, (xmod << 1) - 1, (xmod << 1) - 1, b);
}
fill(b + n + 1, b + lim, 0);
}
I void pow(int *a, int *b, int n, int _k = k){
int lim = 1; while(lim <= n) lim <<= 1;
ln(a, sa, n);
rep(i, 0, n) sa[i] = 1ll * sa[i] * _k % mod;
exp(sa, b, n);
}
}T;
signed main(){
n = read() - 1; k = read();
rep(i, 0, n) g[i] = read();
T.pow(g, f, n);
rep(i, 0, n) printf("%d ", f[i]);
return 0;
}