【POJ3415】 Common Substrings（后缀数组|SAM）

Common Substrings

Description

A substring of a string T is defined as:
T(i, k)=T_iT_i₊₁...T_i+k_-1, 1≤i≤i+k-1≤|T|.
Given two strings A, B and one integer K, we define S, a set of triples (i, j, k):
S = {(i, j, k) | k≥K, A(i, k)=B(j, k)}.
You are to give the value of |S| for specific A, B and K.

Input

The input file contains several blocks of data. For each block, the first line contains one integer K, followed by two lines containing strings A and B, respectively. The input file is ended by K=0.

1 ≤ |A|, |B| ≤ 10⁵
1 ≤ K ≤ min{|A|, |B|}
Characters of A and B are all Latin letters.

Output

For each case, output an integer |S|.

Sample Input
2
aababaa
abaabaa
1
xx
xx
0
Sample Output
22
5

【题意】

　　给两个长度不超过100000的字符串A和B, 求三元组(i, j, k)的数量, 即满足A串从i开始的后缀与B串从j开始的后缀有长度为k的公共前缀, 还要求k不小于某个给你的数K.

【分析】

　　如果i后缀与j后缀的LCP长度为L, 在L不小于K的情况下, 它对答案的贡献为L - K + 1.

　　所以问题就是快速求出两个串的后缀的LCP。

　　就要用到后缀数组，把两个串拼成一个串，中间加一个特殊字符，然后做后缀数组。

　　求出height数组后根据k分组（min大于等于k的分在一组），同一组的LCP一定大于等于k。（height数组的性质啦，这个很基本）

　　接下来就是求出sum｛L-K+1｝，如果直接暴力的话for2遍加一个RMQ也很慢。

　　然后看到大神们说用什么单调栈，哦~~ 单调的优美性质！！！

　　可以分成两步求，先求B串在下，A串在上的ans。再反过来求一遍。

　　过程好像比较难说清楚。

　　反正...对于B串上面的一堆A串，因为LCP是求min，所以LCP必定是单调递增，所以像下图那样做：

　　【POJ3415】 Common Substrings（后缀数组|SAM）

　　开了三个数组，s表示其单位贡献值，cnt表示有多少个A串是这个贡献值，h为从当前位置到底部的cnt*s的和（h便于计算的）。

　　然后就更新和替换。遇到一个B串就把h[tp]累加到ans中即可。　

　　反过来做也是一样的。

　　单调！！单调！！单调哦！！

　　好像很厉害！！

　　主要看代码~~

　　还有，这道题是有大写字母的！！！坑的我RE了巨久！！

　　还有记得long long！！

代码如下：

  1 #include<cstdio>
  2 #include<cstdlib>
  3 #include<cstring>
  4 #include<iostream>
  5 #include<algorithm>
  6 #include<queue>
  7 using namespace std;
  8 #define INF 0xfffffff
  9 #define Maxl 200010
 10 #define LL long long
 11 
 12 int k,la;
 13 char a[Maxl],b[Maxl];
 14 int c[Maxl];
 15 int cl;
 16 
 17 int sa[Maxl],rank[Maxl],Rs[Maxl],wr[Maxl],y[Maxl];
 18 //sa -> 排名第几的是谁
 19 //rank -> i的排名
 20 //Rs数值小于等于i的有多少个
 21 //y -> 第二关键字排名第几的是谁（类似sa）
 22 int height[Maxl];
 23 
 24 void get_sa(int m)
 25 {
 26    memcpy(rank,c,sizeof(rank));
 27     for(int i=0;i<=m;i++) Rs[i]=0;
 28     for(int i=1;i<=cl;i++) Rs[rank[i]]++;
 29     for(int i=1;i<=m;i++) Rs[i]+=Rs[i-1];
 30     for(int i=cl;i>=1;i--) sa[Rs[rank[i]]--]=i;
 31     
 32     int ln=1,p=0;
 33     while(p<cl)
 34     {
 35        int k=0;
 36         for(int i=cl-ln+1;i<=cl;i++) y[++k]=i;
 37         for(int i=1;i<=cl;i++) if(sa[i]>ln) y[++k]=sa[i]-ln;
 38         for(int i=1;i<=cl;i++) wr[i]=rank[y[i]];
 39         
 40         for(int i=0;i<=m;i++) Rs[i]=0;
 41         for(int i=1;i<=cl;i++) Rs[wr[i]]++;
 42         for(int i=1;i<=m;i++) Rs[i]+=Rs[i-1];
 43         for(int i=cl;i>=1;i--) sa[Rs[wr[i]]--]=y[i];
 44         
 45         for(int i=1;i<=cl;i++) wr[i]=rank[i];
 46         for(int i=cl+1;i<=cl+ln;i++) wr[i]=0;
 47         p=1;rank[sa[1]]=1;
 48         for(int i=2;i<=cl;i++)
 49         {
 50             if(wr[sa[i]]!=wr[sa[i-1]]||wr[sa[i]+ln]!=wr[sa[i-1]+ln]) p++;
 51             rank[sa[i]]=p;
 52         }
 53         m=p,ln*=2;
 54     }
 55     sa[0]=rank[0]=0;
 56 }
 57 
 58 void get_he()
 59 {
 60     int kk=0;
 61     for(int i=1;i<=cl;i++)
 62     {
 63         // if(rank[i]==1) break;
 64         int j=sa[rank[i]-1];
 65         if(kk) kk--;
 66         while(c[i+kk]==c[j+kk]&&i+kk<=cl&&j+kk<=cl) kk++;
 67         height[rank[i]]=kk;
 68     }
 69 }
 70 
 71 LL h[Maxl],cnt[Maxl],s[Maxl];
 72 int tp;
 73 void ffind()
 74 {
 75     LL ans=0;tp=0;
 76     s[0]=INF;
 77     //******
 78     for(int i=1;i<=cl-1;i++)
 79     {
 80         if(sa[i]>=la+1) //B串
 81          ans+=h[tp];
 82         if(height[i+1]-k+1<s[tp])//替换
 83         {
 84             LL sum=0;
 85             while(height[i+1]-k+1<=s[tp]&&tp) sum+=cnt[tp--];
 86             s[++tp]=(height[i+1]-k+1),cnt[tp]=sum,h[tp]=h[tp-1]+cnt[tp]*s[tp];
 87         }
 88         else s[++tp]=height[i+1]-k+1,cnt[tp]=0,h[tp]=h[tp-1];
 89         if(sa[i]<=la) cnt[tp]++,h[tp]+=s[tp];//A串
 90         if(height[i+1]<k)
 91         {
 92             tp=0;s[0]=INF;
 93         }
 94     }tp=0;s[tp]=INF;
 95     for(int i=1;i<=cl-1;i++)
 96     {
 97         if(sa[i]<=la) //A串
 98          ans+=h[tp];
 99         if(height[i+1]-k+1<s[tp])//替换
100         {
101             LL sum=0;
102             while(height[i+1]-k+1<=s[tp]&&tp) sum+=cnt[tp--];
103             s[++tp]=(height[i+1]-k+1),cnt[tp]=sum,h[tp]=h[tp-1]+cnt[tp]*s[tp];
104         }
105         else s[++tp]=height[i+1]-k+1,cnt[tp]=0,h[tp]=h[tp-1];
106         if(sa[i]>=la+1) cnt[tp]++,h[tp]+=s[tp];//B串
107         if(height[i+1]<k)
108         {
109             tp=0;s[tp]=INF;
110         }
111     }
112     printf("%I64d\n",ans);
113 }
114 
115 void init()
116 {
117     scanf("%s%s",a,b);
118     int l=strlen(a);cl=0;
119     la=l;
120     for(int i=0;i<l;i++)
121     {
122         if(a[i]>='a'&&a[i]<='z') c[++cl]=a[i]-'a'+1;
123         else c[++cl]=a[i]-'A'+27;
124     }
125     c[++cl]=55;
126     l=strlen(b);
127     for(int i=0;i<l;i++)
128     {
129         if(b[i]>='a'&&b[i]<='z') c[++cl]=b[i]-'a'+1;
130         else c[++cl]=b[i]-'A'+27;
131     }
132 }
133 
134 int main()
135 {
136     while(1)
137     {
138         scanf("%d",&k);
139         if(k==0) break;
140         init();
141         get_sa(55);
142         get_he();
143         ffind();
144     }
145     return 0;
146 }

[POJ3415]