POJ 3415 Common Substrings（后缀数组+单调栈）

Description
给出两个串a和b，求a和b的长度不小于k的公共子串数量
Input
多组用例，每组用例第一行为一整数k，之后为两个字符串表示a和b，两个串的串长不超过10^5，k不超过两串串长，以k=0结束输入
Output
对于每组用例，输出a和b的长度不小于k的公共子串数量
Sample Input
2
aababaa
abaabaa
1
xx
xx
0
Sample Output
22
5
Solution
将两串连在一起（中间用隔离字符隔离以避免匹配越界），求出height数组后，按k分组，问题转化为求每一组的sum(lcp(i,j)-k+1)（第i个后缀属于第一个串，第j个后缀属于第二个串），如果常规去固定一个后缀去枚举另一个后缀，虽然rmq预处理可以使得每次操作是O(1)的，但是总复杂度也是O(n^2)的，所以要优化。考虑第i个后缀，假设其属于第一个串，那么这个值对答案的贡献就是sum(min(height[j]~height[i])-k+1)（j=2,3,…,i-1且第j个后缀属于第二个串），而这个最小值数组是非减的，那么考虑构造一个单调递增的栈，将height-k+1这一列数按顺序入栈，如果当前height[i]-k+1值大于等于栈顶元素则将height[i]-k+1入栈，否则就一直更新栈上部元素，将其值变为height[i]-k+1并累加height[i]-k+1-sta[p]（小于0）的值进入num[id[p]]中，id[p]记录第p个后缀属于哪一个串，这样做的意义就是在第i个后缀进栈时，假设在其之前进栈的元素与其的最长公共前缀都是height[i]-k+1，而多计算的部分用num数组来记录并累加到答案中，这样就可以用接近O(n)的复杂度求出答案，具体实现见代码，注意答案用long long记录
Code

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
using namespace std;
#define maxn 222222
typedef long long ll;
int t1[maxn],t2[maxn],c[maxn],sa[maxn],rank[maxn],height[maxn];
bool cmp(int *r,int a,int b,int l)
{
    return r[a]==r[b]&&r[a+l]==r[b+l];
}
void da(int str[],int n,int m)
{
    n++;
    int i,j,p,*x=t1,*y=t2;
    for(i=0;i<m;i++)c[i]=0;
    for(i=0;i<n;i++)c[x[i]=str[i]]++;
    for(i=1;i<m;i++)c[i]+=c[i-1];
    for(i=n-1;i>=0;i--)sa[--c[x[i]]]=i;
    for(j=1;j<=n;j<<=1)
    {
        p=0;
        for(i=n-j;i<n;i++)y[p++]=i;
        for(i=0;i<n;i++)if(sa[i]>=j)y[p++]=sa[i]-j;
        for(i=0;i<m;i++)c[i]=0;
        for(i=0;i<n;i++)c[x[y[i]]]++;
        for(i=1;i<m;i++)c[i]+=c[i-1];
        for(i=n-1;i>=0;i--)sa[--c[x[y[i]]]]=y[i];
        swap(x,y);
        p=1;x[sa[0]]=0;
        for(i=1;i<n;i++)
        x[sa[i]]=cmp(y,sa[i-1],sa[i],j)?p-1:p++;
        if(p>=n)break;
        m=p;
    }
    int k=0;
    n--;
    for(i=0;i<=n;i++)rank[sa[i]]=i;
        for(i=0;i<n;i++)
        {
            if(k)k--;
            j=sa[rank[i]-1];
            while(str[i+k]==str[j+k])k++;
            height[rank[i]]=k;
        }
}
int k,n,a[maxn],id[maxn],sta[maxn],p;
char s[maxn];
int main()
{
    while(~scanf("%d",&k),k)
    {
        scanf("%s",s);
        int len=strlen(s);
        s[len]='#';
        scanf("%s",s+len+1);
        n=strlen(s);
        for(int i=0;i<n;i++)a[i]=s[i];
        a[n]=0;
        da(a,n,200);
        ll ans=0,num[3]={0};
        memset(id,0,sizeof(id));
        for(int i=2;i<=n;i++)
        {
            if(height[i]<k)
            {
                p=num[1]=num[2]=0;
                continue;
            }
            for(int pp=p;pp&&sta[pp]>height[i]-k+1;pp--)
            {
                num[id[pp]]+=height[i]-k+1-sta[pp];
                sta[pp]=height[i]-k+1;
            }
            sta[++p]=height[i]-k+1;
            if(sa[i-1]<len)id[p]=1;
            else if(sa[i-1]>len)id[p]=2;
            num[id[p]]+=height[i]-k+1;
            if(sa[i]<len)ans+=num[2];
            else if(sa[i]>len)ans+=num[1];
        }
        printf("%lld\n",ans);
    }
}