[BZOJ2179]FFT快速傅立叶
[BZOJ2179]FFT快速傅立叶
题目大意:
求\(a\times b(1\le a,b\le10^{60000})\)。
思路:
FFT模板。
源代码:
#include<cmath>
#include<cstdio>
#include<cctype>
#include<complex>
#include<algorithm>
inline int getint() {
register char ch;
while(!isdigit(ch=getchar()));
register int x=ch^'0';
while(isdigit(ch=getchar())) x=(((x<<2)+x)<<1)+(ch^'0');
return x;
}
inline int getdigit() {
register char ch;
while(!isdigit(ch=getchar()));
return ch^'0';
}
typedef std::complex<double> complex;
const double pi=M_PI;
const int N=131075;
int lim,ans[N];
complex a[N],b[N],c[N],omega[N],iomega[N];
inline void init_omega(const int &n) {
for(register int i=0;i<n;i++) {
omega[i]=(complex){cos(2*pi*i/lim),sin(2*pi*i/lim)};
iomega[i]=conj(omega[i]);
}
}
inline void fft(complex f[],complex w[],const int &n) {
for(register int i=0,j=0;i<n;i++) {
if(i>j) std::swap(f[i],f[j]);
for(register int l=n>>1;(j^=l)<l;l>>=1);
}
for(register int i=2;i<=n;i<<=1) {
const int m=i>>1;
for(register int j=0;j<n;j+=i) {
for(register int k=0;k<m;k++) {
const complex z=f[j+m+k]*w[n/i*k];
f[j+m+k]=f[j+k]-z;
f[j+k]+=z;
}
}
}
}
int main() {
const int n=getint();
for(register int i=0;i<n;i++) a[i]=getdigit();
for(register int i=0;i<n;i++) b[i]=getdigit();
std::reverse(&a[0],&a[n]);
std::reverse(&b[0],&b[n]);
for(lim=1;lim<n;lim<<=1);lim<<=1;
init_omega(lim);
fft(a,omega,lim);
fft(b,omega,lim);
for(register int i=0;i<lim;i++) {
c[i]=a[i]*b[i];
}
fft(c,iomega,lim);
for(register int i=0;i<lim;i++) {
ans[i]=round(c[i].real()/lim);
}
int len=0;
for(register int i=0;i<lim;i++) {
if(ans[i]) len=i;
ans[i+1]+=ans[i]/10;
ans[i]%=10;
}
for(register int i=len;i>=0;i--) {
printf("%d",ans[i]);
}
return 0;
}