后缀数组+莫比乌斯函数
1 #include <stdio.h> 2 #include <string.h> 3 #include<algorithm> 4 using namespace std; 5 const int N = 401000; 6 7 const int MAXN = N; //数列的最大值和个数上界 8 9 struct SuffixArray{ 10 int wa[MAXN]; //用来进行基数排序或临时变量 11 int wb[MAXN]; //用来进行基数排序或临时变量 12 int wv[MAXN]; //用来进行基数排序或临时变量 13 int ws[MAXN]; //用来进行基数排序或临时变量 14 15 int sa[MAXN]; //sa[i]代表排名为i的后缀在原数列起始下标(数列的下标从0开始)，sa[0]肯定等于n，因为标兵为最小的。 16 int rank[MAXN]; //rank[i]代表suffix[i]的排名，rank[n]肯定等于0，理由同上。 17 int height[MAXN]; //height[i]代表排名为i - 1的后缀 和排名为i的后缀 的最长公共连续子序列 的长度。 18 int r[MAXN]; //r[]存放原数列下标从0到n，r[n]为a标兵，是r[]里面最小的. 19 int n; //数列的元素个数，不包括标兵 20 int m; //存放最大值，r[]数组的数都要小于m，用来进行基数排序 21 22 void input(int *val, int len, int Max){//Max要大于r[0..len - 1]，因为内部采用了基数排序 23 for (int i = 0;i < len;i++) 24 r[i] = val[i]; 25 r[len] = 0; //最小值，起标兵作用 26 n = len; 27 m = Max; 28 calSa(); 29 calHeight(); 30 } 31 32 int cmp(int *r, int a, int b, int l){ 33 return (r[a] == r[b] && r[a + l] == r[b + l]); 34 } 35 36 void calSa(){ //求sa数组 37 int i, j, p, *x = wa, *y = wb, *t; 38 for (i = 0;i < m;i++) ws[i] = 0; 39 for (i = 0;i < n + 1;i++) ws[x[i] = r[i]]++; 40 for (i = 1;i < m;i++) ws[i] += ws[i - 1]; 41 for (i = n;i >= 0;i--) sa[--ws[x[i]]] = i; 42 for (j = 1, p = 1;p < n + 1;j *= 2, m = p){ 43 for (p = 0, i = n - j + 1;i < n + 1;i++) y[p++] = i; 44 for (i = 0;i < n + 1;i++) if (sa[i] >= j) y[p++] = sa[i] - j; 45 for (i = 0;i < n + 1;i++) wv[i] = x[y[i]]; 46 for (i = 0;i < m;i++) ws[i] = 0; 47 for (i = 0;i < n + 1;i++) ws[wv[i]]++; 48 for (i = 1;i < m;i++) ws[i] += ws[i - 1]; 49 for (i = n;i >= 0;i--) sa[--ws[wv[i]]] = y[i]; 50 for (t = x, x = y, y = t, p = 1, x[sa[0]] = 0, i = 1; i < n + 1;i++) 51 x[sa[i]] = cmp(y, sa[i - 1], sa[i], j) ? p - 1 : p++; 52 } 53 } 54 55 void calHeight(){ //求rank和height数组 56 int i, j, k = 0; 57 for (i = 1;i <= n;i++) rank[sa[i]] = i; 58 for (i = 0;i < n;height[rank[i++]] = k) 59 for (k?k--:0, j = sa[rank[i]- 1];r[i + k] == r[j + k];k++); 60 } 61 //-求任意两后缀最长公共前缀问题，如果没用到可以去掉---------------------- 62 int Log[MAXN]; 63 int best[20][MAXN]; 64 void initRMQ() { //初始化RMQ 标准RMQ 预处理O(nlgn) 65 Log[0] = -1; 66 for(int i = 1;i <= MAXN;i++){ 67 Log[i]=(i & (i - 1))?Log[i - 1] : Log[i - 1] + 1 ; 68 } 69 for(int i = 1; i <= n ; i ++) best[0][i] = height[i]; 70 for(int i = 1; i <= Log[n] ; i ++) { 71 int limit = n - (1<<i) + 1; 72 for(int j = 1; j <= limit ; j ++) { 73 best[i][j] = (best[i-1][j] > best[i-1][j+(1<<i>>1)]) ? best[i-1][j+(1<<i>>1)] : best[i-1][j]; 74 } 75 } 76 } 77 int lcp(int a,int b) { //询问suffix[a]于suffix[b]的最长公共前缀(标准RMQ 询问O(1)) 使用这个函数之前要先后调用initRMQ()和input() 78 a = rank[a]; b = rank[b]; 79 if(a > b){ 80 int t = a; 81 a = b; 82 b = t; 83 } 84 a ++; 85 int t = Log[b - a + 1]; 86 return (best[t][a] > best[t][b - (1<<t) + 1])? best[t][b - (1<<t) + 1] : best[t][a]; 87 } 88 }SA1,SA2; 89 char s[N]; 90 int len,a[N],b[N],i,j,k,L,R,da,db; 91 int prime[N],mu[N],cnt,vis[N]; 92 long long ans; 93 int main() 94 { 95 mu[1]=1; 96 for (i=2;i<=150000;i++) 97 { 98 if (!vis[i]) 99 { 100 prime[cnt++]=i; 101 mu[i]=-1; 102 } 103 for (j=0;j<cnt&&i*prime[j]<150000;j++) 104 { 105 vis[i*prime[j]]=1; 106 if (i%prime[j]) mu[i*prime[j]]=-mu[i]; 107 else {mu[i*prime[j]]=0;break;} 108 } 109 } 110 int test; 111 scanf("%d",&test); 112 while (test--) 113 { 114 scanf("%s",&s); 115 len=strlen(s); 116 ans=0; 117 for (i=1;i<=len;i++) ans+=i; 118 for (i=0;i<len;i++) 119 { 120 a[i]=s[i]; 121 b[len-i-1]=s[i]; 122 } 123 SA1.input(a,len,200);SA1.initRMQ(); 124 SA2.input(b,len,200);SA2.initRMQ(); 125 for (i=1;i<=len;i++) 126 { 127 j=0; 128 while (j*i<len) 129 { 130 for (k=j+1;k*i<len;k++) 131 if (SA1.lcp(j*i,k*i)<i) break; 132 L=j*i;R=min(len-1,k*i-1); 133 if (R+1!=len) da=min(i,SA1.lcp(L,R+1));else da=0; 134 if (L!=0) db=min(i,SA2.lcp(len-1-R,len-L));else db=0; 135 L=L-db;R=R+da; 136 j=k; 137 for (k=2;k*i<=R-L+1;k++) 138 ans+=mu[k]*(R-L+1-k*i+1); 139 if (R==len-1) break; 140 } 141 //printf(" "); 142 } 143 printf("%I64d ",ans); 144 } 145 }