hdu1402 FFT入门

参考这里：http://www.cnblogs.com/pdev/p/4354705.html

题意：求大数乘法A*B

A和B位数很长。裸高精度时间复杂度是O(nm)，会完蛋

不妨回忆下裸高精度的过程：

其实乘法的那一步很类似前面介绍过的多项式快速乘法诶(⊙▽⊙)

所以就可以用前述方法计算咯，时间复杂度O(nlogn)

我是这样理解的：

每个乘数都是都是一坨时域信号（一个大混合物），然后对乘数分别进行DFT（Discrete Fourier Transform）得到频域信号（一堆纯净物）。

然后对纯净物按类别分别相加，就得到了新信号（这里即乘法结果）对应的频域信号（一堆纯净物）

然后再来一次IDFT（Inverse DFT，逆变换）把频域再转成时域（一个大混合物，即真正的乘法结果）就好啦

  1 //HDU 1402求高精度乘法
  2 #include  <stdio.h>
  3 #include  <string.h>
  4 #include  <iostream>
  5 #include  <algorithm>
  6 #include  <math.h>
  7 using namespace  std;
  8 const double PI = acos(-1.0);
  9 //复数结构体
 10 struct  Complex
 11 {
 12     double x,y;//实部和虚部  x+yi
 13     Complex(double _x = 0.0,double _y = 0.0)
 14     {
 15         x = _x;
 16         y = _y;
 17     }
 18     Complex operator -(const Complex &b)const
 19     {
 20         return  Complex(x-b.x,y-b.y);
 21     }
 22     Complex operator +(const Complex &b)const
 23     {
 24         return  Complex(x+b.x,y+b.y);
 25     }
 26     Complex operator *(const Complex &b)const
 27     {
 28         return  Complex(x*b.x-y*b.y,x*b.y+y*b.x);
 29     }
 30 };
 31 /*
 32 *进行FFT和IFFT前的反转变换。
 33 *位置i和（i二进制反转后位置）互换
 34 * len必须去2的幂
 35 */
 36 void change(Complex y[],int  len)
 37 {
 38     int  i,j,k;
 39     for(i = 1, j = len/2; i <len-1; i++)
 40     {
 41         if(i < j)swap(y[i],y[j]);
 42 //交换互为小标反转的元素，i<j保证交换一次
 43 //i做正常的+1，j左反转类型的+1,始终保持i和j是反转的
 44         k = len/2;
 45         while(j >= k)
 46         {
 47             j -= k;
 48             k /= 2;
 49         }
 50         if(j < k)j += k;
 51     }
 52 }
 53 /*
 54 *做FFT
 55 * len必须为2^k形式，
 56 * on==1时是DFT，on==-1时是IDFT
 57 */
 58 //fft(x,len,1):对向量x做DFT（时域->频域），向量长度为1--len
 59 //fft(x,len,-1):做IDFT（频域->时域）
 60 void fft(Complex y[],int len,int  on)
 61 {
 62     change(y,len);
 63     for(int h = 2; h <= len; h <<= 1)
 64     {
 65         Complex wn(cos(-on*2*PI/h),sin(-on*2*PI/h));
 66         for(int j = 0; j < len; j+=h)
 67         {
 68             Complex  w(1,0);
 69             for(int k = j; k < j+h/2; k++)
 70             {
 71                 Complex u = y[k];
 72                 Complex t = w*y[k+h/2];
 73                 y[k] = u+t;
 74                 y[k+h/2] = u-t;
 75                 w = w*wn;
 76             }
 77         }
 78     }
 79     if(on == -1)
 80         for(int i = 0; i < len; i++)
 81             y[i].x /= len;
 82 }
 83 const int MAXN = 200010;
 84 Complex  x1[MAXN],x2[MAXN];
 85 char  str1[MAXN/2],str2[MAXN/2];
 86 int  sum[MAXN];
 87 int main()
 88 {
 89     while(scanf("%s%s",str1,str2)==2)
 90     {
 91         int len1 = strlen(str1);
 92         int len2 = strlen(str2);
 93         int len = 1;
 94         while(len < len1*2 || len < len2*2)len<<=1;
 95         for(int i = 0; i < len1; i++)
 96             x1[i] =  Complex(str1[len1-1-i]-‘0‘,0);
 97         for(int i = len1; i < len; i++)
 98             x1[i] =  Complex(0,0);
 99         for(int i = 0; i < len2; i++)
100             x2[i] =  Complex(str2[len2-1-i]-‘0‘,0);
101         for(int i = len2; i < len; i++)
102             x2[i] =  Complex(0,0);
103 //x1[i]:x1对应的向量
104 //例如1989就是(9,0)、(8,0)、(9,0)、(1,0)、(0,0)、...
105
106 //求DFT
107         fft(x1,len,1);
108         fft(x2,len,1);
109         for(int i = 0; i < len; i++)
110             x1[i] = x1[i]*x2[i];
111         fft(x1,len,-1);
112         for(int i = 0; i < len; i++)
113             sum[i] = (int)(x1[i].x+0.5);
114         for(int i = 0; i < len; i++)
115         {
116             sum[i+1]+=sum[i]/10;
117             sum[i]%=10;
118         }
119         len = len1+len2-1;
120         while(sum[len] <= 0 && len > 0)len--;
121         for(int i = len; i >= 0; i--)
122             printf("%c",sum[i]+‘0‘);
123         printf("\n");
124     }
125     return  0;
126 }