(http://www.elijahqi.win/2018/02/04/uoj34-%E5%A4%9A%E9%A1%B9%E5%BC%8F%E4%B9%98%E6%B3%95-fft%E6%9D%BF%E5%AD%90%E7%9E%8E%E6%90%9E%E5%AD%A6%E4%B9%A0xx%E8%A2%AB%E8%99%90%E8%AE%B0/)
这是一道模板题。
给你两个多项式,请输出乘起来后的多项式。
输入格式
第一行两个整数
n
n
和
m
m
,分别表示两个多项式的次数。
第二行
n+1
n+1
个整数,表示第一个多项式的
0
0
到
n
n
次项系数。
第三行
m+1
m+1
个整数,表示第二个多项式的
0
0
到
m
m
次项系数。
输出格式
一行
n+m+1
n+m+1
个整数,表示乘起来后的多项式的
0
0
到
n+m
n+m
次项系数。
样例一
input
1 2
1 2
1 2 1
output
1 4 5 2
explanation
(1+2x)⋅(1+2x+
x
2
)=1+4x+5
x
2
+2
x
3
(1+2x)⋅(1+2x+x2)=1+4x+5x2+2x3
。
限制与约定
0≤n,m≤
10
5
0≤n,m≤105
,保证输入中的系数大于等于
0
0
且小于等于
9
9
。
时间限制:
1s
1s
蒟蒻我暂时只会搞多项式乘法?
n-1次多项式A(x)与m-1次多项式B(x)相乘获得l=n+m-1次多项式c(x) 直接朴素做是n*m的复杂度 对于这题 显然不可做 考虑新学一种算法fft
C(x)=A(x)*B(x)
Ci=∑n−1p=0∑m−1q=0[p+q==i]Ap∗Bq
C
i
=
∑
p
=
0
n
−
1
∑
q
=
0
m
−
1
[
p
+
q
==
i
]
A
p
∗
B
q
Ci=∑n−1p=0∑m−1q=0[p+q≡i(mod l)]Ap∗Bq
C
i
=
∑
p
=
0
n
−
1
∑
q
=
0
m
−
1
[
p
+
q
≡
i
(
m
o
d
l
)
]
A
p
∗
B
q
已知一个式子 对于任意v均满足
如果需要知道这个式子首先需要知道复数的应用 例如欧拉公式eiπ=1
e
i
π
=
1
ωl
ω
l
表示l次单位复数根
[v≡0(mod l)]=1l∗∑l−1k=0ωv∗kl
[
v
≡
0
(
m
o
d
l
)
]
=
1
l
∗
∑
k
=
0
l
−
1
ω
l
v
∗
k
对于v≡0(mod l)
v
≡
0
(
m
o
d
l
)
显然可以知道 对于这种情况 是成立的
对于v≢0(mod l)
v
≢
0
(
m
o
d
l
)
的情况 可以知道 对于原式是
1l∗∑l−1i=0ωi∗vl
1
l
∗
∑
i
=
0
l
−
1
ω
l
i
∗
v
对于这个∑
∑
我可以考虑等比数列搞出 化简之后就是1−ωl∗vl1−ωvl=0
1
−
ω
l
l
∗
v
1
−
ω
l
v
=
0
所以可以证明这个式子成立
带入原始式子Ci=∑n−1p=0∑m−1q=0[p+q−i≡0(mod l)]Ap∗Bq
C
i
=
∑
p
=
0
n
−
1
∑
q
=
0
m
−
1
[
p
+
q
−
i
≡
0
(
m
o
d
l
)
]
A
p
∗
B
q
Ci=∑n−1p=0∑m−1q=0(1l∗∑l−1k=0ω(p+q−i)∗kl)Ap∗Bq
C
i
=
∑
p
=
0
n
−
1
∑
q
=
0
m
−
1
(
1
l
∗
∑
k
=
0
l
−
1
ω
l
(
p
+
q
−
i
)
∗
k
)
A
p
∗
B
q
交换枚举对象
Ci=1l∗∑l−1k=0ω−i∗kl(∑n−1p=0ωp∗klAp)∗(∑m−1q=0ωq∗k∗Bq)
C
i
=
1
l
∗
∑
k
=
0
l
−
1
ω
l
−
i
∗
k
(
∑
p
=
0
n
−
1
ω
l
p
∗
k
A
p
)
∗
(
∑
q
=
0
m
−
1
ω
q
∗
k
∗
B
q
)
根据多项式定义后面的两部分可以分别表示为
∑n−1p=0ωp∗klAp=A[ωkl]
∑
p
=
0
n
−
1
ω
l
p
∗
k
A
p
=
A
[
ω
l
k
]
B同理
所以ci=1l∗∑l−1k=0ω−i∗kl∗A[ωkl]∗B[ωkl]
c
i
=
1
l
∗
∑
k
=
0
l
−
1
ω
l
−
i
∗
k
∗
A
[
ω
l
k
]
∗
B
[
ω
l
k
]
设dk=A[ωkl]∗B[ωkl]
d
k
=
A
[
ω
l
k
]
∗
B
[
ω
l
k
]
显然如果我有办法快速求出A的所有点值那我可以采取同样的方法快速求出每个dk的值 相当于只是 带入的值变化了 如何快速球值这就需要我们选取的这个值即可 注意我们选择2^n次复数根J即将n,m同时扩大到和他们最近的一个2^n次幂的地方
如何快速求出A[x]即 ∑n−1i=0xi∗ai
∑
i
=
0
n
−
1
x
i
∗
a
i
显然我可以针对下标划分成奇偶两种形式 设A0表示偶数次的项 A1表示奇数次的项 A0(x)=∑n−1i=0,i+=2xi∗ai
A
0
(
x
)
=
∑
i
=
0
,
i
+
=
2
n
−
1
x
i
∗
a
i
A0(x)=∑⌊n−12⌋i=0x2i∗a2i
A
0
(
x
)
=
∑
i
=
0
⌊
n
−
1
2
⌋
x
2
i
∗
a
2
i
设di=a2i
d
i
=
a
2
i
那么可以表示成设D(x)=∑⌊n−12⌋i=0xi∗di
D
(
x
)
=
∑
i
=
0
⌊
n
−
1
2
⌋
x
i
∗
d
i
那么显然可以知道A0(x)=D(x2)
A
0
(
x
)
=
D
(
x
2
)
同理我们去搞一下奇数次的项
A1(x)=∑n−1i=1,i+=2xi∗ai
A
1
(
x
)
=
∑
i
=
1
,
i
+
=
2
n
−
1
x
i
∗
a
i
A1(x)=x∗∑⌊n−12⌋i=0x2i+1∗a2i+1
A
1
(
x
)
=
x
∗
∑
i
=
0
⌊
n
−
1
2
⌋
x
2
i
+
1
∗
a
2
i
+
1
设ei=a2i+1
e
i
=
a
2
i
+
1
E(x)=∑⌊n−12⌋i=0xi∗e2i
E
(
x
)
=
∑
i
=
0
⌊
n
−
1
2
⌋
x
i
∗
e
2
i
A1(x)=x∗E(x2)
A
1
(
x
)
=
x
∗
E
(
x
2
)
所以A(x)=A1(x)+A0(x)
A
(
x
)
=
A
1
(
x
)
+
A
0
(
x
)
当k<l2
l
2
A(ωkl)=A0(ωkl)+A1(ωkl)
A
(
ω
l
k
)
=
A
0
(
ω
l
k
)
+
A
1
(
ω
l
k
)
A(ωkl)=D(ω2kl)+ωkl∗E(ω2kl)
A
(
ω
l
k
)
=
D
(
ω
l
2
k
)
+
ω
l
k
∗
E
(
ω
l
2
k
)
A(ωkl)=D(ωk2l)+ωkl∗E(ωk2l)
A
(
ω
l
k
)
=
D
(
ω
2
l
k
)
+
ω
l
k
∗
E
(
ω
2
l
k
)
同理获得当k>=l2
k
>=
l
2
此时设定k=k−l2
k
=
k
−
l
2
A(ωk+l2l)=D(ω2k+ll)+ωk+l2l∗E(ω2k+ll)
A
(
ω
l
k
+
l
2
)
=
D
(
ω
l
2
k
+
l
)
+
ω
l
k
+
l
2
∗
E
(
ω
l
2
k
+
l
)
此时和上面相同只不过ωl2l
ω
l
l
2
是-1所以要注意负号
递归版
#include<cmath>
#include<cstdio>
#include<complex>
#define C complex<double>
#define N 300000
#define pi acos(-1)
using namespace std;
inline char gc(){
static char now[1<<16],*S,*T;
if (T==S){T=(S=now)+fread(now,1,1<<16,stdin);if (T==S) return EOF;}
return *S++;
}
inline int read(){
int x=0,f=1;char ch=gc();
while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=gc();}
while(ch<='9'&&ch>='0') x=x*10+ch-'0',ch=gc();
return x*f;
}
C a[N],b[N];
int n,m;
inline void fft(C *x,int n,int f){
if (n==1) return;C l[n>>1],r[n>>1];
for (int i=0;i<n;i+=2) l[i>>1]=x[i],r[i>>1]=x[i+1];
fft(l,n>>1,f);fft(r,n>>1,f);int nn=n>>1;
C wn(cos(2*pi/n),f*sin(2*pi/n)),w(1,0),t;
for (int i=0;i<nn;++i,w*=wn)
t=w*r[i],x[i]=l[i]+t,x[i+nn]=l[i]-t;
}
int main(){
freopen("uoj34.in","r",stdin);
n=read();m=read();
for (int i=0;i<=n;++i) a[i]=read();
for (int i=0;i<=m;++i) b[i]=read();
m=n+m;for (n=1;n<=m;n<<=1) ;
fft(a,n,1);fft(b,n,1);
for (int i=0;i<=n;++i) a[i]*=b[i];
fft(a,n,-1);
for (int i=0;i<=m;++i) printf("%d,int(a[i].real()/n+0.5));
return 0;
}
迭代版
预处理出我的变换顺序 这个R是怎么得到的 有点dp思想 就是1~i-1的已经搞好了那么我现在相当于第i位不好 所以把他取出来 然后填到第一位上即可
#include<cmath>
#include<cstdio>
#include<complex>
#include<algorithm>
#define N 300000
#define pi acos(-1)
#define C complex<double>
using namespace std;
inline char gc(){
static char now[1<<16],*S,*T;
if (T==S){T=(S=now)+fread(now,1,1<<16,stdin);if (T==S) return EOF;}
return *S++;
}
inline int read(){
int x=0,f=1;char ch=gc();
while(ch<'0'||ch>'9') {if (ch=='-') f=-1;ch=gc();}
while(ch<='9'&&ch>='0') x=x*10+ch-'0',ch=gc();
return x*f;
}
C a[N],b[N];
int n,m,R[N];
inline void fft(C *x,int f){
for (int i=0;i<n;++i) if (i<R[i]) swap(x[i],x[R[i]]);
for (int i=1;i<n;i<<=1){
C wn(cos(pi/i),f*sin(pi/i));
for (int j=0;j<n;j+=i<<1){
C w(1,0),t1,t2;
for (int k=0;k<i;++k,w*=wn)
t1=x[j+k],t2=x[j+k+i]*w,x[k+j]=t1+t2,x[k+j+i]=t1-t2;
}
}
}
int main(){
freopen("uoj34.in","r",stdin);
n=read();m=read();
for (int i=0;i<=n;++i) a[i]=read();
for (int i=0;i<=m;++i) b[i]=read();
m=n+m;for (n=1;n<=m;n<<=1);int l=log2(n);
for (int i=1;i<n;++i) R[i]=R[i>>1]>>1|(i&1)<<l-1;
//for (int i=0;i<8;++i) printf("%d,R[i]);puts("");
fft(a,1);fft(b,1);for (int i=0;i<n;++i) a[i]*=b[i];
fft(a,-1);for (int i=0;i<=m;++i) printf("%d,(int)(a[i].real()/n+0.5));
return 0;
}