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

Problem:
给定两个大整数,计算两者的积。
Solution:
利用傅里叶变换,把10,100理解为阶,每一位理解为系数。

//fft
#include<cstdio>
#include<iostream>
#include<sstream>
#include<cstdlib>
#include<cmath>
#include<cctype>
#include<string>
#include<cstring>
#include<algorithm>
#include<stack>
#include<queue>
#include<set>
#include<map>
#include<ctime>
#include<vector>
#include<fstream>
#include<list>
using namespace std;

#define ms(s) memset(s,0,sizeof(s))
typedef unsigned long long ULL;
typedef long long LL;

const double PI = 3.141592653589;
const int INF = 0x3fffffff;
const int maxn = 200010;

struct Complex{
    double a,b;
    Complex(double _a = 0.0, double _b = 0.0){
        a = _a;  b = _b;
    }
    Complex operator - (const Complex& rhs) const{
        return Complex(a-rhs.a, b-rhs.b);
    }
    Complex operator + (const Complex& rhs) const{
        return Complex(a+rhs.a, b+rhs.b);
    }
    Complex operator * (const Complex& rhs) const{
        return Complex(a*rhs.a-b*rhs.b, a*rhs.b+b*rhs.a);
    }
};

void bit_reverse(Complex y[], int len){
    int i = 1, j = len/2,k;
    for(; i < len-1; ++i){
        if(i < j)
            swap(y[i],y[j]);
        k = len/2;
        while(j >= k){
            j -= k;  k /= 2;
        }
        if(j < k)
            j += k;
    }
}

void fft(Complex y[], int len, int on){
    bit_reverse(y,len);
    for(int h = 2; h <= len; h <<= 1){
        Complex wn(cos(-on*2*PI/h),sin(-on*2*PI/h));
        for(int j = 0; j < len; j += h){
            Complex w(1,0);
            for(int k = j; k < j+h/2; ++k){
                Complex u = y[k];
                Complex t = w*y[k+h/2];
                y[k] = u+t;
                y[k+h/2] = u-t;
                w = w*wn;
            }
        }
    }
    if(on == -1){
        for(int i = 0; i < len; ++i)
            y[i].a /= len;
    }
}

char str1[maxn/2];
char str2[maxn/2];
Complex x1[maxn],x2[maxn];
int sum[maxn];

int main(){
//    freopen("F:\\input.txt","r",stdin);
//    freopen("F:\\output.txt","w",stdout);
//    ios::sync_with_stdio(false);

    while(~scanf("%s%s",str1,str2)){
        int len1 = strlen(str1);
        int len2 = strlen(str2);
        int len = 1;
        while(len<len1*2 || len<len2*2)
            len <<= 1;//find the minimum 2^n that >= max(len1,len2)
        for(int i = 0; i < len1; ++i)
            x1[i] = Complex(str1[len1-i-1]-'0',0);
        for(int i = len1; i < len; ++i)
            x1[i] = Complex(0,0);
        for(int i = 0; i < len2; ++i)
            x2[i] = Complex(str2[len2-i-1]-'0',0);
        for(int i = len2; i < len; ++i)
            x2[i] = Complex(0,0);
        fft(x1,len,1);
        fft(x2,len,1);
        for(int i = 0; i < len; ++i){
            x1[i] = x1[i]*x2[i];
        }
        fft(x1,len,-1);
        for(int i = 0; i < len; ++i)
            sum[i] = (int)round(x1[i].a);
        for(int i = 0; i < len; ++i){
            sum[i+1] += sum[i]/10;
            sum[i] %= 10;
        }
        len = len1 + len2 -1;
        while(sum[len]==0 && len>0)
            len--;
        for(int i = len; i >= 0; --i)
            printf("%c",sum[i]+'0');
        printf("\n");
    }
    return 0;
}
//ntt
#include<cstdio>
#include<iostream>
#include<sstream>
#include<cstdlib>
#include<cmath>
#include<cctype>
#include<string>
#include<cstring>
#include<algorithm>
#include<stack>
#include<queue>
#include<set>
#include<map>
#include<ctime>
#include<vector>
#include<fstream>
#include<list>
using namespace std;

#define ms(s) memset(s,0,sizeof(s))
typedef unsigned long long ULL;
typedef long long LL;

const double PI = 3.141592653589;
const int INF = 0x3fffffff;

char aa[50005],bb[50005];

const long long P = 50000000001507329LL; //190734863287 * 2 ^ 18 + 1
const int G = 3;

long long wn[25];

long long mul(long long x, long long y) {
    return (x*y - (long long)(x / (long double)P*y + 1e-3) * P + P) % P;
}

long long qpow(long long x, long long k) {
    long long ret = 1;
    while(k) {
        if(k & 1) ret = mul(ret, x);
        k >>= 1;
        x = mul(x, x);
    }
    return ret;
}

void getwn() {
    for(int i = 1; i <= 18; ++i) {
        int t = 1 << i;
        wn[i] = qpow(G, (P - 1)/t);
    }
}
int len;
void change(long long *y, int len) {
    for(int i = 1, j = len/2; i < len - 1; ++i) {
        if(i < j) swap(y[i], y[j]);
        int k = len/2;
        while(j >= k) {
            j -= k;
            k /= 2;
        }
        j += k;
    }
}

void NTT(long long y[], int on) {
    change(y, len);
    int id = 0;

    for(int h = 2; h <= len; h <<= 1) {
        ++id;
        for(int j = 0; j < len; j += h) {
            long long w = 1;
            for(int k = j; k < j + h / 2; ++k) {
                long long u = y[k];
                long long t = mul(y[k+h/2], w);
                y[k] = u + t;
                if(y[k] >= P) y[k] -= P;
                y[k+h/2] = u - t + P;
                if(y[k+h/2] >= P) y[k+h/2] -= P;
                w = mul(w, wn[id]);
            }
        }
    }
    if(on == -1) {
        for(int i = 1; i < len / 2; ++i) swap(y[i], y[len-i]);
        long long inv = qpow(len, P - 2);
        for(int i = 0; i < len; ++i)
            y[i] = mul(y[i], inv);
    }
}
void Convolution(long long A[],long long B[],int n){
    for(len=1; len<(n<<1); len<<=1);
    for(int i=n; i<len; ++i){
        A[i]=B[i]=0;
    }

    NTT(A,1); NTT(B,1);
    for(int i=0; i<len; ++i){
        A[i]=mul(A[i],B[i]);
    }
    NTT(A,-1);
}

long long A[220000];
long long B[220000];
int main(){
//    freopen("D:\\input.txt","r",stdin);
//    freopen("D:\\output.txt","w",stdout);
//    ios::sync_with_stdio(false);
    getwn();
    int idx;
    int la,lb;
    while(~scanf("%s%s",aa,bb)){
        ms(A);  ms(B);
        la = strlen(aa);  lb = strlen(bb);
        idx = la-1;
        for(int i = 0; idx >= 0; ++i){
            A[i] = aa[idx]-'0';
            idx--;
        }
        idx = lb-1;
        for(int i = 0; idx >= 0; ++i){
            B[i] = bb[idx]-'0';
            idx--;
        }
        if(la >= lb){
            for(int i = lb; i < la; ++i)
                B[i] = 0;
        }
        else{
            for(int i = la; i < lb; ++i)
                A[i] = 0;
        }
        int mx = max(la,lb);
        Convolution(A,B,mx);
         len--;
        for(int i = 0; i <= len; ++i){
            int t = 0;
            if(A[i] > 9){
                t = A[i]/10;
                A[i] %= 10;
                A[i+1] += t;
            }
        }
        while(A[len] == 0LL && len >= 0){
            len--;
        }
        if(len < 0){
            printf("0");
        }
        else{
            while(len >= 0){
                printf("%lld",A[len]);
                len--;
            }
        }
        printf("\n");
    }


    return 0;
}
  • 1
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值