FFT杂谈
因为之前的我已经不想去动他了。QAQ
其实FFT并不难学，大纲是这样：
- 前置知识
- 分治
- 蝴蝶变换
- 学完了，没了
之前那个前面不是是抄你谷多项式奆佬的，后面的各种优化没什么实际意义还疯狂掉精度，基本在实战中不会去碰他，所以那篇文章我可以说是弃飞了。
FFT的前置知识并不多，向量和单位根，向量相乘的两个性质一个用相似，一个用两点距离公式爆算都能证明，也就没什么难度。
会用到的性质：
- 向量相乘模相乘
- 向量相乘辐角相加
单位根的定义百度一下，你就知道。
当然，单位根用到的性质在下面列出来了：
- w n 0 = 1 , w n k = w n k m o d n w_{n}^0=1,w_{n}^k=w_{n}^{k\mod n} wn0=1,wnk=wnkmodn
- { x ∈ C x∈C x∈C| x = w n i ( i ∈ [ 0 , n − 1 ] , i ∈ Z ) x=w_{n}^i(i∈[0,n-1],i∈Z) x=wni(i∈[0,n−1],i∈Z)}的元素个数是 n − 1 n-1 n−1个（其实意思就是 w n i w_{n}^i wni互不相等）
- w n k = w 2 n 2 k w_{n}^k=w_{2n}^{2k} wnk=w2n2k
- w n k = − w n k + n 2 ( n ≡ 0 m o d 2 ) w_{n}^k=-w_{n}^{k+\frac{n}{2}}(n≡0\mod2) wnk=−wnk+2n(n≡0mod2)
- ∑ i = 0 n − 1 ( w n j − k ) i = n ∗ [ j = k ] \sum\limits_{i=0}^{n-1}(w_{n}^{j-k})^i=n*[j=k] i=0∑n−1(wnj−k)i=n∗[j=k]（ [ A ] [A] [A]表示如果 A A A条件成立，则 [ A ] = 1 [A]=1 [A]=1，否则为 0 0 0）
- w n k ∗ w n t = w n k + t w_{n}^k*w_{n}^t=w_{n}^{k+t} wnk∗wnt=wnk+t
只要用到这几个性质，就可以完美的推导出FFT。
这里也不展开，主要就是奇偶分治。
至于蝴蝶迭代，这个我之前的博客也就，在这里就不展开了，代码以前也是有的，我看了一下，我之前的代码长得还行(⊙﹏⊙)。
但是如果这里就这么一点东西，那我还做个P的杂谈啊(╯‵□′)╯︵┻━┻
首先，关于快速傅里叶变换，你可以觉得其很神奇，因为他巧妙的利用了单位根的性质，之前那个IDFT的推导，也很神奇，好像一气呵成似的。
当然，这里介绍一种不同的理解方式（听说卷积跟线性代数有关，但是那个我不会，就不展开了）
其实你会发现，其实DFT和IDFT就是一个 1 ∗ n 1*n 1∗n的矩阵乘一个 n ∗ n n*n n∗n的矩阵的过程，所以理论上只要 I D F T IDFT IDFT和 D F T DFT DFT的矩阵互逆，这样就可以完美的相互抵消。
当然，DFT中的矩阵叫范德蒙德矩阵。
图片来自参考链接的奆佬博客Orz中
所以，这看似运气十足的FFT背后也有其科学的数学道理。
NTT
全称：快速速论变换
请注意，本文中对于一个多项式模一个整数，只是对每一个系数取模而已。
优点：
- 全是整数，不用担心掉精度的问题。
- 可以计算特定模数的卷积。
- 因为不用计算复数，所以乘法少了不少，貌似常数加快了b（￣▽￣）d 。
缺点：
- 不能计算小数，过程涉及到取模，所以一旦数域超过了模数，就难以计算。
- 模数必须是 k ∗ 2 t + 1 k*2^t+1 k∗2t+1的形式，且必须是一个质数。
- 由于模运算比较慢，所以常数经常不如FFT。
你是不是以为NTT十分的妙？
不，其实就是FFT的过程然后把单位根换了一个更加神奇的东西。
上面也列了，我们会利用到单位根的一些性质，但是其实只要我们能在整数中找到一个同样能满足上述性质的东西，就可以完美的换成整数啦。
当然，在一般的整数域我们很难找到满足条件的数字。
但是如果套上一个模，如果这个模有原根，则有个与原根相关的东西刚好可以满足上面。
设模为 p p p， p p p是形式为 k ∗ 2 t + 1 k*2^t+1 k∗2t+1形式的质数，同时原根为 g g g。
（原根定义：在模 p p p的时候，如果对于 n n n， φ ( p ) φ(p) φ(p)是最小的正整数 t t t 满足 n t ≡ 1 n^t≡1 nt≡1，则 n n n 为 p 的 原 根 p的原根 p的原根， φ ( p ) φ(p) φ(p)表示小于 p p p的正整数中与 p p p 互质的数字个数）
则我们可以把 w n i w_{n}^i wni替换成 g p − 1 n ∗ i g^{\frac{p-1}{n}*i} gnp−1∗i（假设 n ∣ p − 1 n|p-1 n∣p−1），用 g n i g_{n}^i gni简化符号。
当然，其是否满足上面的式子，我们来看一下：
-
g
n
0
≡
1
,
g
n
i
≡
g
n
i
m
o
d
n
g_{n}^0≡1,g_{n}^i≡g_{n}^{i\mod n}
gn0≡1,gni≡gnimodn
证明： g n 0 ≡ g 0 = 1 g_{n}^0≡g^0=1 gn0≡g0=1，将 i i i 表示成 k n + t ( 0 ≤ t < n ) kn+t(0≤t<n) kn+t(0≤t<n) 的形式，所以 g n i ≡ g n k n + t ≡ g p − 1 n ∗ ( k n + t ) ≡ g p − 1 n ∗ k n ∗ g p − 1 n ∗ t ≡ g p − 1 n ∗ t ≡ g n t g_{n}^i≡g_{n}^{kn+t}≡g^{\frac{p-1}{n}*(kn+t)}≡g^{\frac{p-1}{n}*kn}*g^{\frac{p-1}{n}*t}≡g^{\frac{p-1}{n}*t}≡g_{n}^t gni≡gnkn+t≡gnp−1∗(kn+t)≡gnp−1∗kn∗gnp−1∗t≡gnp−1∗t≡gnt，证毕，同理，因为都等于 g n t g_{n}^t gnt，所以 g n i ≡ g n i ± n g_{n}^i≡g_{n}^{i±n} gni≡gni±n - {
x
∈
Z
x∈Z
x∈Z|
x
=
g
n
i
(
i
∈
[
0
,
n
−
1
]
,
i
∈
Z
)
x=g_{n}^i(i∈[0,n-1],i∈Z)
x=gni(i∈[0,n−1],i∈Z)}的元素个数是
n
−
1
n-1
n−1个（意思就是
g
n
i
g_{n}^i
gni互不相等）
证明：假设 g n i ≡ g n j ( i > j ) g_{n}^i≡g_{n}^j(i>j) gni≡gnj(i>j)
g p − 1 n ∗ i ≡ g p − 1 n ∗ j g^{\frac{p-1}{n}*i}≡g^{\frac{p-1}{n}}*j gnp−1∗i≡gnp−1∗j
( g p − 1 n ∗ ( i − j ) − 1 ) ∗ g p − 1 n ∗ j ≡ 0 (g^{\frac{p-1}{n}*(i-j)}-1)*g^{\frac{p-1}{n}*j}≡0 (gnp−1∗(i−j)−1)∗gnp−1∗j≡0
那么 p ∣ ( g p − 1 n ∗ ( i − j ) − 1 ) ∗ g p − 1 n ∗ j p|(g^{\frac{p-1}{n}*(i-j)}-1)*g^{\frac{p-1}{n}*j} p∣(gnp−1∗(i−j)−1)∗gnp−1∗j，但是如果 p ∤ g p∤g p∤g，则 p ∤ p k ( k ∈ Z ∗ ) p∤p^k(k∈Z^{*}) p∤pk(k∈Z∗)
所以 p ∣ ( g p − 1 n ∗ ( i − j ) − 1 ) p|(g^{\frac{p-1}{n}*(i-j)}-1) p∣(gnp−1∗(i−j)−1)
所以 g p − 1 n ∗ ( i − j ) ≡ 1 g^{\frac{p-1}{n}*(i-j)}≡1 gnp−1∗(i−j)≡1，这显然是不成立的，因为 p − 1 n ∗ ( i − j ) < p − 1 = φ ( p ) \frac{p-1}{n}*(i-j)<p-1=φ(p) np−1∗(i−j)<p−1=φ(p) -
g
n
i
≡
g
2
n
2
i
g_{n}^i≡g_{2n}^{2i}
gni≡g2n2i
证明： g n i ≡ g p − 1 n ∗ i ≡ g p − 1 2 n ∗ 2 i g_{n}^i≡g^{\frac{p-1}{n}*i}≡g^{\frac{p-1}{2n}*2i} gni≡gnp−1∗i≡g2np−1∗2i。 -
g
n
i
≡
−
g
n
i
+
n
2
(
2
∣
n
)
g_{n}^i≡-g_{n}^{i+\frac{n}{2}}(2|n)
gni≡−gni+2n(2∣n)
证明： g n i + n 2 ≡ g p − 1 n ∗ ( i + n 2 ) ≡ g p − 1 n ∗ i ∗ g p − 1 n ∗ n 2 ≡ g n i ∗ g p − 1 2 = − g n i g_{n}^{i+\frac{n}{2}}≡g^{\frac{p-1}{n}*(i+\frac{n}{2})}≡g^{\frac{p-1}{n}*i}*g^{\frac{p-1}{n}*\frac{n}{2}}≡g_{n}^i*g^{\frac{p-1}{2}}=-g_{n}^i gni+2n≡gnp−1∗(i+2n)≡gnp−1∗i∗gnp−1∗2n≡gni∗g2p−1=−gni
当然，至于为什么 g p − 1 2 ≡ − 1 g^{\frac{p-1}{2}}≡-1 g2p−1≡−1，我们还需要证明一个东西：
对于模质数 p p p， x 2 ≡ 1 x^2≡1 x2≡1的解只有两个 1 , − 1 1,-1 1,−1（或者说是其的同余）。
x 2 ≡ 1 x^2≡1 x2≡1
则 ( x + 1 ) ( x − 1 ) ≡ 0 (x+1)(x-1)≡0 (x+1)(x−1)≡0，那么 p ∣ ( x + 1 ) p|(x+1) p∣(x+1)或者 p ∣ ( x − 1 ) p|(x-1) p∣(x−1)，所以 x ≡ p + 1 ≡ 1 x≡p+1≡1 x≡p+1≡1或 x = p − 1 ≡ − 1 x=p-1≡-1 x=p−1≡−1。
但是因为 g g g是原根，所以 g p − 1 2 g^{\frac{p-1}{2}} g2p−1只能是 − 1 -1 −1。 -
g
n
i
∗
g
n
j
≡
g
n
i
+
j
g_{n}^{i}*g_{n}^j≡g^{i+j}_{n}
gni∗gnj≡gni+j。
5,6点换了一个顺序，海星。
证明：展开即可，没什么好说的 -
∑
i
=
0
n
−
1
(
g
n
j
−
k
)
i
≡
n
∗
[
j
=
k
]
\sum\limits_{i=0}^{n-1}(g_{n}^{j-k})^i≡n*[j=k]
i=0∑n−1(gnj−k)i≡n∗[j=k]
证明：其实完全没有什么必要证明的，反正这个东西只要用相同的性质套到FFT中便可以证明了：
先求 ∑ i = 0 n − 1 ( g n j ) i \sum\limits_{i=0}^{n-1}(g_{n}^j)i i=0∑n−1(gnj)i
当 j ≠ 0 j≠0 j=0时
∑ i = 0 n − 1 ( g n j ) i ≡ g n n j − 1 g n j − 1 ≡ 0 \sum\limits_{i=0}^{n-1}(g_{n}^j)i≡\frac{g_{n}^{nj}-1}{g_{n}^j-1}≡0 i=0∑n−1(gnj)i≡gnj−1gnnj−1≡0。
但是，如果 j = 0 j=0 j=0
那么显然 ∑ i = 0 n − 1 ( g n j ) i ≡ n \sum\limits_{i=0}^{n-1}(g_{n}^j)i≡n i=0∑n−1(gnj)i≡n
那么什么时候 j = 0 j=0 j=0呢？也就是上述式子中的 j = k j=k j=k时成立，证毕。
然后就没有什么了，直接把代码中的单位根换成 g g g，就可以求解了，但是，仅适用于整数卷积。
当然，需要注意的是：
- g n i g_{n}^i gni中 n n n必须整除 p − 1 p-1 p−1，而 F F T FFT FFT中 n = 2 i ( i ∈ [ 1 , r ] ) n=2^i(i∈[1,r]) n=2i(i∈[1,r])，所以 p = k ∗ 2 t + 1 p=k*2^{t}+1 p=k∗2t+1就是这么一个道理。
- 同样的，对于补充到 2 2 2的次幂的多项式长度 l e n ′ = 2 r len'=2^r len′=2r，我们要求 r < = t r<=t r<=t，所以对于原来的多项式长度 l e n len len，我们要求 l e n ≤ 2 t len≤2^t len≤2t，所以 p p p的 2 t 2^t 2t要尽可能的大，那是最好的。
- 其次，值域的上届减下界 + 1 +1 +1要小于 p p p，因为模 p p p最怕的就是在值域中有两个数字同余。
- 对于
g
n
i
(
i
<
0
)
g_{n}^i(i<0)
gni(i<0)，此时
g
n
i
≡
1
g
p
−
1
n
∗
(
−
i
)
≡
(
1
g
)
p
−
1
n
∗
(
−
i
)
g_{n}^i≡\frac{1}{g^{\frac{p-1}{n}*(-i)}}≡(\frac{1}{g})^{\frac{p-1}{n}*(-i)}
gni≡gnp−1∗(−i)1≡(g1)np−1∗(−i)，所以只要求出
g
g
g的逆就可以计算
g
n
i
(
i
<
0
)
g_{n}^i(i<0)
gni(i<0)了，当然，
g
g
g的逆也是一个原根，证明如下：
g ′ i ≡ ( g i ) ′ g'^{i}≡(g^i)' g′i≡(gi)′（其中 x ′ x' x′表示 x x x的逆）
而 1 1 1的逆是其本身，且因为一个数不存在两个不同余的逆，所以很容易知道 g ′ g' g′也是原根。 - 在IDFT中应当用 g n − i g_{n}^{-i} gn−i替换 w n i w_{n}^i wni，应该不会只有我一个憨憨用 − g n i -g_{n}^i −gni去代替吧QAQ。
- 在最后输出除长度 n n n时，需要注意的是，不能直接除 p p p，即使值域中的每个数字都小于 p p p，因为虽然值域的长度小于 p p p，但是不代表要求值域的长度 ∗ n *n ∗n以后还小于 p p p，所以还是乘 n ′ n' n′ 的更加保险，当然，如果你这值域长度小的可怜，当我没说。
- 没什么好提示了，还不快去🐎代码？
当然，这里提一下常用的一个模数： 998244353 = 119 ∗ 2 23 + 1 998244353=119*2^{23}+1 998244353=119∗223+1， 2 23 = 8388608 2^{23}=8388608 223=8388608，他的原根是 3 3 3。
一般模数是这个就不用担心长度不够的问题了，因为如果长度还大于 8388608 8388608 8388608，那多半是你做法的问题。
时间复杂度： O ( n log n ) O(n\log n) O(nlogn)
#include<cstdio>
#include<cstring>
#include<cstdlib>
#define N 3100000
using namespace std;
typedef long long LL;
const LL mod=998244353;
inline LL ksm(LL x,LL y)
{
LL ans=1;
while(y)
{
if(y&1)ans=(ans*x)%mod;
x=x*x%mod;y>>=1;
}
return ans;
}
template<class T>
inline void zwap(T &x,T &y){x^=y^=x^=y;}
LL n,m;
LL limit,l;
LL a[N],b[N],r[N]/*蝴蝶变换*/;//用来相乘的多项式
LL g[N],g1=3,g2=ksm(3,mod-2),gg;
void NTT(LL *A,LL type)
{
type==1?gg=g1:gg=g2;
for(int i=1;i<limit;i++)
{
if(i<r[i])zwap(A[i],A[r[i]]);
}
for(LL lon=1;lon<limit;lon<<=1)//表示由长度为2^i合并到2^i+1
{
LL llon=lon<<1;
g[0]=1;g[1]=ksm(gg,(mod-1)/llon);
for(LL j=2;j<lon;j++)g[j]=g[j-1]*g[1]%mod;
for(LL L=0;L<limit;L+=llon)
{
LL mid=L+lon;
for(LL k=0;k<lon;k++)
{
LL x=A[k+L],y=A[k+mid]*g[k]%mod;
A[k+L]=(x+y)%mod;A[k+mid]=(x-y+mod)%mod;
}
}
}
}
int main()
{
scanf("%lld%lld",&n,&m);
for(LL i=0;i<=n;i++)scanf("%lld",&a[i]);
for(LL i=0;i<=m;i++)scanf("%lld",&b[i]);
limit=1;l=0;
while(limit<=n+m)limit<<=1,l++;
r[0]=0;for(int i=1;i<limit;i++)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
NTT(a,1);NTT(b,1);
for(LL i=0;i<limit;i++)a[i]=a[i]*b[i]%mod;
NTT(a,-1);
LL ed=n+m,fuck=ksm(limit,mod-2);
for(LL i=0;i<=ed;i++)printf("%lld ",a[i]*fuck%mod);//我之前就犯了第6点错误QAQ
printf("\n");
return 0;
}
任意模数NTT
前排说明，没有代码，没有代码！
当然，这道题目用MTT也是可以轻松过掉的
当模数为任意数字的时候，我们便不再考虑用使用题目中的模数来跑NTT了，而应该考虑一种更加神仙的东西，即：中国剩余定理合并。
即我们考虑用大质数作为模数，扩大值域，然后直接计算出其准确的值再进行模，讲真其实MTT又何尝不是呢？
考虑三个大质数 a , b , c a,b,c a,b,c，通过中国剩余定理我们知道，我们可以把值域 [ 0 , a − 1 ] , [ 0 , b − 1 ] , [ 0 , c − 1 ] [0,a-1],[0,b-1],[0,c-1] [0,a−1],[0,b−1],[0,c−1]合并成 [ 0 , a b c − 1 ] [0,abc-1] [0,abc−1]。
当然，对于 a , b , c a,b,c a,b,c的选取，要注意 a 2 , b 2 , c 2 a^2,b^2,c^2 a2,b2,c2都要在long long范围内，否则同样会爆掉，这可就功亏一篑了啊。
当然，对于 a , b , c a,b,c a,b,c的选取，有个经典组合：998244353,1004535809,469762049，而且这三个模数的原根都是 3 3 3，也方便写代码。
然后考虑一下这道题目的值域为多少：
大概是在 1 0 9 ∗ 1 0 9 ∗ 1 0 5 = 1 0 23 10^9*10^9*10^5=10^{23} 109∗109∗105=1023级别。
所以如果三个大质数的值域在 1 0 24 10^{24} 1024级别是绝对不会爆炸的。
而三个经典质数的值域可以去到：
4.7106432275 × 1 0 26 4.7106432275×10^{26} 4.7106432275×1026
完全没有任何担忧。
至于合并，其实在合并的时候也有技巧可以不去使用快速乘，又是一个黑科技呢。
首先，已知：
x
≡
x
1
m
o
d
a
x≡x_1\mod a
x≡x1moda
x
≡
x
2
m
o
d
b
x≡x_2\mod b
x≡x2modb
x
≡
x
3
m
o
d
c
x≡x_3\mod c
x≡x3modc
前两个直接合并，当然，合并方法可以是CRT合并，可以是EXGCD合并，当然，也可以用下面的奇技淫巧。
x
≡
x
1
+
k
1
a
m
o
d
b
x≡x_1+k_1a\mod b
x≡x1+k1amodb
x
2
≡
x
1
+
k
1
a
m
o
d
b
x_2≡x_1+k_1a\mod b
x2≡x1+k1amodb
所以 k 1 ≡ x 2 − x 1 a m o d b k1≡\frac{x2-x1}{a}\mod b k1≡ax2−x1modb
求逆元即可，轻松合并，时间复杂度都是一样的，别想着偷懒。
假定 x 4 ≡ x 1 + k 1 a m o d a b x_4≡x_1+k_1a\mod ab x4≡x1+k1amodab
可以发现， x 4 x_4 x4在 m o d a \mod a moda和 m o d b \mod b modb的时候都满足对应的要求。
当然，第二次就是 a b ab ab和 c c c的合并。
同理：
x
≡
x
4
+
k
2
a
b
m
o
d
c
x≡x_4+k_2ab\mod c
x≡x4+k2abmodc
x
3
≡
x
4
+
k
2
a
b
m
o
d
c
x_3≡x_4+k_2ab\mod c
x3≡x4+k2abmodc
所以 k 2 ≡ x 3 − x 4 a b m o d c k2≡\frac{x3-x4}{ab}\mod c k2≡abx3−x4modc
然后同理便可以求出 k 2 k2 k2，所以 x = x 4 + k 2 a b x=x_{4}+k_2ab x=x4+k2ab，代回去检验也会发现没有任何的问题。
而且，可以发现 x ∈ [ 0 , a b c ) , x 4 ∈ [ 0 , a b ) , k 2 ∈ [ 0 , c ) x∈[0,abc),x_{4}∈[0,ab),k_2∈[0,c) x∈[0,abc),x4∈[0,ab),k2∈[0,c)，除了 x x x都没有爆long long，可以完美的模其给定的素数 p p p完成计算，太漂亮。
看完了还不赶紧A了模板题？
🐎代码去！
我就不🐎了。
分治FFT
还是没有🐎代码。
讲真我觉得这玩意实际上就是FFT放在了分治上，甚至你用NTT放在分治上有时候效果都比这个好。
其实就是分治，然后计算左边对右边的贡献的时候用一下FFT统计贡献就行了。
时间复杂度： T ( n ) = n log n + 2 T ( n ) T(n)=n\log{n}+2T(n) T(n)=nlogn+2T(n)
这个复杂度啊，那不就是 O ( n log 2 n ) O(n\log^2{n}) O(nlog2n)
上模板，🐎！
当然，在这里顺便提一下这个模板的另外一个做法。
我的牛顿迭代做法要求复合函数，时间复杂度爆炸了啊QAQ。
多项式求逆
转载自https://www.luogu.com.cn/blog/cain12/solution-p4721
时间复杂度：
O
(
n
log
n
)
O(n\log{n})
O(nlogn)
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
#define md 998244353LL
#define inv(x) power(x,md-2)
#define max_n 1000000
using namespace std;
int n,N,i;
long long num[max_n],ans[max_n];
//const double pi=acos(-1);
//
//struct complex
//{
// double x,y;
// complex operator +(complex a)
// {
// return (complex){x+a.x,y+a.y};
// }
// complex operator -(complex a)
// {
// return (complex){x-a.x,y-a.y};
// }
// complex operator *(complex a)
// {
// return (complex){x*a.x-y*a.y,x*a.y+y*a.x};
// }
// complex operator /(double a)
// {
// return (complex){x/a,y/a};
// }
//};
long long power(long long a,long long n)
{
long long ans=1;
while(n)
{
if(n&1)
ans=ans*a%md;
a=a*a%md;
n>>=1;
}
return ans;
}
int re[max_n];
void build_re(int N)
{
int x,now=N>>1,len=0;
while(now)
len++,now>>=1;
for(int i=0;i<N;i++)
{
x=i;now=0;
for(int j=0;j<len;j++)
now=(now<<1)|(x&1),x>>=1;
re[i]=now;
}
}
//void FFT(complex a[],int len,int opt)
//{
// for(int i=0;i<len;i++)
// if(i<re[i])
// swap(a[i],a[re[i]]);
// complex wn,w,x;
// for(int i=2;i<=len;i<<=1)
// for(int j=(wn=(complex){cos(pi*2/i),opt*sin(pi*2/i)},0);j<len;j+=i)
// for(int k=(w=(complex){1,0},0);k<i>>1;k++,w=w*wn)
// x=w*a[j+k+(i>>1)],a[j+k+(i>>1)]=a[j+k]-x,a[j+k]=a[j+k]+x;
// if(opt==-1)
// for(int i=0;i<len;i++)
// a[i]=a[i]/len;
//}
void NTT(long long a[],int len,int opt)
{
for(int i=0;i<len;i++)
if(i<re[i])
swap(a[i],a[re[i]]);
long long wn,w,x;
for(int i=2;i<=len;i<<=1)
for(int j=(wn=power(3,(md-1)/i),(opt==-1?wn=power(wn,md-2):0),0);j<len;j+=i)
for(int k=(w=1,0);k<i>>1;k++,w=w*wn%md)
x=a[j+k+(i>>1)]*w%md,a[j+k+(i>>1)]=(a[j+k]-x+md)%md,a[j+k]=(a[j+k]+x)%md;
if(opt==-1)
{
long long inv=power(len,md-2);
for(int i=0;i<len;i++)
a[i]=a[i]*inv%md;
}
}
//void poly_mul(long long a[],long long b[],long long tar[],int len)
//{
// static long long A[4*max_n],B[4*max_n];
// memcpy(A,a,sizeof(long long)*len);memcpy(B,b,sizeof(long long)*len);
// build_re(len);
// NTT(A,len,1);NTT(B,len,1);
// for(int i=0;i<len;i++)
// tar[i]=A[i]*B[i]%md;
// NTT(tar,len,-1);
//}
void poly_inv(long long a[],long long tar[],int len)
{
static long long now[4*max_n],ret[4*max_n];
memset(now,0,sizeof(long long)*len);memset(ret,0,sizeof(long long)*len);
ret[0]=inv(a[0]);
for(int i=2;i<=len;i<<=1)
{
build_re(i<<1);
memcpy(now,a,sizeof(long long)*i);
NTT(ret,i<<1,1);NTT(now,i<<1,1);
for(int j=0;j<i<<1;j++)
ret[j]=ret[j]*(2LL-now[j]*ret[j]%md+md)%md;
NTT(ret,i<<1,-1);
memset(ret+i,0,sizeof(long long)*i);
}
memcpy(tar,ret,sizeof(long long)*len);
}
//void poly_dir(long long a[],long long tar[],int len)
//{
// for(int i=0;i<len-1;i++)
// tar[i]=a[i+1]*(i+1)%md;
//}
//void poly_int(long long a[],long long tar[],int len)
//{
// for(int i=len;i>=1;i--)
// tar[i]=a[i-1]*inv(i)%md;
// tar[0]=0;
//}
//void poly_ln(long long a[],long long tar[],int len)
//{
// static long long now[4*max_n],ret[4*max_n];
// memset(now,0,sizeof(long long)*len);memset(ret,0,sizeof(long long)*len);
// poly_inv(a,ret,len);
// poly_dir(a,now,len);
// poly_mul(now,ret,ret,len<<1);
// poly_int(ret,ret,len);
// memcpy(tar,ret,sizeof(long long)*len);
//}
//void poly_exp(long long a[],long long tar[],int len)
//{
// static long long now[4*max_n],ret[4*max_n],ln[4*max_n];
// memset(now,0,sizeof(long long)*len);memset(ret,0,sizeof(long long)*len);memset(ln,0,sizeof(long long)*len);
// ret[0]=1;
// for(int i=2;i<=len;i<<=1)
// {
// poly_ln(ret,now,i);
// for(int j=0;j<i;j++)
// now[j]=(a[j]-now[j]+(j==0)+md)%md;
// poly_mul(now,ret,ret,i<<1);
// }
// memcpy(tar,ret,sizeof(long long)*len);
//}
char Getchar()
{
return getchar();
static char buff[1000000],*p,*end=p;
if(p==end)
end=buff+fread(p=buff,1,1000000,stdin);
return *(p++);
}
template<typename T>void read(T &x)
{
static char rc;
x=0;rc=Getchar();
while(!isdigit(rc))
rc=Getchar();
while(isdigit(rc))
x=x*10+rc-'0',rc=Getchar();
}
int main()
{
read(n);num[0]=1;
for(i=1;i<n;i++)
read(num[i]),num[i]=-num[i];
N=1;
while(N<n)
N<<=1;
poly_inv(num,ans,N);
for(i=0;i<n;i++)
cout<<ans[i]<<' ';
return 0;
}
附原根表
转自https://www.cnblogs.com/Guess2/p/8422205.html
常用素数：
P = 1004535809 ====> pr = 3
P = 998244353 =====> pr = 3
//（g 是mod(r*2^k+1)的原根）
素数 r k g
3 1 1 2
5 1 2 2
17 1 4 3
97 3 5 5
193 3 6 5
257 1 8 3
7681 15 9 17
12289 3 12 11
40961 5 13 3
65537 1 16 3
786433 3 18 10
5767169 11 19 3
7340033 7 20 3
23068673 11 21 3
104857601 25 22 3
167772161 5 25 3
469762049 7 26 3
1004535809 479 21 3
2013265921 15 27 31
2281701377 17 27 3
3221225473 3 30 5
75161927681 35 31 3
77309411329 9 33 7
206158430209 3 36 22
2061584302081 15 37 7
2748779069441 5 39 3
6597069766657 3 41 5
39582418599937 9 42 5
79164837199873 9 43 5
263882790666241 15 44 7
1231453023109121 35 45 3
1337006139375617 19 46 3
3799912185593857 27 47 5
4222124650659841 15 48 19
7881299347898369 7 50 6
31525197391593473 7 52 3
180143985094819841 5 55 6
1945555039024054273 27 56 5
4179340454199820289 29 57 3
可以根据自己的癖好选择
参考资料
这个师兄真的写的很不错啊QAQ，直接看懂了啊（虽然FFT是在你谷学的）：https://blog.csdn.net/a_forever_dream/article/details/106520376
奆佬博客Orz：https://www.luogu.com.cn/blog/command-block/wei-yun-suan-juan-ji-yu-ji-kuo-zhan
1072
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
