[关闭]
@zzzc18 2017-12-03T05:57:35.000000Z 字数 3207 阅读 1071

模板

FFT


见miskcoo博客

FFT

  1. #include<cmath>
  2. #include<cstdio>
  3. #include<cstring>
  4. #include<algorithm>
  5. using namespace std;
  6. namespace Fast_Fourier_Transform{
  7. const int MAXN = 300000+9;
  8. const double PI = acos(-1.0);
  9. int size,bit_length;
  10. int loc[MAXN];
  11. struct C{
  12. double x,y;
  13. C(double _x=0,double _y=0):x(_x),y(_y){}
  14. };
  15. C operator + (const C &A,const C &B){
  16. return C(A.x+B.x,A.y+B.y);
  17. }
  18. C operator - (const C &A,const C &B){
  19. return C(A.x-B.x,A.y-B.y);
  20. }
  21. C operator * (const C &A,const C &B){
  22. return C(A.x*B.x-A.y*B.y,A.x*B.y+A.y*B.x);
  23. }
  24. C A1[MAXN],B1[MAXN];
  25. void INIT(int len){
  26. for(size=1,bit_length=-1;size<len;size<<=1,bit_length++);
  27. int now=0;
  28. for(int i=0;i<size;i++){
  29. loc[i]=now;
  30. for(int j=1<<bit_length;(now^=j)<j;j>>=1);
  31. }
  32. }
  33. void DFT(int *A,C *re){
  34. for(int i=0;i<size;i++)
  35. re[i].x=A[loc[i]];
  36. for(int k=2;k<=size;k<<=1){
  37. int len=k>>1;
  38. C Wn(cos(2.0*PI/(1.0*k)),sin(2.0*PI/(1.0*k)));
  39. for(int i=0;i<size;i+=k){
  40. C W(1,0);
  41. for(int j=0;j<len;j++,W=W*Wn){
  42. C u=re[i+j],v=W*re[i+j+len];
  43. re[i+j]=u+v;
  44. re[i+j+len]=u-v;
  45. }
  46. }
  47. }
  48. }
  49. void IDFT(C *A,C *re){
  50. for(int i=0;i<size;i++)
  51. re[i]=A[loc[i]];
  52. for(int k=2;k<=size;k<<=1){
  53. int len=k>>1;
  54. C Wn(cos(-2.0*PI/(1.0*k)),sin(-2.0*PI/(1.0*k)));
  55. for(int i=0;i<size;i+=k){
  56. C W(1,0);
  57. for(int j=0;j<len;j++,W=W*Wn){
  58. C u=re[i+j],v=W*re[i+j+len];
  59. re[i+j]=u+v;
  60. re[i+j+len]=u-v;
  61. }
  62. }
  63. }
  64. for(int i=0;i<size;i++)
  65. re[i].x/=1.0*size;
  66. }
  67. void FFT(int *A,int *B,int *ans){
  68. DFT(A,A1);DFT(B,B1);
  69. for(int i=0;i<size;i++)
  70. A1[i]=A1[i]*B1[i];
  71. IDFT(A1,B1);
  72. for(int i=0;i<size;i++)
  73. ans[i]=(int)(B1[i].x+0.5);
  74. }
  75. }
  76. using namespace Fast_Fourier_Transform;
  77. int num1[MAXN],num2[MAXN];
  78. int ANS[MAXN];
  79. int main(){
  80. int n,m;
  81. scanf("%d%d",&n,&m);
  82. INIT(n+m+3);
  83. for(int i=0;i<=n;i++)
  84. scanf("%d",&num1[i]);
  85. for(int i=0;i<=m;i++)
  86. scanf("%d",&num2[i]);
  87. FFT(num1,num2,ANS);
  88. for(int i=0;i<=m+n;i++)
  89. printf("%d ",ANS[i]);
  90. puts("");
  91. return 0;
  92. }

NTT

  1. #include<cstdio>
  2. #include<cstring>
  3. #include<algorithm>
  4. using namespace std;
  5. const int MAXN = 300000;
  6. class Fast_Number_Theory_Transform{
  7. private:
  8. const int MOD,G;
  9. int size,bit_length;
  10. int loc[MAXN];
  11. int A1[MAXN],B1[MAXN];
  12. int qpow(int x,int k){
  13. int mul=x;
  14. int re=1;
  15. for(int i=1;i<=k;i<<=1){
  16. if(i&k)
  17. re=1LL*re*mul%MOD;
  18. mul=1LL*mul*mul%MOD;
  19. }
  20. return re;
  21. }
  22. int inv(int x){
  23. return qpow(x,MOD-2);
  24. }
  25. void DFT(int *A,int *re){
  26. for(int i=0;i<size;i++)
  27. re[i]=A[loc[i]];
  28. for(int k=2;k<=size;k<<=1){
  29. int len=k>>1;
  30. int Wn=qpow(G,(MOD-1)/k);
  31. for(int i=0;i<size;i+=k){
  32. int W=1;
  33. for(int j=0;j<len;j++){
  34. int u=re[i+j],v=1LL*re[i+j+len]*W%MOD;
  35. re[i+j]=(u+v)%MOD;
  36. re[i+j+len]=(u-v+MOD)%MOD;
  37. W=1LL*W*Wn%MOD;
  38. }
  39. }
  40. }
  41. }
  42. void IDFT(int *A,int *re){
  43. for(int i=0;i<size;i++)
  44. re[i]=A[loc[i]];
  45. for(int k=2;k<=size;k<<=1){
  46. int len=k>>1;
  47. int Wn=inv(qpow(G,(MOD-1)/k));
  48. for(int i=0;i<size;i+=k){
  49. int W=1;
  50. for(int j=0;j<len;j++){
  51. int u=re[i+j],v=1LL*re[i+j+len]*W%MOD;
  52. re[i+j]=(u+v)%MOD;
  53. re[i+j+len]=(u-v+MOD)%MOD;
  54. W=1LL*W*Wn%MOD;
  55. }
  56. }
  57. }
  58. int tmp=inv(size);
  59. for(int i=0;i<size;i++)
  60. re[i]=1LL*re[i]*tmp%MOD;
  61. }
  62. public:
  63. Fast_Number_Theory_Transform(int _x,int _y):MOD(_x),G(_y){}
  64. void Transform(int *A,int *B,int *ans){
  65. DFT(A,A1);DFT(B,B1);
  66. for(int i=0;i<size;i++)
  67. A1[i]=1LL*A1[i]*B1[i]%MOD;
  68. IDFT(A1,ans);
  69. }
  70. void init(int x){
  71. for(size=1,bit_length=-1;size<x;size<<=1,bit_length++);
  72. int now=0;
  73. for(int i=0;i<size;i++){
  74. loc[i]=now;
  75. for(int j=1<<bit_length;(now^=j)<j;j>>=1);
  76. }
  77. }
  78. }NTT((479<<21)+1,3);
  79. int a[MAXN],b[MAXN];
  80. int ans[MAXN];
  81. int main(){
  82. int n,m;
  83. scanf("%d%d",&n,&m);
  84. for(int i=0;i<=n;i++)
  85. scanf("%d",&a[i]);
  86. for(int i=0;i<=m;i++)
  87. scanf("%d",&b[i]);
  88. NTT.init(n+m+1);
  89. NTT.Transform(a,b,ans);
  90. for(int i=0;i<=m+n;i++)
  91. printf("%d ",ans[i]);
  92. puts("");
  93. return 0;
  94. }
添加新批注
在作者公开此批注前,只有你和作者可见。
回复批注