标签:
没什么好说的,今天有考了FFT(虽然不用FFT也能过)但是确实有忘了怎么写FFT了,于是乎只有重新写一遍FFT模板练一下手了。第一部分普通FFT,第二部分数论FFT,记一下模数2^23*7*17+1
#include<iostream> #include<cstdio> #include<cstring> #include<algorithm> #include<cmath> using namespace std; typedef double real; namespace task1 { const int maxn = 110000; const real pi = 3.1415926535897932; struct Complex { real x,y; Complex(){} Complex(real x,real y=0):x(x),y(y){} }; Complex operator + (Complex c1,Complex c2) { return Complex(c1.x+c2.x,c1.y+c2.y); } Complex operator - (Complex c1,Complex c2) { return Complex(c1.x-c2.x,c1.y-c2.y); } Complex operator * (Complex c1,Complex c2) { return Complex(c1.x*c2.x-c1.y*c2.y,c1.y*c2.x+c1.x*c2.y); } Complex operator / (Complex c1,real k) { return Complex(c1.x/k,c1.y/k); } char s1[maxn],s2[maxn]; Complex g1[maxn],g2[maxn],rr[maxn]; int lst[maxn]; void FFT(Complex g[],int l,int p) { int x; for (int i=0;i<l;i++) { x=0; for (int j=1,k=(l>>1);j<l;j<<=1,k>>=1) if (i&j)x+=k; if (x<i) swap(g[x],g[i]); } Complex wt,w,tmp; for (int i=1;i<l;i<<=1) { wt=Complex(cos(pi/i*p),sin(pi/i*p)); for (int j=0;j<l;j+=(i<<1)) { Complex w=Complex(1,0); for (int k=0;k<i;k++) { tmp=g[j+k]; g[j+k]=g[j+k]+g[i+j+k]*w; g[i+j+k]=tmp-g[i+j+k]*w; w=w*wt; } } } } void main() { int l1,l2; scanf("%s",s1); scanf("%s",s2); l1=(int)strlen(s1); l2=(int)strlen(s2); for (int i=0;i<l1;i++) g1[(l1-i-1)/4]=g1[(l1-i-1)/4].x*10+s1[i]-‘0‘; for (int i=0;i<l2;i++) g2[(l2-i-1)/4]=g2[(l2-i-1)/4].x*10+s2[i]-‘0‘; int l=max(l1,l2)/4+1; while (l!=(l&(-l)))l-=l&(-l); l<<=2; FFT(g1,l,1); FFT(g2,l,1); for (int i=0;i<l;i++)rr[i]=g1[i]*g2[i]; FFT(rr,l,-1); for (int i=0;i<l;i++) { rr[i]=rr[i]/l; lst[i]=(int)(rr[i].x+0.5); } for (int i=0;i<l;i++) { lst[i+1]+=lst[i]/10000; lst[i]%=10000; } while (l>1 && !lst[l-1])l--; printf("%d",lst[l-1]); for (int i=l-2;i>=0;i--) printf("%04d",lst[i]); printf("\n"); } } namespace task2 { typedef long long qword; const int maxn = 10000; const qword p = 998244353; char s1[maxn],s2[maxn]; qword g1[maxn],g2[maxn];; qword rr[maxn]; qword pow_mod(qword x,qword y) { qword ret=1; while (y) { if (y&1)ret=ret*x%p; x=x*x%p; y>>=1; } return ret; } void FFT(qword g[],int l,int f) { int x=0; for (int i=1;i<l;i++) { x=0; for (int j=1,k=(l>>1);j<l;j<<=1,k>>=1) if (i&j)x+=k; if (i<x) swap(g[i],g[x]); } qword w,wt,tmp; for (int i=1,t=1;i<l;i<<=1,t++) { wt=pow_mod(3,7*17*(1<<(23-t)))%p; if (f==-1) //wt=-wt; wt=pow_mod(wt,p-2); for (int j=0;j<l;j+=(i<<1)) { w=1; for (int k=0;k<i;k++) { tmp=g[j+k]; g[j+k]=(g[j+k]+g[i+j+k]*w)%p; g[i+j+k]=(tmp-g[i+j+k]*w)%p; w=w*wt%p; } } } } void main() { int l1,l2; // for (int i=0;i<5;i++) // printf("[1/%d*pi]:%d\n",(1<<i),7*17*(1<<(23-i))); scanf("%s",s1); scanf("%s",s2); l1=(int)strlen(s1); l2=(int)strlen(s2); for (int i=0;i<l1;i++) g1[i]=s1[l1-i-1]-‘0‘; for (int i=0;i<l2;i++) g2[i]=s2[l2-i-1]-‘0‘; int l=max(l1,l2); while (l!=(l&(-l)))l-=l&(-l); l<<=2; FFT(g1,l,1); FFT(g2,l,1); for (int i=0;i<l;i++)rr[i]=g1[i]*g2[i]%p; FFT(rr,l,-1); qword invl=pow_mod(l,p-2); for (int i=0;i<l;i++)rr[i]=(rr[i]*invl%p+p)%p; for (int i=0;i<l;i++) { rr[i+1]+=rr[i]/10; rr[i]%=10; } while (l>1 && ! rr[l-1])l--; for (int i=l-1;i>=0;i--) printf("%d",(int)rr[i]); } } int main() { freopen("input.txt","r",stdin); task1::main(); task2::main(); }
标签:
原文地址:http://www.cnblogs.com/mhy12345/p/4531821.html