BZOJ4259 残缺的字符串 FFT
首先是csy博客的题解,但csy的并不十分优秀。
事实上,我们可以求出f[i] = sigma(A[j] * inv[B[i - j]]),当f[i] <= n时,匹配,否则不匹配。这样就只需要做一次FFT或者NTT即可。
1 #include<cstdio> 2 #include<cmath> 3 #include<algorithm> 4 #define N 1048576 5 using namespace std; 6 char sa[N],sb[N];int n,m,i,j,k,a[N],b[N],ans,q[N]; 7 struct comp{ 8 double r,i;comp(double _r=0,double _i=0){r=_r;i=_i;} 9 comp operator+(const comp&x){return comp(r+x.r,i+x.i);} 10 comp operator-(const comp&x){return comp(r-x.r,i-x.i);} 11 comp operator*(const comp&x){return comp(r*x.r-i*x.i,r*x.i+i*x.r);} 12 }A[N],B[N],C[N]; 13 const double pi=acos(-1.0); 14 void FFT(comp*a,int n,int t){ 15 for(int i=1,j=0;i<n-1;i++){ 16 for(int s=n;j^=s>>=1,~j&s;); 17 if(i<j)swap(a[i],a[j]); 18 } 19 for(int d=0;(1<<d)<n;d++){ 20 int m=1<<d,m2=m<<1; 21 double o=pi/m*t;comp _w(cos(o),sin(o)); 22 for(int i=0;i<n;i+=m2){ 23 comp w(1,0); 24 for(int j=0;j<m;j++){ 25 comp &A=a[i+j+m],&B=a[i+j],t=w*A; 26 A=B-t;B=B+t;w=w*_w; 27 } 28 } 29 } 30 if(t==-1)for(int i=0;i<n;i++)a[i].r/=n; 31 } 32 int main(){ 33 scanf("%d%d%s%s",&m,&n,sa,sb); 34 for(i=0,j=m-1;i<j;i++,j--)swap(sa[i],sa[j]); 35 for(i=0;i<m;i++)if(sa[i]!='*')a[i]=sa[i]-'a'+1; 36 for(i=0;i<n;i++)if(sb[i]!='*')b[i]=sb[i]-'a'+1; 37 for(k=1;k<n+m;k<<=1); 38 for(i=0;i<k;i++)A[i]=comp(a[i]*a[i]*a[i],0),B[i]=comp(b[i],0); 39 for(FFT(A,k,1),FFT(B,k,1),i=0;i<k;i++)C[i]=C[i]+A[i]*B[i]; 40 for(i=0;i<k;i++)A[i]=comp(a[i],0),B[i]=comp(b[i]*b[i]*b[i],0); 41 for(FFT(A,k,1),FFT(B,k,1),i=0;i<k;i++)C[i]=C[i]+A[i]*B[i]; 42 for(i=0;i<k;i++)A[i]=comp(a[i]*a[i],0),B[i]=comp(b[i]*b[i],0); 43 for(FFT(A,k,1),FFT(B,k,1),i=0;i<k;i++)C[i]=C[i]-A[i]*B[i]*comp(2,0); 44 FFT(C,k,-1); 45 for(i=m-1;i<n;i++)if(C[i].r<0.5)q[++ans]=i-m+2; 46 for(printf("%d\n",ans),i=1;i<ans;i++)printf("%d ",q[i]); 47 if(ans)printf("%d",q[ans]); 48 return 0; 49 }
NTT版本代码
1 #include<bits/stdc++.h> 2 using namespace std; 3 4 const int MOD = 998244353; 5 const int N = 1048576; 6 const int M = 998244353; 7 8 int WM[N + 2], IWM[N + 2]; 9 vector<int> IW[N + 2], W[N + 2]; 10 11 int Pw(int a, int b, int p) 12 { 13 int v = 1; 14 for(; b; b >>= 1, a = 1ll * a * a % p) 15 if (b & 1) 16 v = 1ll * v * a % p; 17 return v; 18 } 19 20 void ycl() 21 { 22 for (int m = 2; m <= N; m <<= 1) 23 { 24 WM[m] = Pw(3, (M - 1) / m, M), IWM[m] = Pw(3, (M - 1) / m * (m - 1), M); 25 int o = 1; 26 W[m].push_back(o); 27 for (int i = 1; i < m; ++i) 28 o = 1ll * o * WM[m] % M, W[m].push_back(o); 29 o = 1; 30 IW[m].push_back(o); 31 for (int i = 1; i < m; ++i) 32 o = 1ll * o * IWM[m] % M, IW[m].push_back(o); 33 } 34 } 35 36 void NTT(int *a, int n, int f = 1) 37 { 38 int i, j, k, m, w, u, v; 39 for (i = n >> 1, j = 1; j < n - 1; ++j) 40 { 41 if (i > j) 42 swap(a[i], a[j]); 43 for (k = n >> 1; k <= i; k >>= 1) 44 i ^= k; 45 i ^= k; 46 } 47 for (m = 2; m <= n; m <<= 1) 48 for (i = 0; i < n; i += m) 49 for (j = i; j < i + (m >> 1); ++j) 50 if (a[j] || a[j + (m >> 1)]) 51 { 52 u = a[j]; 53 v = 1ll * (f == 1 ? W[m][j - i] : IW[m][j - i]) * a[j + (m >> 1)] % M; 54 if ((a[j] = u + v) >= M) 55 a[j] -= M; 56 if ((a[j + (m >> 1)] = u - v) < 0) 57 a[j + (m >> 1)] += M; 58 } 59 if (f == -1) 60 for (w = Pw(n, M - 2, M), i = 0; i < n; ++i) 61 a[i] = 1ll * a[i] * w % M; 62 } 63 64 char sa[N], sb[N], st[N]; 65 int a[5][N], b[5][N], c[N], q[N]; 66 67 int main() 68 { 69 ycl(); 70 int m, n; 71 scanf("%d%d", &m, &n); 72 scanf("%s%s", sa, sb); 73 for (int i = 0, j = m - 1; i < j; ++i, --j) 74 swap(sa[i],sa[j]); 75 for (int i = 0; i < m; ++i) 76 a[1][i] = sa[i] - 'a' + 1; 77 for (int i = 0; i < n; ++i) 78 b[1][i] = sb[i] - 'a' + 1; 79 for (int i = 0; i < m; ++i) 80 if(sa[i] == '*') 81 a[1][i] = 0; 82 for (int i = 0; i < n; ++i) 83 if(sb[i] == '*') 84 b[1][i] = 0; 85 int K; 86 for (K = 1; K < n + m; K <<= 1); 87 for (int i = 2; i <= 3; ++i) 88 { 89 for (int j = 0; j < m; ++j) 90 a[i][j] = a[i - 1][j] * a[1][j]; 91 for (int j = 0; j < n; ++j) 92 b[i][j] = b[i - 1][j] * b[1][j]; 93 } 94 for (int i = 1; i < 4; ++ i) 95 { 96 NTT(a[i], K); 97 NTT(b[i], K); 98 } 99 for (int i = 0; i < K; ++ i) 100 { 101 c[i] = (c[i] + 1LL * a[3][i] * b[1][i] % MOD) % MOD; 102 c[i] = (c[i] + MOD - 2LL * a[2][i] * b[2][i] % MOD) % MOD; 103 c[i] = (c[i] + 1LL * a[1][i] * b[3][i] % MOD) % MOD; 104 } 105 NTT(c, K, -1); 106 107 int ans(0); 108 for (int i = m - 1; i < n; ++i) 109 if (!c[i]) 110 q[++ans] = i - m + 2; 111 printf("%d\n", ans); 112 for (int i = 1; i < ans; ++i) 113 printf("%d ", q[i]); 114 if (ans) 115 printf("%d\n", q[ans]); 116 117 return 0; 118 }