FFT
- 用处：
快速计算多项式的卷积。
- 做法：
我们知道多项式有两种表示：系数表示与点值表示。系数表示就是最常见的表示法，就是拿一个向量来表示多项式的各项系数。而点值表示则是拿几个点
（
x
i
，
y
i
）
（x_i，y_i）
（xi，yi）来表示。把系数表示的多项式转成点值的过程叫做求值，反之叫插值。下面比较一下这两种方法。
对于多项式加法，点值与系数表示都是
O
（
n
+
m
）
O（n+m）
O（n+m）（
n
n
n，
m
m
m分别是原多项式的次数），但是点值的话必须选的
x
x
x都一样。
对多项式乘法，系数是
O
（
n
m
）
O（nm）
O（nm）的，而点值还是
O
（
n
+
m
）
O（n+m）
O（n+m）的（只需要把对应点乘上去）。由此可见，如果在算乘法时，表示法是点值的，那么就可以很快算出。而如果用最朴素的一个点一个点算，还是
O
（
n
m
）
O（nm）
O（nm）的，而
F
F
T
FFT
FFT则是快速的将多项式转成点值，然后再快速将点值转成系数的方法，其时间复杂度为
O
（
n
l
o
g
n
）
O（nlogn）
O（nlogn）的。其中，将系数换成点值的过程叫做
D
F
T
DFT
DFT，反之叫做
I
D
F
T
IDFT
IDFT。
- DFT的实现：
首先，设
w
n
k
=
e
k
2
π
i
n
w_n^k=e^{k \frac{2 \pi i }{n}}
wnk=ekn2πi，即
n
n
n次方程的
k
k
k次单位根。我们有如下几个极其重要的性质：
（1）相消引理
对于任何整数
n
≥
0
,
k
≥
0
,
d
>
0
n\geq 0,k\geq 0,d>0
n≥0,k≥0,d>0，有
w
d
n
d
k
=
w
n
k
w_{dn}^{dk}=w_n^k
wdndk=wnk
（2）折半引理
如果
n
>
0
n>0
n>0且为偶数，则
(
w
n
k
+
n
2
)
2
=
(
w
n
k
)
2
(w_n{k+\frac{n}{2}})^2 = (w_n^k)^2
(wnk+2n)2=(wnk)2
（2）求和引理
任意整数
n
≥
1
n\geq 1
n≥1和不能被
n
n
n整除的非零整数
k
k
k，有
∑
i
=
0
n
−
1
(
w
n
k
)
i
=
0
\sum _{i=0}^{n-1}(w_n^k)^i=0
∑i=0n−1(wnk)i=0
那么有了这三条引理，就可以进行
F
F
T
FFT
FFT了。
F
F
T
FFT
FFT的方法是：将原来的数列奇偶分类得到两个新的
A
0
A0
A0与
A
1
A1
A1。其中
A
0
A0
A0放偶数位置，A1放奇数位置的，这样子，我们可以得到一个式子：
A
(
w
n
k
)
=
A
0
(
(
w
n
k
)
2
)
+
(
w
n
k
)
A
1
(
(
w
n
k
)
2
)
A(w_n^k)=A_0((w_n^k)^2)+(w_n^k)A_1((w_n^k)^2)
A(wnk)=A0((wnk)2)+(wnk)A1((wnk)2)，这样子可以写出最基本的
F
F
T
FFT
FFT。
就是我们对于一个
n
=
2
q
n=2^q
n=2q，对于一个
i
<
n
2
i<\frac{n}{2}
i<2n,有
A
(
w
n
i
)
=
A
0
(
(
w
n
i
)
2
)
+
w
n
i
A
1
(
(
w
n
i
)
2
)
A(w_n^i)=A_0((w_n^i)^2)+w_n^iA_1((w_n^i)^2)
A(wni)=A0((wni)2)+wniA1((wni)2)
A
(
w
n
i
+
n
2
)
=
A
0
(
(
w
n
i
+
n
2
)
2
)
−
w
n
i
A
1
(
(
w
n
i
+
n
2
)
2
)
A(w_n^{i+\frac{n}{2}})=A_0((w_n^{i+\frac{n}{2}})^2)-w_n^iA_1((w_n^{i+\frac{n}{2}})^2)
A(wni+2n)=A0((wni+2n)2)−wniA1((wni+2n)2)
所以只用算
n
2
\frac{n}{2}
2n次就可以了。
我们还有一种叫做蝴蝶操作的方法降低复杂度。这个是迭代版FFT用到的，具体就是只用一个数组计算。具体的，我也不是很明白还是看代码。
代码
（递归版，常数大，空间大，好记忆）
void FFT(com *A,int len,int fl)
{
if (len==1) return;
com A0[len],A1[len];
for (int i=0;i<(len/2);i++)
{
A0[i] = A[2*i];
A1[i] = A[2*i+1];
}
FFT(A0,len/2,fl);FFT(A1,len/2,fl);
com wn(cos(2*pi/len*fl),sin(2*pi/len*fl)) , w(1,0);
for (int i=0;i<(len/2);i++)
{
A[i] = A0[i] + w * A1[i];
A[i+(len/2)] = A0[i] - w * A1[i];
w = w * wn;
}
}
（非递归版，常数，空间都小，难记忆）
void FFT2(com c[],int len,int fl)
{
int i,j,k,clen;
for (i=(len>>1),j=1;j<len;j++)
{
if (i<j) swap(c[i],c[j]);
for (k=(len>>1);i&k;k>>=1) i^=k;
i^=k;
}
for (clen=2;clen<=len;clen<<=1)
{
com wn=com(cos(2.0*pi/clen*fl),sin(2.0*pi/clen*fl));
for (j=0;j<len;j+=clen)
{
com w=com(1,0);
for (k=j;k<j+(clen>>1);k++)
{
com tmp1=c[k],tmp2=c[k+(clen>>1)]*w;
c[k]=tmp1+tmp2;c[k+(clen>>1)]=tmp1-tmp2;
w=w*wn;
}
}
}
if (fl==-1)
{
for (i=0;i<len;i++) c[i]=c[i]/len;
}
}
- IDFT
之前不是传进了一个
f
l
fl
fl吗？当
f
l
fl
fl是
1
1
1的时候，就是
F
F
T
FFT
FFT，是
−
1
-1
−1的时候，就是
I
D
F
T
IDFT
IDFT。至于怎么证明就是这样，你相信就是了,可以参见其他dalao的博客。最后算出来的结果还要除以一个
l
e
n
len
len。（这个在非递归版中已经除了，但在递归版中写不了，所以要在外面那个地方再除一次）。
NTT
N
T
T
NTT
NTT与
F
F
T
FFT
FFT的方法类似，但是我们这次要选择原根。那么原根满足之前的那三条引理吗？可以证明是满足的，所以可以直接搞。
这里补一下代码：
const ll p = 998244353;
const ll g = 3;
const ll invg = 332748118;
void NTT(ll *a,int len,int fl)
{
int i,j,k;
for (i=(len>>1),j=1 ; j<len ; j++)
{
if (i < j) swap(a[i],a[j]);
for (k = (len>>1) ;i & k; k>>=1) i ^= k;
i ^= k;
}
for (int clen=2;clen<=len;clen*=2)
{
ll wn = qpow(( fl == 1 ) ? g : invg , (mod-1) / clen );
for (int j=0;j<len;j+=clen)
{
ll w = 1;
for (int k=j;k<j+(clen>>1);k++)
{
ll x = a[k] , y = a[k+(clen>>1)] * w %mod;
a[k] = (x + y)%mod; a[k+(clen>>1)] = (x - y+mod)%mod;
w = (w * wn)%mod;
}
}
}
if (fl == -1)
{
ll inv = qpow(len,mod-2);
for (int i=0;i<len;i++) a[i] = a[i]*inv%mod;
}
}
上面的膜数以及其原根有许多取法，这里给的是最常见的
998244353
998244353
998244353与
3
3
3。
2018.3.10补：
有时候我们会看见这种题目：
C
=
A
×
B
C = A \times B
C=A×B
求
C
（
1
）
C（1）
C（1），那么就是求
C
C
C的系数之和。这种问题，我们当然可以用
F
F
T
和
N
T
T
FFT和NTT
FFT和NTT来做，但是我们有复杂度更低的做法，就是考虑
A
i
A_i
Ai会与
B
B
B的那些相乘，那么明显是
B
B
B的全部项，所以答案就是
A
A
A的系数和乘上
B
B
B的系数和。
FWT
详细的在这里：FWT。
我们之前的
F
F
T
FFT
FFT解决的多项式卷积可以写成这样的形式
C
[
k
]
=
∑
i
=
0
k
A
[
i
]
B
[
k
−
i
]
C[k] = \sum _{i=0} ^{k}A[i] B[k-i]
C[k]=i=0∑kA[i]B[k−i]
从下标的角度来看，就是我们要卷积的下标之和为新的下标。那么，这种下标就是加法意义的下标。那么
F
W
T
FWT
FWT解决便是集合意义下的下标，就是
C
[
k
]
=
∑
i
⊕
j
=
k
k
A
[
i
]
B
[
j
]
C[k] = \sum _{i⊕j=k} ^{k}A[i] B[j]
C[k]=i⊕j=k∑kA[i]B[j]
这里
⊕
⊕
⊕符号不是指异或，而是一种二进制运算。
那么怎么做呢？我先说一个简单的，与运算
a
n
d
and
and吧。
首先我们要明白，这个东西
F
W
T
FWT
FWT其实和
F
F
T
FFT
FFT差不多，都是把原式变成点值之类的东西，再把点值乘起来，然后再还原的大体思路。那么显而易见的，这种变化要满足
A
′
×
B
′
=
C
′
A'\times B'=C'
A′×B′=C′，其中
A
’
,
B
’
,
C
’
A’,B’,C’
A’,B’,C’ 是变化后的式子。显而易见要满足这一点，不然你把变化后的两个式子乘起来都不等于原式了。
那么我们就设这种变化为
F
W
T
FWT
FWT, 其逆变换为
U
F
W
T
UFWT
UFWT。那么我就说与运算怎么搞吧。
也不知道是哪个神人想出的变换，使得
A
′
[
i
]
=
∑
i
∧
j
=
i
A
[
j
]
(
∧
这
个
符
号
是
与
，
因
为
打
不
出
a
n
d
)
A'[i] = \sum _{i ∧ j=i}A[j] (∧这个符号是与，因为打不出and)
A′[i]=i∧j=i∑A[j](∧这个符号是与，因为打不出and)，莫名其妙吧，但是仔细看一下，好像是满足
A
′
×
B
′
=
C
′
A'\times B'=C'
A′×B′=C′这个东西的。那么我们考虑怎么求这个
A
′
A'
A′。
类似于
F
F
T
FFT
FFT，我们考虑分治做，同样我们设一个
A
0
′
,
A
1
′
A'_0,A'_1
A0′,A1′，但是与
F
F
T
FFT
FFT中不同，我们设
A
0
′
A'_0
A0′是
A
′
A'
A′的前一半。设
A
′
A'
A′的长度为
2
k
2^k
2k，那么这一半的下标二进制下的最前面一位第
k
k
k位是为
0
0
0的，然后
A
1
′
A'_1
A1′就是后面一半，那么最前面一位第
k
k
k位是为
1
1
1的。这样有啥好处？就是我们就可以知道，如果算出了
A
0
′
,
A
1
′
A'_0,A'_1
A0′,A1′，那么我们可以知道
A
′
A'
A′的后一部分必须还得是后一部分，还是
A
1
′
A'_1
A1′，而前一部分因为是与运算，所以
(
0
∧
1
)
(
0
∧
0
)
(0∧1)(0∧0)
(0∧1)(0∧0)的结果都是
0
0
0，所以前一部分是
A
0
′
+
A
1
′
A'_0+A'_1
A0′+A1′。
那么这样子就构造完毕了，至于逆运算，也很简单。你都知道了
A
′
0
,
A
′
1
A'0,A'1
A′0,A′1与
A
0
,
A
1
A0,A1
A0,A1的关系了，那么之前是
A
0
=
A
0
′
+
A
1
′
,
A
1
=
A
1
′
A_0 = A'_0+A'_1,A_1 = A'_1
A0=A0′+A1′,A1=A1′，那么我们可以把式子反过来，变成
A
0
′
=
A
0
−
A
1
,
A
1
′
=
A
1
A'_0=A_0-A_1,A'_1=A_1
A0′=A0−A1,A1′=A1，就是简单的反过来。
好这样子的话，就把与运算解决了，那么其他的常见运算或运算与异或运算读者可以自行去做，下面贴出结论：
或运算的
A
0
=
A
0
′
,
A
1
=
A
0
′
+
A
1
′
A_0 = A'_0 , A_1 = A'_0 + A'_1
A0=A0′,A1=A0′+A1′
A
0
′
=
A
0
,
A
1
′
=
A
0
+
A
1
A'_0 = A_0 , A'_1 = A_0 + A_1
A0′=A0,A1′=A0+A1
##鬼畜的异或运算的
A
0
=
A
0
′
+
A
1
′
,
A
1
=
A
0
′
−
A
1
′
A_0 = A'_0 + A'_1, A_1 = A'_0 - A'_1
A0=A0′+A1′,A1=A0′−A1′
A
0
′
=
A
0
+
A
1
2
,
A
1
′
=
A
0
−
A
1
2
A'_0 = \frac{A_0 + A_1 }{ 2}, A'_1 = \frac{A_0 - A_1}{2}
A0′=2A0+A1,A1′=2A0−A1
下面是代码：
void FWT(ll *c,ll len)
{
ll i,j,k,clen;
for (i=(len>>1),j=1 ; j<len ; j++)
{
if (i < j) swap(c[i],c[j]);
for (k = (len>>1) ;i & k; k>>=1) i ^= k;
i ^= k;
}
for (clen = 2;clen<= len;clen<<=1)
{
for (j = 0;j < len;j += clen)
{
for ( k = j;k < j + (clen>>1) ;k++)
{
ll x = c[k] , y = c[k + (clen>>1) ] ;
// xor c[k] = x + y , c[k + (clen>>1)] = x-y;
// and c[k] = x + y , c[k + (clen>>1)] = c[k + (clen>>1)]
// or c[k] = c[k] c[k + (clen>>1)] = x+y
c[k] = x + y ; c[k + (clen>>1) ] = x - y ;
}
}
}
}
void UFWT(ll *c,ll len)
{
ll i,j,k,clen;
for (i=(len>>1),j=1 ; j<len ; j++)
{
if (i < j) swap(c[i],c[j]);
for (k = (len>>1) ;i & k; k>>=1) i ^= k;
i ^= k;
}
for (clen = 2;clen<= len;clen<<=1)
{
for (j = 0;j < len;j += clen)
{
for ( k = j;k < j + (clen>>1) ;k++)
{
ll x = c[k] , y = c[k + (clen>>1) ] ;
// xor c[k] = x + y / 2, c[k + (clen>>1)] = x-y / 2;
// and c[k] = x - y , c[k + (clen>>1)] = c[k + (clen>>1)]
// or c[k] = c[k] c[k + (clen>>1)] = y-x
c[k] = (x + y)*rev ; c[k + (clen>>1) ] = (x - y)*rev ;//rev 是2的逆元
}
}
}
}
细心的读者肯定发现了，这不就是刚才的
F
F
T
FFT
FFT的版改了一下运算部分吗？没错，就是这样，所以对于
F
W
T
FWT
FWT，只用记住那几个运算的变换，然后在套上
F
F
T
FFT
FFT的模板就可以了。
这里提几个细节，就是传
l
e
n
len
len的时候，一定要传
2
k
2^k
2k这种
l
e
n
len
len，而不是
2
k
−
1
2^k-1
2k−1这种看似正确的
l
e
n
len
len，因为你的
l
e
n
len
len的意义是有多长，所以下标为
0
0
0的那个也要算。
然后，前面那个二进制分组可以去掉，因为
F
W
T
FWT
FWT相当于是每一位独立算的，可以不用二进制分组。
还有代码中那个
r
e
v
rev
rev，它的意思是
2
2
2在你要模的数的意义下的逆元。
K进制下的FWT
void FWT(ll *c,ll len) {
ll i , j , clen;
for (clen = 1;clen < len;clen *= k) {
for (i = 0;i < len;i += clen * k) {
for (j = i;j < i + clen;i++) {
ll a[j] , a[j + len] , a[j + len * 2] , a[j + len * 3]...
a[j] = w(0 * 0) a[j] + w(0 * 1) a[j + len] + w(0 * 2) a[j + len * 2]...;
a[j + len] = w(1 * 0) a[j] + w(1 * 1) a[j + len] + w(1 * 2) a[j + len * 2]...;
}
}
}
}
扫一扫
541
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
