• 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
    ",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
    ", ans);
    112     for (int i = 1; i < ans; ++i)
    113         printf("%d ", q[i]);
    114     if (ans)
    115         printf("%d
    ", q[ans]);
    116 
    117     return 0;
    118 }
  • 相关阅读:
    EMF介绍系列(一、EMF与MDA)
    EMF介绍系列(四、枚举类型、自定义类型和Map)
    使用osgi.util.NLS简化资源文件访问
    2012 定制化产品探讨(周金根).pdf
    敏捷个人理念与模型PPT及今年唯一一次的公开线上课堂
    生活:父与子三亚行
    与北邮学子交流成长,敏捷个人总体介绍 PPT
    敏捷个人教你如何制作2012生活看板
    敏捷个人架构图 V1.3
    敏捷个人微刊封面及敏捷个人使命和加入社区方式
  • 原文地址:https://www.cnblogs.com/aseer/p/8968537.html
Copyright © 2020-2023  润新知