POJ 3693 Maximum repetition substring（后缀数组+RMQ）

Maximum repetition substring

Time Limit: 1000MS	Memory Limit: 65536K
Total Submissions: 11322	Accepted: 3487

Description
The repetitionnumber of a string is defined as the maximum numberR such that thestring can be partitioned into
R same consecutive substrings. Forexample, the repetition number of "ababab" is 3 and "ababa"is 1.
Given a stringcontaining lowercase letters, you are to find a substring of it with maximumrepetition number.
Input
The input consistsof multiple test cases. Each test case contains exactly one line, which

gives a non-empty string consisting of lowercase letters. The length of thestring will not be greater than 100,000.
The last test caseis followed by a line containing a '#'.
Output
For each testcase, print a line containing the test case number( beginning with 1) followedby the substring of maximum repetition number. If there are multiple substringsof maximum repetition number, print the lexicographically
smallest one.
Sample Input
ccabababc
daabbccaa
#
Sample Output
Case 1: ababab
Case 2: aa
Source
2008 Asia Hefei Regional Contest Online by USTC

        此题也是一道非常之巧妙的后缀数组题目。
        大致题意，就是求一个字符串的最多连续出现的子串，即一个子串，可以分成最多个相同的循环节的组合。
        对于本题，我一开始有一个想法，根据之前做题累积的经验套路，我们可以枚举重复次数，然后再枚举一个开始节点，求区间rmq，如果大于0，那么就是一个可能的解。但是显然这样子有重大的问题。首先，忽视了这个“连续”，求区间rmq可以求出出现至少k次的最大长度，但是这个出现k次并不一定是连续的。举个例子说，ababac，有后缀ababac、abac、ac，这连续区间内，height的最小值为1也即a出现了3次，但是这三个a出现的缺失间隔出现不是连续出现的。其次，这个算法是O(N^2)的，在时间上也行不通。

       所以得考虑别的方法。看了著名论文的解释和网上许多dalao的博客之后，才算是彻底理解了。首先，就是要枚举循环节的长度，然后根据观察，如果确实存在这样长度的循环节，那么在s[0]、s[i]、s[2*i]、……、s[j*i]中，至少有一对相邻的位置的字符（注意是单个字符）相同，具体证明不好给出。接下来的步骤，举例子来说明。例如：1234234234，当长度l枚举到3的时候，s[3]==s[6]=='4'，我们再计算lcp(3,6)，这个根据定理是区间 ( Rank[3] , Rank[6] ]
的height的最小值。用ST表计算rmq可以得到lcp，这里lcp==4。显然从3这个位置开始，重复次数是lcp/l +1 == 2，也即对应后面两个'234234'。但是很不幸，'234'这个循环节实际上出现了3次，因为我们初始位置3并不一定是循环节的首位，所以可能会少算一些。为了解决这个问题，我们可以相应往前挪，即补全这段长度。起始位置L=3-(l-lcp%l)，这个好理解，我们在看看 
( Rank[L] , Rank[L+l] ]  的lcp是否大于等于之前求的lcp，相等则说明这前面也是一个循环节，重复次数加一。我们求出重复次数最多的循环节长度，作为备选解存入数组中。

        因为会有多解，而且要求输出字典序最小，所以我们继续利用后缀数组。这次是利用sa[]，由于sa[]就是按照字典序排的，所以我们对于数组中保存的所有备选循环节长度进行尝试，找到第一个符合条件的解输出即可。具体见代码：

#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
#define LL long long
#define N 100010

using namespace std;

struct Suffix_Array
{
int sa
,Rank
,h
,n;
int xx
,yy
,c
; char *s;
bool cmp(int *s,int x,int y,int k)
{return (s[x]==s[y])&&(s[x+k]==s[y+k]);}
void ins(char *str) {s=str;n=strlen(s)+1;}

void DA()
{
memset(c,0,sizeof(c));
int *x=xx,*y=yy,m=130,*t,i;
for(i=0;i<n;i++) x[i]=s[i];
for(i=0;i<n;i++) c[x[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(int k=1,tot=0;tot<n;k<<=1,m=tot)
{
memset(c,0,sizeof(c));
for(i=0;i<n;i++) c[x[i]]++;
for(i=1;i<m;i++) c[i]+=c[i-1];
for(i=n-k,tot=0;i<n;i++) y[tot++]=i;
for(i=0;i<n;i++) if (sa[i]>=k) y[tot++]=sa[i]-k;
for(i=n-1;i>=0;i--) sa[--c[x[y[i]]]]=y[i];
for(i=tot=1,t=x,x=y,y=t,x[sa[0]]=0;i<n;i++)
x[sa[i]]=cmp(y,sa[i-1],sa[i],k)?tot-1:tot++;
}
}

void cal_height()
{
int i,j,k=0;
for(i=1;i<n;i++) Rank[sa[i]]=i;
for(i=0;i<n-1;h[Rank[i++]]=k)
for(k?k--:0,j=sa[Rank[i]-1];s[i+k]==s[j+k];k++);
}

} SA;

struct ST
{
int dp
[18];

void init(int n)
{
int m=floor(log(n+0.0)/log(2.0));
for(int i=1;i<=n;i++) dp[i][0]=SA.h[i];
for(int i=1;i<=m;i++)
for(int j=n;j;j--)
{
dp[j][i]=dp[j][i-1];
if(j+(1<<(i-1))<=n) dp[j][i]=min(dp[j][i],dp[j+(1<<(i-1))][i-1]);
}
}

int rmq(int l,int r)
{
int m=floor(log(r-l+1.0)/log(2.0));
return min(dp[l][m],dp[r-(1<<m)+1][m]);
}

} ST;

int n,k,ans
;
char s
;

int lcp(int l,int r)
{
l=SA.Rank[l];
r=SA.Rank[r];
if (l>r) swap(l,r);
return  ST.rmq(l+1,r);
}

int main()
{
int T_T=0;
while(~scanf("%s",s))
{
int maxlen=0,tot;
if (s[0]=='#') break;
SA.ins(s);SA.DA();
SA.cal_height();
ST.init(n=strlen(s));
for(int l=1;l<=n/2;l++)
{
for(int i=1;i+l<=n;i+=l)
{
int LCP=lcp(i,i+l);
int tmp=LCP/l+1,res=LCP%l,L;
if (res&&(L=i-l+res)>0&&lcp(L,L+l)>=LCP) tmp++;
if (tmp>maxlen) maxlen=tmp,ans[tot=1]=l;
else if (tmp==maxlen) ans[++tot]=l;
}
}
int flag=1,b,e;
for(int i=1;i<=n&&flag;i++)
{
for(int j=1;j<=tot;j++)
{
int l=ans[j],sa=SA.sa[i];
if (sa+l<=n&&lcp(sa,sa+l)>=(maxlen-1)*l) b=sa,e=b+maxlen*l,flag=0;
}
}
printf("Case %d: ",++T_T);
for(int i=b;i<e;i++)
printf("%c",s[i]);
puts("");
}
return 0;
}