预备知识：FFT/NTT
-
多项式的逆
给定一个多项式 F(x)F(x)F(x)，请求出一个多项式 G(x)G(x)G(x)，满足 F(x)∗G(x)≡1(mod xn)F(x)*G(x) \equiv 1(mod\ x^n)F(x)∗G(x)≡1(mod xn)。
系数对 998244353998244353998244353 取模，1≤n≤1051≤n≤10^51≤n≤105首先将多项式的长度拓展至222的次幂，然后我们要求的是
G(x)∗F(x)≡1 (mod xn)G(x)*F(x) \equiv 1\ (mod\ x^n)G(x)∗F(x)≡1 (mod xn)假设已经求出了
H(x)∗F(x)≡1 (mod x⌈n2⌉)H(x)*F(x) \equiv 1\ (mod\ x^{⌈\frac n2⌉})H(x)∗F(x)≡1 (mod x⌈2n⌉)又因为有
G(x)∗F(x)≡1 (mod x⌈n2⌉)G(x)*F(x) \equiv 1\ (mod\ x^{⌈\frac n2⌉})G(x)∗F(x)≡1 (mod x⌈2n⌉)两式相减有
(G(x)−H(x))∗F(x)≡0 (mod x⌈n2⌉)(G(x)-H(x))*F(x) \equiv 0\ (mod\ x^{⌈\frac n2⌉})(G(x)−H(x))∗F(x)≡0 (mod x⌈2n⌉)G(x)−H(x)≡0 (mod x⌈n2⌉)G(x)-H(x) \equiv 0\ (mod\ x^{⌈\frac n2⌉})G(x)−H(x)≡0 (mod x⌈2n⌉)因为左边得到的多项式前⌈n2⌉{⌈\frac n2⌉}⌈2n⌉项都是000，平方前nnn项都是000，所以有
(G(x)−H(x))2≡0 (mod xn)(G(x)-H(x))^2 \equiv 0\ (mod\ x^{n})(G(x)−H(x))2≡0 (mod xn)G2(x)−2G(x)H(x)+H2(x)≡0 (mod xn)G^2(x)-2G(x)H(x)+H^2(x)\equiv 0\ (mod\ x^n)G2(x)−2G(x)H(x)+H2(x)≡0 (mod xn)两边同时乘以F(x)F(x)F(x)
G(x)−2H(x)+H2(x)F(x)≡0 (mod xn)G(x)-2H(x)+H^2(x)F(x)\equiv 0\ (mod\ x^n)G(x)−2H(x)+H2(x)F(x)≡0 (mod xn)G(x)≡2H(x)−H2(x)F(x) (mod xn)G(x)\equiv2H(x)-H^2(x)F(x)\ (mod\ x^n)G(x)≡2H(x)−H2(x)F(x) (mod xn)
我们就得到了递推式，时间复杂度如下T(n)=T(n/2)+O(nlog2n)=O(nlog2n)T(n)=T(n/2)+O(nlog_2n)=O(nlog_2n)T(n)=T(n/2)+O(nlog2n)=O(nlog2n) -
多项式开根
给定一个 n−1n−1n−1 次多项式 F(x)F(x)F(x)，求一个在 mod xnmod\ x^nmod xn 意义下的多项式 G(x)G(x)G(x)，使得 G2(x)≡F(x) (mod xn)G^2(x) \equiv F(x) \ (mod\ x^n)G2(x)≡F(x) (mod xn)（多项式的系数在 mod 998244353mod\ 998244353mod 998244353 的意义下进行运算，1≤n≤1051≤n≤10^51≤n≤105）
同样假设求出了H2(x)≡F(x) (mod x⌈n2⌉)H^2(x)\equiv F(x)\ (mod\ x^{⌈\frac n2⌉})H2(x)≡F(x) (mod x⌈2n⌉)且有
G2(x)≡F(x) (mod x⌈n2⌉)G^2(x)\equiv F(x)\ (mod\ x^{⌈\frac n2⌉})G2(x)≡F(x) (mod x⌈2n⌉)所以
G(x)−H(x)≡0 (mod x⌈n2⌉)G(x)-H(x) \equiv 0\ (mod\ x^{⌈\frac n2⌉})G(x)−H(x)≡0 (mod x⌈2n⌉)(G(x)−H(x))2≡0 (mod xn)(G(x)-H(x))^2 \equiv 0\ (mod\ x^{n})(G(x)−H(x))2≡0 (mod xn)G2(x)−2G(x)H(x)+H2(x)≡0 (mod xn)G^2(x)-2G(x)H(x)+H^2(x)\equiv 0\ (mod\ x^n)G2(x)−2G(x)H(x)+H2(x)≡0 (mod xn)F(x)−2G(x)H(x)+H2(x)≡0 (mod xn)F(x)-2G(x)H(x)+H^2(x)\equiv 0\ (mod\ x^n)F(x)−2G(x)H(x)+H2(x)≡0 (mod xn)
所以2G(x)≡F(x)H(x)+H(x) (mod xn)2G(x)\equiv \frac{F(x)}{H(x)}+H(x)\ (mod\ x^n)2G(x)≡H(x)F(x)+H(x) (mod xn)
此处除法转化为乘上多项式的逆，就能递推了，时间复杂度也是O(nlog2n)O(nlog_2n)O(nlog2n) -
多项式求自然对数
给出 n−1n−1n−1 次多项式 F(x)F(x)F(x)，求一个 mod xnmod\ x^nmod xn 下的多项式 G(x)G(x)G(x)，满足 G(x)≡lnF(x) (mod xn)G(x) \equiv \ln F(x)\ (mod\ x^n)G(x)≡lnF(x) (mod xn).
多项式的系数在 mod 998244353mod\ 998244353mod 998244353 的意义下进行运算，1≤n≤1051≤n≤10^51≤n≤105
G(x)=ln F(x)=∫F′(x)F(x)dx=∫F′(x)F−1(x)dx\begin{aligned}
G(x)&=ln\ F(x)\\
&=\int \frac{F'(x)}{F(x)}dx\\
&=\int F'(x)F^{-1}(x)dx
\end{aligned}G(x)=ln F(x)=∫F(x)F′(x)dx=∫F′(x)F−1(x)dx所以G=∫F′F−1G=\int F'F^{-1}G=∫F′F−1，直接多项式求逆+多项式求导积分就行了
求逆复杂度O(nlog2n)O(nlog_2n)O(nlog2n)，求导积分复杂度O(n)O(n)O(n)，总时间复杂度仍是O(nlog2n)O(nlog_2n)O(nlog2n) -
多项式exp
多项式exp是用牛顿迭代做的,我也不会证明牛顿迭代,只会背公式.
本来牛顿迭代公式是这个样子的xn+1=xn−f(xn)f′(xn)x_{n+1}=x_n-\frac{f(x_n)}{f'(x_n)}xn+1=xn−f′(xn)f(xn)
对于多项式的牛顿迭代长这样un+1=un−f(u,v)fu(u,v)u_{n+1}=u_n-\frac{f(u,v)}{f_u(u,v)}un+1=un−fu(u,v)f(u,v)- 其中f(u,v)f(u,v)f(u,v)表示已知多项式为vvv,未知多项式为uuu的多项式方程,就拿exp为例,原方程为u≡ev (mod xn)u\equiv e^v\ (mod\ x^n)u≡ev (mod xn)两边取对数为ln u≡v (mod xn)ln u−v≡0 (mod xn)\begin{aligned}ln\ u\equiv v\ (mod\ x^n)\\ln\ u-v\equiv 0\ (mod\ x^n)\end{aligned}ln u≡v (mod xn)ln u−v≡0 (mod xn)那么f(u,v)=ln u−vf(u,v)=ln\ u-vf(u,v)=ln u−v.
- fu(u,v)f_u(u,v)fu(u,v)表示的是f(u,v)f(u,v)f(u,v)对uuu取导.注意是对uuu,不是对多项式里的xxx.
举个栗子,ln vln\ vln v对vvv求导是1v\frac 1vv1,但是对xxx求导是v′v\frac{v'}{v}vv′(此处v′v'v′是对xxx求导)
那么我们就可以推推柿子un+1≡un−(ln un−v(ln un−v)′) (mod xn)≡un−(ln un−v1un) (mod xn)≡un−(ln un−v)un (mod xn)≡un(1−ln un+v) (mod xn)\begin{aligned}u_{n+1}&\equiv u_n-\left(\frac{ln\ u_n-v}{(ln\ u_n-v)'}\right)\ &(mod\ x^n)\\&\equiv u_n-\left(\frac{ln\ u_n-v}{\frac1{u_n}}\right)\ &(mod\ x^n)\\&\equiv u_n-(ln\ u_n-v)u_n\ &(mod\ x^n)\\&\equiv u_n(1-ln\ u_n+v)\ &(mod\ x^n)\end{aligned}un+1≡un−((ln un−v)′ln un−v) ≡un−(un1ln un−v) ≡un−(ln un−v)un ≡un(1−ln un+v) (mod xn)(mod xn)(mod xn)(mod xn)然后就套用前面的模板迭代做了
TIP:TIP:TIP:其实多项式开根也可以用这个方法推出来,有兴趣可以去试试看.
代码
#include <cmath>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;
namespace READ {
inline void read(int &num) {
char ch; int flg=1;
while((ch=getchar())<'0'||ch>'9')if(ch=='-')flg=-flg;
for(num=0;ch>='0'&&ch<='9';num=num*10+ch-'0',ch=getchar());
num*=flg;
}
}
using namespace READ;
const int mod = 998244353, G = 3, inv2 = 499122177;
const int MAXN = 1<<18; // 注意MAXN要大于2*n
int rev[MAXN], Wn[MAXN+1][2], inv[MAXN+1];
inline void reset(int len) { for(int i = 0; i < len; ++i) rev[i] = (rev[i>>1]>>1)|((i&1)?(len>>1):0); }
inline int qmul(int a, int b) {
int res = 1;
while(b) {
if(b&1) res = 1ll * res * a % mod;
a = 1ll * a * a % mod; b>>=1;
}
return res;
}
inline int add(int x, int y) { return x+y>mod ? x+y-mod : x+y; }
namespace Poly {
int A[MAXN],B[MAXN],C[MAXN],D[MAXN],E[MAXN],F[MAXN],GG[MAXN],H[MAXN];//每个不同函数开不同的临时数组,否则可能重复使用
inline void clear(int *b, int n) { for(int i = 0; i < n; ++i) b[i] = 0; } //多清零
inline void NTT(int *arr, int n, int flg) {
int wn, w, A0, wA1, i, j, k;
for(i = 0; i < n; ++i) if(i < rev[i]) swap(arr[i], arr[rev[i]]);
for(i = 2; i <= n; i<<=1) {
wn = Wn[i][(~flg)?0:1];
for(j = 0, w = 1; j < n; j += i, w = 1)
for(k = j; k < j + i/2; ++k, w = 1ll * w * wn % mod) {
A0 = arr[k], wA1 = 1ll * w * arr[k + i/2] % mod;
arr[k] = add(A0, wA1);
arr[k + i/2] = add(A0, -wA1+mod);
}
}
if(!(~flg)) for(i = 0; i < n; ++i) arr[i] = 1ll * arr[i] * inv[n] % mod;
}
inline void INV(int *a, int *b, int n) {
clear(b, n<<1), b[0] = qmul(a[0], mod-2);
clear(A, n<<1), clear(B, n<<1);
int len, lim;
for(len = 2; len < (n<<1); len<<=1) {
lim = len<<1;
for(int i = 0; i < len; ++i) A[i] = a[i], B[i] = b[i];
reset(lim); NTT(A, lim, 1); NTT(B, lim, 1);
for(int i = 0; i < lim; ++i) b[i] = ((2ll - 1ll * A[i] * B[i] % mod) * B[i] % mod + mod) % mod;
NTT(b, lim, -1);
for(int i = len; i < lim; ++i) b[i] = 0;
}
for(int i = 0; i < len; ++i) A[i] = B[i] = 0;
for(int i = n; i < len; ++i) b[i] = 0;
}
inline void SQRT(int *a, int *b, int n) {
clear(b, n<<1); b[0] = int(sqrt(a[0])+0.5);
int len, lim, *A = GG, *B = H;
clear(A, n<<1), clear(B, n<<1);
for(len = 2; len < (n<<1); len<<=1) {
lim = len<<1;
for(int i = 0; i < len; ++i) A[i] = a[i];
INV(b, B, len);
reset(lim); NTT(A, lim, 1); NTT(B, lim, 1);
for(int i = 0; i < lim; ++i) A[i] = 1ll * A[i] * B[i] % mod;
NTT(A, lim, -1);
for(int i = 0; i < len; ++i) b[i] = 1ll * (b[i] + A[i]) % mod * inv2 % mod;
for(int i = len; i < lim; ++i) b[i] = 0;
}
for(int i = 0; i < len; ++i) A[i] = B[i] = 0;
for(int i = n; i < len; ++i) b[i] = 0;
}
inline void INT(int *a, int *b, int n) {
clear(b, n<<1);
for(int i = n-1; i; --i)
b[i] = 1ll * a[i-1] * qmul(i, mod-2) % mod;
b[0] = 0;
}
inline void DIF(int *a, int *b, int n) {
clear(b, n<<1);
for(int i = 1; i < n; ++i)
b[i-1] = 1ll * a[i] * i % mod;
b[n-1] = 0;
}
inline void LN(int *a, int *b, int n) {
int *A = E, *B = F, len = 1;
clear(A, n<<1), clear(B, n<<1);
INV(a, A, n); DIF(a, B, n);
while(len < (n<<1)) len<<=1;
reset(len); NTT(A, len, 1); NTT(B, len, 1);
for(int i = 0; i < len; ++i) A[i] = 1ll * A[i] * B[i] % mod;
NTT(A, len, -1); INT(A, b, n);
}
inline void EXP(int *a, int *b, int n) {
clear(b, n<<1), b[0] = 1;
int len, lim, *A = C, *B = D;
clear(A, n<<1), clear(B, n<<1);
for(len = 2; len < (n<<1); len<<=1) {
lim = len<<1;
LN(b, A, len);
for(int i = 0; i < len; ++i) A[i] = ((i==0)-A[i]+a[i]+mod) % mod, B[i] = b[i];
for(int i = len; i < lim; ++i) A[i] = B[i] = 0;
reset(lim); NTT(A, lim, 1); NTT(B, lim, 1);
for(int i = 0; i < lim; ++i) A[i] = 1ll * A[i] * B[i] % mod;
NTT(A, lim, -1);
for(int i = 0; i < len; ++i) b[i] = A[i];
}
for(int i = n; i < len; ++i) b[i] = 0;
}
}
using namespace Poly;
inline void Pre() { //预处理原根及逆元
for(int i = 1; i <= MAXN; i<<=1)
Wn[i][0] = qmul(G, (mod-1)/i),
Wn[i][1] = qmul(G, (mod-1)-(mod-1)/i),
inv[i] = qmul(i, mod-2);
}
int main () {
Pre();
}
Upd:Upd:Upd:现在看,发现我以前sb了…不用开不同的数组来取地址,直接static定义就行了…
带余除法（待更…）
多项式多点求值（待更…）