多项式操作
1 多项式求逆
求 对于 \(A(x)\) 求 \(B(x)\) 使得 \(A(x) * B(x) \equiv 1 \ (mod\ n)\)
\[AB \equiv 1 (mod\ n)
\]
\[AB - 1\equiv 0(mod\ n)
\]
\[(AB - 1) ^ 2 \equiv 0 (mod\ 2n)
\]
因为
\[F(x) = \sum\limits_{i = 0}^ {n - 1} \sum\limits_{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\) 答案就是它的根
还是牛顿迭代一下
\[\begin{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;
}
