hdu 1402 A * B Problem Plus 快速傅里叶变换

A * B Problem Plus

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 65536/32768 K (Java/Others)
Total Submission(s): 14598    Accepted Submission(s): 2741

Problem Description

Calculate A * B.

 

Input

Each line will contain two integers A and B. Process to end of file.

Note: the length of each integer will not exceed 50000.

 

Output

For each case, output A * B in one line.

 

Sample Input

  
  

   
   1
2
1000
2
  
  

 

Sample Output

  
  

   
   2
2000
  
  

 

Author

DOOM III

 

题目：

求A*B的的结果，长度为5000，如果直接模拟乘法，需要5000*5000的复杂度，会超时。

解析：

一个数anan-1................a0是由n+1个数字构成的

可以看出an*10^n + an-1*10^n-1,......,+a0*1把令x = 10

则  变为f(x) = an*x^n+........+a0*x^0

两个数相乘就可以看成系数不同的f(x)相乘。

变成 ----------  A(x)*B(x)       ------------- 然后ax^n * bx^m = a*b*x^(n+m)可以看成这两个位置上的数字相乘。得到a*b，是没有进位的。所以最后的结果要进行进位处理

代码如下：fft是模板题，想要了解详情看算法导论第七部分，30章。但是需要一点线性代数，复数，数值分析的知识。

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
using namespace std;
struct Complex{
    double real;
    double imag;
    Complex (){}
    Complex(double the){
        real = cos(the);
        imag = sin(the);
    }
    Complex(double a,double b){
        real = a;
        imag = b;
    }
};
Complex operator + (Complex a,Complex b){
    return Complex(a.real+b.real,a.imag+b.imag);
}
Complex operator - (Complex a,Complex b){
    return Complex(a.real-b.real,a.imag-b.imag);
}
Complex operator * (Complex a,Complex b){
    return Complex(a.real*b.real-a.imag*b.imag,a.imag*b.real+a.real*b.imag);
}

#define maxn 140000
double pi2 = 2*acos(-1.0);
void fft(Complex* A,int len, int ref){
    //A[rev[k]] = ak 系数向量转置
    int bitn = log2(len);
    int i,j,k;
    for(i = 0;i < len; i++){
        k = 0;
        for( j = 0;j < bitn; j++){
            k = (k<<1);
            if(i&(1<<j))
                k |= 1;
        }
        if(k > i)
            swap(A[i],A[k]);
    }
    //fft计算得到点值
    Complex wm,w,t,u;
    for(int m = 2,f=1; m <= len; m<<=1,f<<=1){
        wm = Complex(ref*pi2/m);
        for( k = 0;k < len; k += m){
            w = Complex(1.0,0);
            for( j = k; j < k+f; j++){
                t = w*A[j+f];
                u = A[j];
                A[j] = u+t;
                A[j+f] = u-t;
                w = w*wm;
            }
        }
    }
    if(ref == -1){
        for( i = 0;i < len; i++){
            A[i].real = A[i].real/len;
        }
    }
}
char word1[maxn];
char word2[maxn];
Complex a1[maxn];
Complex a2[maxn];
int result[maxn];
int main(){
    while(gets(word1)){
        gets(word2);
        int len1 = strlen(word1);
        int len2 = strlen(word2);
        int len3 = 1;
//计算的长度必须为2的n次幂
        while(len3 < len1+len2)
            len3 <<= 1;
        //计算第一个多项式在个点的值
        for(int i = len1; i < len3; i++)
            a1[i] = Complex(0,0);
        for(int i = 0;i < len1; i++)
            a1[i] = Complex(word1[len1-i-1]-'0',0);
        fft(a1,len3,1);
        //计算第二个多项式在个点的值
        for(int i = len2; i < len3; i++)
            a2[i] = Complex(0,0);
        for(int i = 0;i < len2; i++)
            a2[i] = Complex(word2[len2-i-1]-'0',0);
        fft(a2,len3,1);
        //计算两个多项式在各个点的乘积
        for(int i = 0;i < len3; i++)
            a1[i] = a1[i]*a2[i];
        //逆向fft求多项式的系数
        fft(a1,len3,-1);
        //每个数的实部就是多项式的系数，由于精度问题需要+0.5
        int flag = 0;
        for(int i = 0;i < len3; i++){
            result[i] = int(a1[i].real+0.5);
        }
        //十进制数每位都是0-9的数字，需要进位
        for(int i = 0;i < len3; i++){
            result[i+1] += result[i]/10;
            result[i] %= 10;
        }
        for(flag = len3-1; flag > 0 && result[flag]==0;flag--);
        for(;flag>=0;flag--){
            putchar(result[flag]+'0');
        }
        printf("\n");
    }
    return 0;
}