FFT快速傅立叶

时间:2023-03-09 16:36:15
FFT快速傅立叶

Description

给出两个n位10进制整数x和y,你需要计算x*y。

Input

第一行一个正整数n。第二行描述一个位数为n的正整数x。第三行描述一个位数为n的正整数y。

Output

输出一行,即x*y的结果。

Sample Input

1
3
4

Sample Output

12

数据范围:
n<=60000

HINT

Source

FFT裸题

具体见算导

code:

 #include<cstdio>
#include<iostream>
#include<cmath>
#include<cstring>
#include<algorithm>
#define maxn 131072
#define pi 3.14159265358979323846
using namespace std;
char ch;
int m,n,ans[maxn];
bool ok;
struct comp{
double rea,ima;
void clear(){rea=ima=;}
comp operator +(const comp &x){return (comp){rea+x.rea,ima+x.ima};}
comp operator -(const comp &x){return (comp){rea-x.rea,ima-x.ima};}
comp operator *(const comp &x){return (comp){rea*x.rea-ima*x.ima,rea*x.ima+ima*x.rea};}
}a[maxn],b[maxn],c[maxn],tmp[maxn],w,wn;
void read(int &x){
for (ok=,ch=getchar();!isdigit(ch);ch=getchar()) if (ch=='-') ok=;
for (x=;isdigit(ch);x=x*+ch-'',ch=getchar());
if (ok) x=-x;
}
void read(comp *a){
for (int i=m-;i>=;i--){
for (ch=getchar();!isdigit(ch);ch=getchar());
a[i].rea=ch-'';
}
}
void write(int *a){
int i;
for (i=n-;i>=&&!a[i];i--);
if (i<){puts("");return;}
for (;i>=;i--) printf("%d",(int)a[i]);
puts("");
}
void fft(comp *a,int st,int siz,int step,int op){
if (siz==) return;
fft(a,st,siz>>,step<<,op),fft(a,st+step,siz>>,step<<,op);
int x=st,x1=st,x2=st+step;
w=(comp){,},wn=(comp){cos(op**pi/siz),sin(op**pi/siz)};
for (int i=;i<(siz>>);i++,x+=step,x1+=(step<<),x2+=(step<<),w=w*wn)
tmp[x]=a[x1]+(w*a[x2]),tmp[x+(siz>>)*step]=a[x1]-(w*a[x2]);
for (int i=st;siz;i+=step,siz--) a[i]=tmp[i];
}
int main(){
read(m),n=;
while (n<(m<<)) n<<=;
read(a),read(b);
fft(a,,n,,),fft(b,,n,,);
for (int i=;i<n;i++) c[i]=a[i]*b[i];
fft(c,,n,,-);
for (int i=;i<n;i++) ans[i]=(int)round(c[i].rea/(1.0*n));
for (int i=;i<n;i++) ans[i+]+=ans[i]/,ans[i]%=;
write(ans);
system("pause");
return ;
}