第一次学FFT,先膜拜一下法法塔大神ORZ
关于FFT的话,有一篇博文特别赞http://z55250825.blog.163.com/blog/static/150230809201431274653644/
他后面还有关于高精度和jsoi2014 力的题解写的特别好
其次算导讲的真的不错
不过这篇博文讲得更算导差不多了ORZ
直接上代码吧
尼玛重载运算符老写错QAQ
好吧突然发现以前有一点错误,然后插了别人的代码来check,后来自己的就没了= = sorry
CODE:
#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cmath>
#include<cstring>
#define pi acos(-1.0)
using namespace std;
#define maxn 400010
struct com{
double r,i;
com(double real=0.0,double imag=0.0){
r=real,i=imag;
}
com operator + (const com x){
return com(x.r+r,i+x.i);
}
com operator - (const com x){
return com(r-x.r,i-x.i);
}
com operator * (const com x){
return com(x.r*r-x.i*i,r*x.i+i*x.r);
}
}x1[maxn*2],x2[maxn*2];
bool bo[maxn*2];
int bitrev(com *a,int l){
memset(bo,0,sizeof(bo));//O(logn)的在后面
for(int i=1;i<l-1;i++){
int x=i,y=0;
int m=(int)log2(l)+0.1;
if(bo[x])continue;
while(m--){
y<<=1;
y|=(x&1);
x>>=1;
}
bo[i]=bo[y]=1;
swap(a[i],a[y]);
}
return 0;
}
void fft(com *y,int l,double on) // FFT O(nlogn)
// 其中on==1时为DFT,on==-1为IDFT
{
register int h,i,j,k;
com u,t;
bitrev(y,l); // 调用反转置换
for(h=2;h<=l;h<<=1) // 控制层数
{
// 初始化单位复根
com wn(cos(on*2*pi/h),sin(on*2*pi/h));
for(j=0;j<l;j+=h) // 控制起始下标
{
com w(1,0); // 初始化螺旋因子
for(k=j;k<j+h/2;k++) // 配对
{
u=y[k];
t=w*y[k+h/2];
y[k]=u+t;
y[k+h/2]=u-t;
w=w*wn; // 更新螺旋因子
} // 据说上面的操作叫蝴蝶操作…
}
}
if(on==-1) for(i=0;i<l;i++) y[i].r/=l; // IDFT
}
char a[maxn],b[maxn];
int sum[maxn*2];
int main(){
int l1,l2,l;
scanf("%s%s",a,b);
l1=strlen(a);l2=strlen(b);
l=1;
while (l<l1*2||l<l2*2) {
l<<=1;
}
for (int j=0;j<l1;j++) {
x1[j].r=a[l1-j-1]-'0';
x1[j].i=0.0;
}
for (int i=l1;i<l;i++) x1[i]=com(0.0,0.0);
for (int j=0;j<l2;j++) {
x2[j].r=b[l2-j-1]-'0';
x2[j].i=0.0;
}
for (int i=l2;i<l;i++) x2[i]=com(0.0,0.0);
fft(x1,l,1);
fft(x2,l,1);
for (int i=0;i<l;i++) {
x1[i]=x1[i]*x2[i];
}
fft(x1,l,-1);
for (int i=0;i<l;i++) sum[i]=x1[i].r+0.5;
for (int i=0;i<l;i++) {
sum[i+1]+=sum[i]/10;
sum[i]%=10;
}
l=l1+l2-1;
while(sum[l]<=0 && l>0) l--; // 检索最高位
for(int i=l;i>=0;i--) putchar(sum[i]+'0'); // 倒序输出
return 0;
}