标签：
FFT基本操作。。。讲解请自己看大学信号转置系列。。。
15-5-30更新：改成结构体的，跪烂王学长啊啊啊啊机智的static。。。
1 #include<iostream> 2 #include<cstdio> 3 #include<cmath> 4 #include<algorithm> 5 #include<queue> 6 #include<cstring> 7 #define PAU putchar(‘ ‘) 8 #define ENT putchar(‘\n‘) 9 #pragma comment(linker,"/STACK:10240000,10240000") 10 using namespace std; 11 const double PI=acos(-1.0); 12 const int maxn=800000+10; 13 struct FFT{ 14 struct cox{ 15 double r,i;cox(double _r = 0.0,double _i = 0.0){r=_r;i=_i;} 16 cox operator +(const cox &b){return cox(r+b.r,i+b.i);} 17 cox operator -(const cox &b){return cox(r-b.r,i-b.i);} 18 cox operator *(const cox &b){return cox(r*b.r-i*b.i,r*b.i+i*b.r);} 19 }f[maxn];int len,lenx; 20 void init(int*s,int L,int Len){ 21 len=L;lenx=Len; 22 for(int i=0;i<L;i++) f[i]=cox(s[L-1-i],0); 23 for(int i=L;i<Len;i++) f[i]=cox(0,0);return; 24 } 25 void change(){ 26 for(int i=1,j=lenx>>1;i<lenx-1;i++){ 27 if(i<j)swap(f[i],f[j]); 28 int k=lenx>>1; 29 while(j>=k) j-=k,k>>=1; 30 if(j<k) j+=k; 31 } return; 32 } 33 void cal(int tp){ 34 change(); 35 for(int h=2;h<=lenx;h<<=1){ 36 double tr=-tp*2*PI/h; 37 cox wn(cos(tr),sin(tr)); 38 for(int j=0;j<lenx;j+=h){ 39 cox w(1,0); 40 for(int k=j;k<j+(h>>1);k++){ 41 cox u=f[k],t=w*f[k+(h>>1)]; 42 f[k]=u+t;f[k+(h>>1)]=u-t;w=w*wn; 43 } 44 } 45 } if(tp==-1) for(int i=0;i<lenx;i++) f[i].r/=lenx;return; 46 } 47 }; 48 void mul(int*s1,int*s2,int L1,int L2,int&L,int*ans){ 49 L=1;while(L<L1<<1||L<L2<<1) L<<=1; 50 static FFT a,b;a.init(s1,L1,L);b.init(s2,L2,L);a.cal(1);b.cal(1); 51 for(int i=0;i<L;i++) a.f[i]=a.f[i]*b.f[i];a.cal(-1); 52 for(int i=0;i<L;i++) ans[i]=(int){a.f[i].r+0.5};return; 53 } 54 int s1[maxn>>1],s2[maxn>>1],ans[maxn],L1,L2,L; 55 void init(){ 56 char ch;int tot; 57 do ch=getchar(); while(!isdigit(ch));tot=0; 58 while(isdigit(ch)){s1[tot++]=ch-‘0‘;ch=getchar();}L1=tot; 59 do ch=getchar(); while(!isdigit(ch));tot=0; 60 while(isdigit(ch)){s2[tot++]=ch-‘0‘;ch=getchar();}L2=tot; 61 mul(s1,s2,L1,L2,L,ans); 62 return; 63 } 64 void work(){ 65 for(int i=0;i<L;i++){ 66 ans[i+1]+=ans[i]/10;ans[i]%=10; 67 } L=L1+L2-1;return; 68 } 69 void print(){ 70 while(ans[L]<=0&&L>0) L--; 71 for(int i=L;i>=0;i--) putchar(ans[i]+‘0‘); 72 return; 73 } 74 int main(){init();work();print();return 0;}
原来写的丑哭了TAT
1 #include<iostream> 2 #include<cstdio> 3 #include<cmath> 4 #include<algorithm> 5 #include<queue> 6 #include<cstring> 7 #define PAU putchar(‘ ‘) 8 #define ENT putchar(‘\n‘) 9 using namespace std; 10 const double PI=acos(-1.0); 11 struct complex{ 12 double r,i; 13 complex(double _r = 0.0,double _i = 0.0){r=_r;i=_i;} 14 complex operator +(const complex &b){return complex(r+b.r,i+b.i);} 15 complex operator -(const complex &b){return complex(r-b.r,i-b.i);} 16 complex operator *(const complex &b){return complex(r*b.r-i*b.i,r*b.i+i*b.r);} 17 }; 18 void change(complex y[],int len){ 19 int i,j,k; 20 for(i=1,j=len/2;i<len-1;i++){ 21 if(i<j)swap(y[i],y[j]); 22 k=len/2; 23 while(j>=k) j-=k,k/=2; 24 if(j<k) j+=k; 25 } return; 26 } 27 void fft(complex y[],int len,int on){ 28 change(y,len); 29 for(int h=2;h<=len;h<<=1){ 30 complex wn(cos(-on*2*PI/h),sin(-on*2*PI/h)); 31 for(int j=0;j<len;j+=h){ 32 complex w(1,0); 33 for(int k=j;k<j+h/2;k++){ 34 complex u=y[k],t=w*y[k+h/2]; 35 y[k]=u+t; 36 y[k+h/2]=u-t; 37 w=w*wn; 38 } 39 } 40 } 41 if(on==-1) for(int i=0;i<len;i++) y[i].r/=len; 42 return; 43 } 44 const int MAXN=2000+10; 45 complex x1[MAXN],x2[MAXN]; 46 char str1[MAXN>>1],str2[MAXN>>1]; 47 int sum[MAXN]; 48 inline int read(){ 49 int x=0,sig=1;char ch=getchar(); 50 while(!isdigit(ch)){if(ch==‘-‘)sig=-1;ch=getchar();} 51 while(isdigit(ch))x=10*x+ch-‘0‘,ch=getchar(); 52 return x*=sig; 53 } 54 inline void write(int x){ 55 if(x==0){putchar(‘0‘);return;}if(x<0)putchar(‘-‘),x=-x; 56 int len=0,buf[15];while(x)buf[len++]=x%10,x/=10; 57 for(int i=len-1;i>=0;i--)putchar(buf[i]+‘0‘);return; 58 } 59 int len,len1,len2; 60 void init(){ 61 scanf("%s%s",str1,str2); 62 len1=strlen(str1),len2=strlen(str2),len=1; 63 return; 64 } 65 void work(){ 66 while(len<len1<<1||len<len2<<1) len<<=1; 67 for(int i=0;i<len1;i++) x1[i]=complex(str1[len1-1-i]-‘0‘,0); 68 for(int i=len1;i<len;i++) x1[i]=complex(0,0); 69 for(int i=0;i<len2;i++) x2[i]=complex(str2[len2-1-i]-‘0‘,0); 70 for(int i=len2;i<len;i++) x2[i]=complex(0,0); 71 fft(x1,len,1);fft(x2,len,1); 72 for(int i=0;i<len;i++) x1[i]=x1[i]*x2[i]; 73 fft(x1,len,-1); 74 for(int i=0;i<len;i++) sum[i]=(int)(x1[i].r+0.5); 75 for(int i=0;i<len;i++){ 76 sum[i+1]+=sum[i]/10; 77 sum[i]%=10; 78 } len=len1+len2-1;return; 79 } 80 void print(){ 81 while(sum[len]<=0&&len>0) len--; 82 for(int i=len;i>=0;i--) putchar(sum[i]+‘0‘); 83 return; 84 } 85 int main(){init();work();print();return 0;}
标签：
原文地址：http://www.cnblogs.com/chxer/p/4507345.html
