FFT:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 |
struct Complex { double real, image; Complex( double _real, double _image) { real = _real; image = _image; } Complex(){} }; Complex operator + (const Complex &c1, const Complex &c2) { return Complex(c1.real + c2.real, c1.image + c2.image); } Complex operator - (const Complex &c1, const Complex &c2) { return Complex(c1.real - c2.real, c1.image - c2.image); } Complex operator * (const Complex &c1, const Complex &c2) { return Complex(c1.real*c2.real - c1.image*c2.image, c1.real*c2.image + c1.image*c2.real); } int rev(int id, int len) { int ret = 0; for(int i = 0; (1 << i) < len; i++) { ret <<= 1; if(id & (1 << i)) ret |= 1; } return ret; } Complex A[1 << 19]; void FFT(Complex *a, int len, int DFT) { for(int i = 0; i < len; i++) A[rev(i, len)] = a[i]; for(int s = 1; (1 << s) <= len; s++) { int m = (1 << s); Complex wm = Complex(cos(DFT*2*PI/m), sin(DFT*2*PI/m)); for(int k = 0; k < len; k += m) { Complex w = Complex(1, 0); for(int j = 0; j < (m >> 1); j++) { Complex t = w*A[k + j + (m >> 1)]; Complex u = A[k + j]; A[k + j] = u + t; A[k + j + (m >> 1)] = u - t; w = w*wm; } } } if(DFT == -1) for(int i = 0; i < len; i++) A[i].real /= len, A[i].image /= len; for(int i = 0; i < len; i++) a[i] = A[i]; return; } |
NTT
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 |
#include <iostream> #include <cstdio> #include <cstring> #include <algorithm> #include <cmath> using namespace std; const int maxn=300020; const int P=880803841; const int G=26; int 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; } int qpow(int x,int k,int p) { int ret=1; while(k) { if (k&1) ret=1LL*ret*x%p; k>>=1; x=1LL*x*x%p; } return ret; } void getwn() { for (int i=1;i<=21;i++) { int t=1<<i; wn[i]=qpow(G,(P-1)/t,P); } } void change(int *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(int *y,int len,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) { int w=1; for (int k=j;k<j+h/2;k++) { int u=y[k]; int t=1LL*y[k+h/2]*w%P; 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=1LL*w*wn[id]%P; } } } if (on==-1) { for (int i=1;i<len/2;i++) { swap(y[i],y[len-i]); } int inv=qpow(len,P-2,P); for (int i=1;i<len;i++) { y[i]=1LL*y[i]*inv%P; } } } int main() { getwn(); return 0; } |