HDU

题意:计算A*B,A,B均为长度小于50000的整数。

这是FFT在大整数相乘中的一个应用,我本来想用NTT做的,但NTT由于取模很可能取炸,所以base必须设得很小,而且效率也比不上FFT。

A和B的存储均用long long,在计算乘积的时候转化成double,计算完成后再转回来即可。

测得base在精度允许范围内最多能开到10000。

把平方和快速幂的函数也写上了,可以当模板用~

  1 #include<bits/stdc++.h>
  2 using namespace std;
  3 typedef long long ll;
  4 typedef double db;
  5 const int N=1e5+10,inf=0x3f3f3f3f;
  6 const db pi=acos(-1);
  7 const ll P[]= {1,10,100,1000,10000};
  8 struct cpl {
  9     db x,y;
 10     cpl operator-(cpl& b) {return {x-b.x,y-b.y};}
 11     cpl operator+(cpl& b) {return {x+b.x,y+b.y};}
 12     cpl operator*(cpl& b) {return {x*b.x-y*b.y,x*b.y+y*b.x};}
 13 } A[N],B[N];
 14 void change(cpl* a,int n) {
 15     for(int i=1,j=n>>1,k; i<n-1; ++i) {
 16         if(i<j)swap(a[i],a[j]);
 17         k=n>>1;
 18         while(j>=k)j-=k,k>>=1;
 19         j+=k;
 20     }
 21 }
 22 void FFT(cpl* a,int n,int f) {
 23     change(a,n);
 24     for(int k=1; k<n; k<<=1) {
 25         cpl wn= {cos(pi*f/k),sin(pi*f/k)};
 26         for(int i=0; i<n; i+=k<<1) {
 27             cpl w{1,0};
 28             for(int j=i; j<i+k; ++j) {
 29                 cpl x=a[j],y=w*a[j+k];
 30                 a[j]=x+y,a[j+k]=x-y;
 31                 w=w*wn;
 32             }
 33         }
 34     }
 35     if(!~f)for(int i=0; i<n; ++i)a[i].x/=n,a[i].y/=n;
 36 }
 37 struct Bigint {
 38     static const int N=1e5+10;
 39     static const int bit=4;
 40     static const ll base=pow(10,bit);
 41     int n;
 42     ll a[N];
 43     ll& operator[](int x) {return a[x];}
 44     void init(ll x) {
 45         memset(a,0,sizeof a);
 46         a[0]=x,n=0;
 47     }
 48     void init(char* s) {
 49         memset(a,0,sizeof a);
 50         int m=strlen(s);
 51         for(int i=0; i<m; ++i)a[(m-1-i)/bit]+=(s[i]-'0')*P[(m-1-i)%bit];
 52         n=(m-1)/bit;
 53     }
 54     void Sqr() {
 55         int m;
 56         for(m=1; m<=n*2; m<<=1);
 57         for(int i=0; i<m; ++i)A[i]= {a[i],0};
 58         FFT(A,m,1);
 59         for(int i=0; i<m; ++i)A[i]=A[i]*A[i];
 60         FFT(A,m,-1);
 61         for(int i=0; i<m; ++i)a[i]=A[i].x+0.5;
 62         for(int i=0; i<m; ++i) {
 63             if(a[i]>=base) {
 64                 a[i+1]+=a[i]/base;
 65                 a[i]%=base;
 66                 if(i==m-1)++m;
 67             }
 68         }
 69         for(n=m-1; n>1&&!a[n]; --n);
 70     }
 71     void Mul(Bigint& b) {
 72         int m;
 73         for(m=1; m<=n*2||m<=b.n*2; m<<=1);
 74         for(int i=0; i<m; ++i)A[i]= {a[i],0},B[i]= {b[i],0};
 75         FFT(A,m,1),FFT(B,m,1);
 76         for(int i=0; i<m; ++i)A[i]=A[i]*B[i];
 77         FFT(A,m,-1);
 78         for(int i=0; i<m; ++i)a[i]=A[i].x+0.5;
 79         for(int i=0; i<m; ++i) {
 80             if(a[i]>=base) {
 81                 a[i+1]+=a[i]/base;
 82                 a[i]%=base;
 83                 if(i==m-1)++m;
 84             }
 85         }
 86         for(n=m-1; n>1&&!a[n]; --n);
 87     }
 88     void Pow(int p) {
 89         Bigint b=*this;
 90         this->init(1);
 91         for(; p; p>>=1,b.Sqr())if(p&1)this->Mul(b);
 92     }
 93     void pr() {
 94         for(int i=n; i>=0; --i) {
 95             if(i==n)printf("%lld",a[i]);
 96             else printf("%04lld",a[i]);
 97         }
 98         puts("");
 99     }
100 } a,b;
101 char s1[N],s2[N];
102 
103 int main() {
104     while(scanf("%s%s",s1,s2)==2) {
105         a.init(s1),b.init(s2);
106         a.Mul(b);
107         a.pr();
108     }
109     return 0;
110 }
原文地址:https://www.cnblogs.com/asdfsag/p/10492859.html