快速傅里叶变换学习笔记（总结向）

初次接触请点击：一小时学会快速傅里叶变换（Fast Fourier Transform）

关键字：多项式系数表达法与点值表达法（有关插值可行性的证明），复数，单位根，n次单位向量，折半引理，消去引理，

　　蝴蝶操作（为减小常数、空间）

复数类型可以自己实现，也可以使用C++自带STL库中的complex ，详细见：C++STL complex吃书使用指南

小注意点：

　　1、最终的n为2的整数次幂（否则单位根数组下标会有分数的类似情况，很不好处理。数组下标不能为实型）。

　　2、因为运算多次设计复数，要注意各个运算符链接的数的类型。

　　3、注意精度，对一个double型数a四舍五入取整，可直接（int）（a+0.5），这是精度最高的一种方法。

　　4、最后答案的形式多为整数，注意这时整数的规模范围为做乘法前的多项式中单个的系数规模的平方再乘项数（总之多半情况会很大，int型容易爆）

题目：P1919 【模板】A*B Problem升级版（FFT快速傅里叶）

  1 #include<complex>
  2 #include<iostream>
  3 #include<cstdio>
  4 #include<cmath>
  5 #include<cstring>
  6 //  *标记为FFT重点代码
  7 using namespace std;
  8 
  9 const int N=524288;
 10 
 11 const double PI=acos(-1);
 12 
 13 
 14 double a[N],b[N];
 15 
 16 long long ans[N];
 17 
 18 char num[1000005];
 19 int w1,w2,n;
 20 
 21 complex<double>omege[N],reomege[N],d[N],c[N];
 22 
 23 void init()
 24 {
 25     scanf("%s",num);
 26     int l=strlen(num);
 27     while(l>4)
 28     {
 29         a[w1]=(num[l-4]-'0')*1000+(num[l-3]-'0')*100+(num[l-2]-'0')*10+num[l-1]-'0';
 30         l-=4;
 31         c[w1]=a[w1];
 32         w1++; 
 33     }
 34     if(l)
 35     {
 36         for(int i=0;i<l;++i)
 37             a[w1]=a[w1]*10+num[i]-'0';
 38         c[w1]=a[w1];
 39         w1++;
 40     }
 41     scanf("%s",num);
 42     l=strlen(num);
 43     while(l>4)
 44     {
 45         b[w2]=(num[l-4]-'0')*1000+(num[l-3]-'0')*100+(num[l-2]-'0')*10+num[l-1]-'0';
 46         d[w2]=b[w2];
 47         l-=4;
 48         w2++;
 49     }
 50     if(l)
 51     {
 52         for(int i=0;i<l;++i)
 53             b[w2]=b[w2]*10+num[i]-'0';
 54         d[w2]=b[w2];
 55         w2++;
 56     }
 57     n=w1+w2-1;
 58     n=1<<((int)ceil(log(n)/log(2)));
 59     for(int i=0;i<n;++i)//*
 60     {
 61         omege[i]=complex<double>(cos(2*PI/n*i),sin(2*PI/n*i));
 62         reomege[i]=conj(omege[i]);
 63     }
 64 }
 65 
 66 inline void FFT(complex<double> *a,complex<double> *omege)//*
 67 {
 68     for(int i=0,j=0;i<n;++i)
 69     {
 70         if(i>j) swap(a[i],a[j]);
 71         for(int l=n>>1;(j^=l)<l;l>>=1);//反向加1
 72     }
 73     complex<double> k;
 74     for(int l=2;l<=n;l<<=1)
 75     {
 76         for(int i=0;i<n;i+=l)
 77             for(int j=i;j<i+l/2;++j)
 78             {
 79                 k=omege[n/l*(j-i)]*a[j+l/2];
 80                 a[j+l/2]=a[j]-k;
 81                 a[j]=a[j]+k;
 82             }
 83     }
 84 }
 85 
 86 int main()
 87 {
 88     init();
 89     FFT(c,omege);
 90     FFT(d,omege);
 91     for(int i=0;i<n;++i)//*
 92         c[i]*=d[i];
 93     FFT(c,reomege);
 94     for(int i=0;i<n;++i)//*
 95         c[i]/=n;
 96     int l=w1+w2-1;
 97     for(int i=0;i<l;++i)
 98     {
 99         ans[i]+=c[i].real()+0.5;
100         if(ans[i]>=10000)
101         {
102             ans[i+1]+=ans[i]/10000;
103             ans[i]%=10000;
104         }
105     }
106     while(!ans[l]&&l)
107         --l;
108     printf("%d",ans[l--]);
109     while(l>=0)
110         printf("%04d",ans[l--]);
111     return 0;
112 }