「笔记」后缀数组

目录

写在前面

投了日报过了。
@ComeIntoPower 审日报的时候我正好在线，其实这篇一开始是没有过的。
草
心死了，犇犇都发出去了，伤心了 10min 又来看了一眼发现过了（（（
于是顺便来扩充了一下内容。

网上部分题解直接对着优化后面目全非的代码开讲。
*这太野蛮了*
这里主要参考了 OI-wiki 上的讲解。

一些约定

  1. (left| sum ight|)：字符集大小。
  2. (S[i:j])：由字符串 (S) 中 (S_isim S_j) 构成的子串。
  3. (S_1<S_2)：字符串 (S_1) 的字典序 (<S_2)。
  4. 后缀：从某个位置 (i) 开始，到串末尾结束的子串，后缀 (i) 等价于子串 (S[i:n])。每一个后缀都与一个 (1sim n) 的整数一一映射。

SA

SA 的定义

字符串 (S) 的后缀数组定义为两个数组 (sa,rk)。
(sa) 储存 (S) 的所有后缀 按字典序排序后的起始下标，满足 (S[sa_{i-1}:n]<S[sa_i:n])。
(rk) 储存 (S) 的所有后缀的排名。

举例：这里有一个可爱的字符串 (S= exttt{yuyuko})。
有 ( exttt{k}< exttt{o}< exttt{u}< exttt{y})，则它的后缀数组 (sa = [5,6,4,2,3,1])，(rk = [6,4,5,3,1,2])。具体地，有：

排名 1 2 3 4 5 6
下标 (5) (6) (4) (2) (3) (1)
后缀 ( exttt{ko}) ( exttt{o}) ( exttt{uko}) ( exttt{uyuko}) ( exttt{yuko}) ( exttt{yuyuko})

不同后缀的排名必然不同（长度不等），(rk) 中不会有重复值出现。

倍增法构造

考虑字符串大小的从前向后的比较过程，可以先将所有长度为 (2^k) 的子串进行排序，通过相邻子串合并并比较大小，求出所有长度为 (2^{k+1}) 的子串的大小关系。

对于 (S[i:i+2^k-1]) 和 (S[j:j+2^k-1])，分别将它们裂开，分成两长度为 (2^{k-1}) 的串。设 (A_i = S[i:i+2^{k-1}-1])，(B_i = S[i+2^{k-1}:i+2^k-1])。
考虑字典序排序的过程，则 (S[i:i+2^k-1] <S[j:j+2^k-1]) 的条件为：

[[A_i<A_j] operatorname{or} [A_i=A_j operatorname{and} B_i<B_j] ]

考虑每一次倍增时，都使用 sort 按双关键字 (A_i) 和 (B_i) 进行排序，每次倍增都进行依次排序，时间复杂度为 (O(nlog^2 n))。

P3809 【模板】后缀排序
以下是一份简单易懂的代码。

这里定义了两个数组：
(sa_i)：倍增中 排名为 (i) 的长度为 (2^{k-1}) 的子串。
(rk_i)：倍增过程中子串 (S[i:i+2^k-1]) 的排名，
显然它们互为反函数，(sa_{rk_i}=rk_{sa_i} = i)。

初始化 (rk_i = S_i)，即 (S_i) 的 ( ext{ASCII}) 值。
虽然这样不满足值域在 ([1,n]) 内，但体现了大小关系，可用于更新。 (rk) 的值之后还会更新。

//知识点：SA
/*
By:Luckyblock
*/
#include <algorithm>
#include <cctype>
#include <cstdio>
#include <cstring>
#define LL long long
const int kN = 1e6 + 10;
//=============================================================
char s[kN];
int n, m, w, sa[kN], rk[kN << 1], oldrk[kN << 1];
//rk[i]: 倍增过程中子串[i:i+2^k-1]的排名，
//sa[i] 排名为i的子串 [i:i+2^k-1]，
//它们互为反函数。
//存在越界风险，如果不写特判，rk 和 oldrk 要开 2 倍空间。
//=============================================================
inline int read() {
  int f = 1, w = 0;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
void Chkmax(int &fir_, int sec_) {
  if (sec_ > fir_) fir_ = sec_;
}
void Chkmin(int &fir_, int sec_) {
  if (sec_ < fir_) fir_ = sec_;
}
bool CMP(int fir_, int sec_) {
  if (rk[fir_] != rk[sec_]) return rk[fir_] < rk[sec_];
  return rk[fir_ + w] < rk[sec_ + w];
}
int main() {
  scanf("%s", s + 1);
  n = strlen(s + 1);
  m = std::max(n, 300);
  //初始化 rk 和 sa。
  //观察下面的代码可知，每次倍增时都会根据 rk 来更新 sa，则仅须保证 sa 数组是一个 1~n 的排列即可。
  for (int i = 1; i <= n; ++ i) sa[i] = i, rk[i] = s[i];
  //倍增过程。w 是已经推出的子串长度，数值上等于上述分析中的 2^{k-1}。
  //注意此处的 sa 数组存的并不是后缀的排名，存的是倍增过程中指定长度子串的排名。
  for (w = 1; w < n; w <<= 1) {
    std::sort(sa + 1, sa + n + 1, CMP);
    for (int i = 1; i <= n; ++ i) oldrk[i] = rk[i];
    for (int i = 1, p = 0; i <= n; ++ i) {
      if (oldrk[sa[i]] == oldrk[sa[i - 1]] && //判断两个子串是否相等。
          oldrk[sa[i] + w] == oldrk[sa[i - 1] + w]) { //越界风险，2倍空间
        rk[sa[i]] = p;
      } else {
        rk[sa[i]] = ++p;
      }
    }
  }
  for (int i = 1; i <= n; ++ i) printf("%d ", sa[i]);
  return 0;
}

计数排序 与 基数排序

优化上述算法的前置知识。

可以参考：OI-wiki 计数排序 and OI-wiki 基数排序。

计数排序是一种与桶排序类似的排序方法。
将长度为 (n) 的数列 (a) 排序后放入 (b) 的代码如下， 其中 (w) 为值域，即 (max{a_i})。

int a[kMaxn], b[kMaxn], cnt[kMaxw];
for (int i = 1; i <= n; ++ i) ++ cnt[a[i]];
for (int i = 1; i <= w; ++ i) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; -- i) b[cnt[a[i]] --] = a[i];

其中，在对 (cnt) 求前缀和后， (cnt_i) 为不大于 (i) 的数的数量，即为 (i) 的排名。
因此在下一步中，可以根据排名赋值。
复杂度为 (O(n+w))，值域与 (n) 同阶时复杂度比较优秀。

个人认为基数排序只是一种思想，并不算一种独立的排序方法。
它仅仅是将 (k) 个排序关键字分开，按优先级升序依次考虑，从而实现多比较字的排序。实际每次排序还是靠其他排序方法实现。常常与计数排序相结合。

优化

请确保完全理解上述朴素实现后再阅读下文。

发现后缀数组值域即为 (n)，又是多关键字排序，考虑基数排序。
上面已经给出一个用于比较的式子：

[[A_i<A_j] operatorname{or} [A_i=A_j operatorname{and} B_i<B_j] ]

倍增过程中 (A_i,B_i) 大小关系已知，先将 (B_i) 作为第二关键字排序，再将 (A_i) 作为第一关键字排序，两次计数排序实现即可。
单次计数排序复杂度 (O(n + w))（(w) 为值域，显然最大与 (n) 同阶），总时间复杂度变为 (O(nlog n))。

实现时将所有排序替换为基数排序即可。注意初始化。

//知识点：SA
/*
By:Luckyblock
I love Marisa;
But Marisa has died;
*/
#include <cstdio>
#include <ctype.h>
#include <cstring>
#include <algorithm>
#define ll long long
const int kMaxn = 1e6 + 10;
//=============================================================
char S[kMaxn];
//rk[i]: 倍增过程中子串[i:i+2^k-1]的排名，
//sa[i] 排名为i的子串 [i:i+2^k-1]，
//它们互为反函数。
//存在越界风险，如果不写特判，rk 和 oldrk 要开 2 倍空间。
int n, m, sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1];
int id[kMaxn], cnt[kMaxn]; //用于计数排序的两个 temp 数组
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
//=============================================================
int main() {
  scanf("%s", S + 1);
  n = strlen(S + 1);
  m = std :: max(n, 300); //值域大小
  
  //初始化 rk 和 sa
  for (int i = 1; i <= n; ++ i) ++ cnt[rk[i] = S[i]];
  for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
  for (int i = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;

  //倍增过程。w 是已经推出的子串长度，数值上等于上述分析中的 2^{k-1}。
  //注意此处的 sa 数组存的并不是后缀的排名，存的是倍增过程中指定长度子串的排名。
  for (int w = 1; w < n; w <<= 1) {
    //按照后半截 rk[i+w] 作为第二关键字排序。
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) id[i] = i;
    for (int i = 1; i <= n; ++ i) ++ cnt[rk[id[i] + w]]; //越界风险，2倍空间
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rk[id[i] + w]] --] = id[i];

    //按照前半截 rk[i] 作为第一关键字排序。
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) id[i] = sa[i];
    for (int i = 1; i <= n; ++ i) ++ cnt[rk[id[i]]];
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rk[id[i]]] --] = id[i];

    //更新 rk 数组，可以滚动数组一下，但是可读性会比较差(
    for (int i = 1; i <= n; ++ i) oldrk[i] = rk[i];
    for (int p = 0, i = 1; i <= n; ++ i) {
      if (oldrk[sa[i]] == oldrk[sa[i - 1]] &&  //判断两个子串是否相等。
          oldrk[sa[i] + w] == oldrk[sa[i - 1] + w]) { //越界风险，2倍空间
        rk[sa[i]] = p;
      } else {
        rk[sa[i]] = ++ p;
      }
    }
  }
  for (int i = 1; i <= n; ++ i) printf("%d ", sa[i]);
  return 0;
}

有一点小问题，排后半截时会枚举到 (id_i+w > n) 怎么办？
考虑实际意义，出现此情况，表示该子串后半截为空。
空串字典序最小，考虑直接把 (rk) 开成两倍空间，则 (rk_i=0 (i>n)) 恒成立。防止了越界，也处理了空串的字典序。

再优化

两次计排太慢啦！ 观察对后半截排序时的特殊性质：

考虑更新前的 (sa_i) 的含义：排名为 (i) 的长度为 (2^{k-1}) 的子串 (S[sa_i, sa_i + 2^{k-1}])。
在本次排序中，(S[sa_i, sa_i + 2^{k-1}]) 是长度为 (2^k) 的子串 (S[sa_{i}-2^{k-1}:sa_i+2^{k-1}]) 的后半截，(sa_i) 的排名将作为排序的关键字。
(S[sa_i, sa_i + 2^{k-1}]) 的排名为 (i)，则第一次计排后 (S[sa_{i}-2^{k-1}:sa_i+2^{k-1}]) 的排名必为 (i)。考虑直接赋值，那么原来的第一次计排就可以写成这样：

int p = 0;
for (int i = n; i > n - w; -- i) id[++ p] = i; //后半截为空的串
for (int i = 1; i <= n; ++ i) { //根据后半截，直接推整个串的排名
  if (sa[i] > w) id[++ p] = sa[i] - w;
}

注意后半截为空串的情况，这样的串排名相同且最小。

以及一些奇怪的常数优化：

  • 减小值域。 值域大小 (m) 与排序复杂度有关，其最小值应为 (rk) 的最大值，更新 (rk) 时更新 (m) 即可。
  • 减少数组嵌套的使用，减少不连续内存访问。 在第二次计数排序时，将 (rk_{id_i}) 存下来。
  • 用 cmp 函数判断两个子串是否相同。同样是减少不连续内存访问，详见代码。
//知识点：SA
/*
By:Luckyblock
I love Marisa;
*/
#include <cstdio>
#include <ctype.h>
#include <cstring>
#include <algorithm>
#define ll long long
const int kMaxn = 1e6 + 10;
//=============================================================
char S[kMaxn];
int n, m, sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1];
int id[kMaxn], cnt[kMaxn], rkid[kMaxn];
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
bool cmp(int x, int y, int w) { //判断两个子串是否相等。
  return oldrk[x] == oldrk[y] && 
         oldrk[x + w] == oldrk[y + w]; 
}
//=============================================================
int main() {
  scanf("%s", S + 1);
  n = strlen(S + 1);
  m = std :: max(n, 300); //值域大小
  
  //初始化 sa数组
  for (int i = 1; i <= n; ++ i) ++ cnt[rk[i] = S[i]];
  for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
  for (int i = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;

  //倍增过程。 
  //此处 w = 2^{k-1}，是已经推出的子串长度。
  //注意此处的 sa 数组存的并不是后缀的排名，
  //存的是指定长度子串的排名。
  for (int p, w = 1; w < n; w <<= 1) {
    //按照后半截 rk[i+w] 作为第二关键字排序。
    p = 0;
    for (int i = n; i > n - w; -- i) id[++ p] = i; //后半截为空的串
    for (int i = 1; i <= n; ++ i) { //根据后半截，直接推整个串的排名
      if (sa[i] > w) id[++ p] = sa[i] - w;
    }

    //按照前半截 rk[i] 作为第一关键字排序。
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) ++ cnt[rkid[i] = rk[id[i]]];
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rkid[i]] --] = id[i];

    //更新 rk 数组。
    //这里可以滚动数组一下，但是可读性会比较差（
    //卡常可以写一下。
    std ::swap(rk, oldrk);
    m = 0; //直接更新值域 m
    for (int i = 1; i <= n; ++ i) {
      rk[sa[i]] = (m += (cmp(sa[i], sa[i - 1], w) ^ 1));
    }
  }
  for (int i = 1; i <= n; ++ i) printf("%d ", sa[i]);
  return 0;
}

线性构建

在大多数题目中，常数较小的倍增是完全够用的。走某些特殊题目中可以使用 DC3/SA-IS 算法实现线性构建后缀数组。
具体做法可以参考：诱导排序与 SA-IS 算法 与 DC3：[2009]后缀数组——处理字符串的有力工具 by. 罗穗骞。

LCP 问题

特别鸣谢：论文爷！后缀数组-许智磊

一些定义

(operatorname{lcp}(S,T)) 定义为字符串 (S) 和 (T) 的最长公共前缀 (Longest common prefix)， 即最大的 (lle min{left| S ight|,left| T ight|})，满足 (S_i=T_i(1le ile l))。
在许多后缀数组相关问题中，都需要它的帮助。

下文以 (operatorname{lcp}(i,j)) 表示后缀 (i)，(j) 的最长公共前缀，并延续上文中一些概念：(sa_i)：排名为 (i) 的后缀，(rk_i)：后缀 (i) 的排名。
并将会用 (sa_i) 直接代表排名为 (i) 的后缀，即 (sa_i = S[sa_i:n])。

定义一些新的概念。
(operatorname{height}_i) 表示排名为 (i-1) 和 (i) 的两后缀的最长公共前缀。

[operatorname{height}_i = operatorname{lcp}(sa_{i-1},sa_i) ]

(h_i) 表示后缀 (i) 和排名在 (i) 之前一位的后缀的最长公共前缀。

[h_i=operatorname{height}_{rk_i} = operatorname{lcp}(sa_{rk_i-1}, sa_{rk_i})= operatorname{lcp}(i, sa_{rk_i -1}) ]

引理：LCP Lemma

[forall 1le i<j<kle n, \,operatorname{lcp}(sa_i,sa_k) = min{operatorname{lcp}(sa_i,sa_j), operatorname{lcp}(sa_j,sa_k)} ]

此引理是证明其他引理的基础。
证明，设 (p = min{operatorname{lcp}(sa_i,sa_j), operatorname{lcp}(sa_j,sa_k)})，则有：

[operatorname{lcp}(sa_i,sa_j)ge p,\, operatorname{lcp}(sa_j,sa_k)ge p ]

则 (sa_i[1:p] = sa_j[1:p] = sa_k[1:p])，可得 (operatorname{lcp}(sa_i,sa_k)ge p)。

再考虑反证法，设 (operatorname{lcp}(sa_i,sa_k) =q > p)。则 (sa_i[1:q]=sa_k[1:q])，有 (sa_i[p+1]=sa_k[p+1])。
对 (p) 的取值分类讨论：

  1. (p=operatorname{lcp}(sa_i,sa_j) < operatorname{lcp}(sa_j,sa_k))：则有 (sa_i[p+1] < sa_j[p+1] = sa_k[p+1])。
  2. (p=operatorname{lcp}(sa_j,sa_k) < operatorname{lcp}(sa_i,sa_j))：则有 (sa_i[p+1] = sa_j[p+1] < sa_k[p+1])。
  3. (p=operatorname{lcp}(sa_j,sa_k) = operatorname{lcp}(sa_i,sa_j))：则有 (sa_i[p+1] < sa_j[p+1] < sa_k[p+1])。

(sa_i[p+1]<sa_k[p+1]) 恒成立，与已知矛盾，则 (operatorname{lcp}(sa_i,sa_k)le p)。
综合上述两个结论，得证引理成立。

引理：LCP Theorem

[forall 1le i < jle n,\, operatorname{lcp}(sa_i,sa_j) = min_{k=i+1}^j{operatorname{height_k}} ]

由 LCP Lemma，可知显然成立。

根据这个优美的式子，求解任意两个后缀的 (operatorname{lcp}) 变为求解 (operatorname{height}) 的区间最值问题。
可通过 st 表 实现 (O(nlog n)) 预处理，(O(1)) 查询。
问题只剩下如何快速求 (operatorname{height}) 了。

推论：LCP Corollary

[operatorname{lcp}(sa_i,sa_j) ge operatorname{lcp}(sa_i, sa_k)\, (ile j<k) ]

表示排名不相邻的两个后缀的 (operatorname{lcp}) 不超过它们之间任何相邻元素的 (operatorname{lcp})。
证明由引理 LCP Lemma 显然可得。
但是涛哥钦定我写一下证明，那我就不胜惶恐地写了（

类似 LCP Lemma，考虑反证法。设 (operatorname{lcp}(sa_i,sa_j)< operatorname{lcp}(sa_i, sa_k))，则有下图：

Lb

考虑字典序比较的过程。若 (sa_i < sa_j)，则有 (sa_i[{operatorname{lcp}(sa_i,sa_j)+1}] <sa_j[{operatorname{lcp}(sa_i,sa_j) + 1}])。
即图中的字符 (x<y)。

此时考虑比较 (sa_j) 与 (sa_k) 的字典序。由图，显然有 (operatorname{lcp}(sa_j,sa_k) = operatorname{lcp}(sa_i,sa_j))。而 (operatorname{lcp}(sa_i,sa_k) > operatorname{lcp}(sa_i,sa_j))，则 (sa_k[{operatorname{lcp}(sa_j,sa_k)+1}] = x)。
又 (x<y)，可得 (sa_k) 的字典序小于 (sa_j)。

与已知矛盾，反证原结论成立。

引理

[forall 1le ile n,\, h_ige h_{i-1}-1 ]

[h_i=operatorname{height}_{rk_i} = operatorname{lcp}(sa_{rk_i-1}, sa_{rk_i})= operatorname{lcp}(i, sa_{rk_i -1}) ]

用来快速计算 (operatorname{height}) 的引理，个人喜欢叫它不完全单调性。
证明考虑数学归纳。首先当 (h_{i-1}le 1) 时，结论显然成立，因为 (h_i ge 0)。

当 (h_{i-1}>1) 时，设 (u = i, \, v = sa_{rk_i-1})，有 (h_i = operatorname{lcp}(u,v))。同时，设 (u' = i-1, \, v' = sa_{rk_{i-1}-1})，有 (h_{i-1} = operatorname{lcp}(u',v'))。
由 (h_{i-1} = operatorname{lcp}(u',v')>1)，则 (u',v') 必有公共前缀。

考虑 删去 (u',v') 的 第一个 字符，设其分别变成 (x,y)。显然 (operatorname{lcp}(x,y) = h_{i-1}-1)，且仍满足字典序 (y<x)。
(u' = i-1)，则删去第一个字符后，(x) 等于后缀 (i)。
则对于他们在 (sa) 中的排名，有 (y<x=i=u)。

又 (sa) 中，(v) 在 (u) 前一位置，则有 (yle v<u)。根据 LCP Corollary，有：

[h_i = operatorname{lcp}(u,v)ge operatorname{lcp}(u,y) = operatorname{lcp}(x,y) = h_{i-1}-1 ]

得证。

快速求 height

由 定义 (h_i = operatorname{height}_{rk_i})，只需快速求出 (h)，便可 (O(n)) 地获得 (operatorname{height})。
由引理已知 (forall 1le ile n,\, h_ige h_{i-1}-1)。
(h_i=operatorname{lcp}(i, sa_{rk_i -1})) 具有不完全单调性，考虑正序枚举 (i) 进行递推。

当 (rk_i=1) 时， (sa_{rk_i-1}) 不存在，特判 (h_i=0)。
当 (i=1)，暴力比较出 (operatorname{lcp}(i,sa_{rk_i-1}))，比较次数 (<n)。
若上述情况均不满足，由引理知，(h_i=operatorname{lcp}(i,sa_{rk_i-1})ge h_{i-1}-1)，两后缀前 (h_{i-1}-1) 位相同。
可从第 (h_{i-1}) 位开始比较两后缀计算出 (h_i)，比较次数 (=h_i-h_{i-1}+2)

代码中并没有专门开 (h) 数组，其中(h_i = k)：

void GetHeight() {
  for (int i = 1, k = 0; i <= n; ++ i) {
    if (rk[i] == 1) k = 0;
    else {
      if (k > 0) k --;
      int j = sa[rk[i] - 1];
      while (i + k <= n && j + k <= n && 
             S[i + k] == S[j + k]) {
        ++ k;
      }
    }
    height[rk[i]] = k;
  }
}

(kle n)，最多减 (n) 次，则最多会在比较中加 (2n) 次。总复杂度为 (O(n)) 级别。

例题

「JSOI2007」字符加密

无法简述的题面。

断环成链，把字符串复制一遍扔到后面，跑 SA 即可。
板子背诵检查，可以练下手感。

SP705 SUBST1 - New Distinct Substrings

(T) 组数据，每次给定一个字符串 (s)，求该字符串本质不同的子串数量。
两个子串本质不同，当且仅当两个子串长度不等，或长度相等但有任意一位不同。
(1le Tle 1le|s|le 5 imes 10^4)。
280ms，1.46GB。

一种想法是用所有子串的个数 (frac{n(n+1)}{2}) 减去重复子串的个数，显然重复的串一定出现在某两个后缀的公共前缀部分。

考虑加入 (sa_i) 后，新增的本质不同的子串的数量，显然即 (operatorname{length}(sa_i) - operatorname{length}(operatorname{lcp}(sa_i, sa_{i-1})))，代表不作为之前加入的后缀的前缀的，(sa_i) 的前缀的数量。最终答案即：

[frac{n(n+1)}{2} - sum_{i = 2}^{n}operatorname{height}_i ]

SA 简单实现即可，总复杂度 (O(nlog n))。

如果想了解直观的证明解释可以阅读这篇文章：「笔记」后缀树。

//知识点：SA 
/*
By:Luckyblock
*/
#include <algorithm>
#include <cctype>
#include <cstdio>
#include <cstring>
#define LL long long
const int kN = 1e5 + 10;
//=============================================================
char s[kN];
int n, m, sa[kN], rk[kN << 1], oldrk[kN << 1], height[kN];
int id[kN], cnt[kN], rkid[kN];
//=============================================================
inline int read() {
  int f = 1, w = 0;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) {
    w = (w << 3) + (w << 1) + (ch ^ '0');
  }
  return f * w;
}
void Chkmax(int &fir_, int sec_) {
  if (sec_ > fir_) fir_ = sec_;
}
void Chkmin(int &fir_, int sec_) {
  if (sec_ < fir_) fir_ = sec_;
}
bool cmp(int x_, int y_, int w_) {
  return oldrk[x_] == oldrk[y_] && 
         oldrk[x_ + w_] == oldrk[y_ + w_];
}
void GetHeight() {
  for (int i = 1, k = 0; i <= n; ++ i) {
    if (rk[i] == 1) k = 0;
    else {
      if (k > 0) -- k;
      int j = sa[rk[i] - 1];
      while (i + k <= n && j + k <=n && 
             s[i + k] == s[j + k]) {
               ++ k;
      }
    }
    height[rk[i]] = k;
  }
}
void SuffixSort() {
  scanf("%s", s + 1);
  m = 300;
  n = strlen(s + 1);
  
  memset(cnt, 0, sizeof (cnt));
  for (int i = 1; i <= n; ++ i) cnt[rk[i] = s[i]] ++;
  for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
  for (int i = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;
  
  for (int p, w = 1; w < n; w <<= 1) {
    p = 0;
    for (int i = n; i > n - w; -- i) id[++ p] = i;
    for (int i = 1; i <= n; ++ i) {
      if (sa[i] > w) id[++ p] = sa[i] - w;
    }
    
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) cnt[rkid[i] = rk[id[i]]] ++;
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rkid[i]] --] = id[i];
    
    m = 0;
    memcpy(oldrk, rk, sizeof (rk));
    for (int i = 1; i <= n; ++ i) {
      m += (cmp(sa[i], sa[i - 1], w) ^ 1);
      rk[sa[i]] = m;
    }
  }
  GetHeight();
}
//=============================================================
int main() {
  int T = read();
  while (T --) {
    SuffixSort();
    LL ans = 1ll * n * (n + 1) / 2ll; 
    for (int i = 1; i <= n; ++ i) ans -= height[i];
    printf("%lld
", ans);
  }
  return 0;
}

SP1811 LCS - Longest Common Substring

给定两字符串 (S_1, S_2)，求它们的最长公共子串长度。
(|S_1|,|S_2|le 2.5 imes 10^5)。
294ms，1.46GB。

套路地把两个字符串连起来，答案即：

[max_{1le ile |S_1| < jle |S_1+S_2|}operatorname{lcp}(i,j) ]

显然答案即为满足 (sa_i,sa_{i-1}) 分属不同字符串 的 (operatorname{height}_{i}) 的最大值。
正确性非常显然，留与读者自证。这里给出一种证明，可以参考这里：「双串最长公共子串」。

//知识点：SA
/*
By:Luckyblock
*/
#include <algorithm>
#include <cstdio>
#include <cstring>
#include <ctype.h>
#define ll long long
const int kMaxn = 5e5 + 10;
//=============================================================
char S[kMaxn];
int n1, n, m, ans, cnt[kMaxn], id[kMaxn], rkid[kMaxn];
int sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1], height[kMaxn];
int MaxHeight[kMaxn][20], Log2[kMaxn];
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
void GetMax(int &fir, int sec) {
  if (sec > fir) fir = sec;
}
bool cmp(int x, int y, int w) { //判断两个子串是否相等。
  return oldrk[x] == oldrk[y] && 
         oldrk[x + w] == oldrk[y + w]; 
}
void GetHeight() {
  for (int i = 1, k = 0; i <= n; ++ i) {
    if (rk[i] == 1) k = 0;
    else {
      if (k > 0) k --;
      int j = sa[rk[i] - 1];
      while (i + k <= n && j + k <= n && 
             S[i + k] == S[j + k]) {
        ++ k;
      }
    }
    height[rk[i]] = k;
  }
}
void SuffixSort() {
  m = 300;
  for (int i = 1; i <= n; ++ i) ++ cnt[rk[i] = S[i]];
  for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
  for (int i = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;
  for (int p, w = 1; w < n; w <<= 1) {
    p = 0;
    for (int i = n; i > n - w; -- i) id[++ p] = i;
    for (int i = 1; i <= n; ++ i) {
      if (sa[i] > w) id[++ p] = sa[i] - w;
    }
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) ++ cnt[(rkid[i] = rk[id[i]])];
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rkid[i]] --] = id[i];
    std ::swap(rk, oldrk);
    m = 0;
    for (int i = 1; i <= n; ++ i) {
      m += (cmp(sa[i], sa[i - 1], w) ^ 1);
      rk[sa[i]] = m;
    }
  }
  GetHeight();
}
bool Judge(int x, int y) {
  return (x <= n1 && y > n1) || (x > n1 && y < n1);
}
//=============================================================
int main() {
  scanf("%s", S + 1); n1 = strlen(S + 1);
  S[n1 + 1] = 'z' + 1;
  scanf("%s", S + n1 + 1 + 1); n = strlen(S + 1);
  SuffixSort();
  for (int i = 2; i <= n; ++ i) {
    if (Judge(sa[i - 1], sa[i])) GetMax(ans, height[i]);
  }
  printf("%d", ans);
  return 0;
}

「HAOI2016」找相同字符

给定两字符串 (S_1, S_2)，求出在两字符串中各取一个子串，使得这两个子串相同的方案数。
两方案不同当且仅当这两个子串中有一个位置不同。
(1le |S_1|, |S_2|le 2 imes 10^5)。
1S，256MB。

考察对 (operatorname{lcp}) 单调性的理解。
(S_1) 后面加个终止符，(S_2) 串扔到 (S_1) 后面，跑 SA。
显然答案即后半段的后缀，与前半段的后缀的所有 (operatorname{lcp}) 之和。

考虑按字典序枚举后半段的后缀，设当前枚举到的后缀为 (sa_i)。
先 仅考虑 字典序 (<sa_i) 的 前半段的后缀 (sa_j (j<i)) 的贡献。其对 (sa_i) 的贡献为 (operatorname{lcp}(sa_i, sa_j))。

由 (operatorname{lcp}) 的单调性，对于 最小 的大于 (sa_i) 的后半段的后缀 (sa_k(k>i))，有 (operatorname{lcp}(sa_{k}, sa_j)le operatorname{lcp}(sa_i,sa_j))，考虑贡献的变化情况。

若 (operatorname{lcp}(sa_{k}, sa_j)< operatorname{lcp}(sa_i,sa_j))，则 (sa_j) 对 (sa_k) 的贡献应变为：

[operatorname{lcp}(sa_k, sa_j) = min{operatorname{lcp}(sa_i,sa_j), minlimits_{l=i+1}^{k}{operatorname{height}_l}} ]

此外，若存在 (sa_l, lin (i,k)) 为 前半段的后缀 时，作出贡献的元素增加。

考虑在枚举后缀的过程中，用权值线段树维护 字典序小于 (sa_i) 的 前半段 的后缀 (sa_j (j<i)) 的不同长度的 (operatorname{lcp}) 的 数量。
上述两操作，即为区间赋值 与 单点插入。

再按字典序倒序枚举后缀，计算字典序 (>sa_i) 的 前半段的后缀的贡献。
分析很屑，代码有详细注释。

总复杂度 (O(nlog n))。线段树写法是自己 YY 的，比较无脑，也可以单调栈简单维护，复杂度也为 (O(nlog n)) 级别。
此外还有优美的广义 SAM 写法，可以参考：「HAOI2016」找相同字符。

//知识点：SA，线段树
/*
By:Luckyblock 
*/
#include <cstdio>
#include <ctype.h>
#include <cstring>
#include <algorithm>
#define ll long long
#define lson (now_<<1)
#define rson (now_<<1|1)
const int kMaxn = 4e5 + 10;
//=============================================================
char S[kMaxn];
int n1, n, m, sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1], height[kMaxn];
int id[kMaxn], cnt[kMaxn], rkid[kMaxn];
ll ans, size[kMaxn << 2], sum[kMaxn << 2]; //size 维护数量，sum 维护 lcp 之和。
bool tag[kMaxn << 2];
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
bool cmp(int x, int y, int w) { //判断两个子串是否相等。
  return oldrk[x] == oldrk[y] && 
         oldrk[x + w] == oldrk[y + w]; 
}
void GetHeight() {
  for (int i = 1, k = 0; i <= n; ++ i) {
    if (rk[i] == 1) k = 0;
    else {
      if (k > 0) k --;
      int j = sa[rk[i] - 1];
      while (i + k <= n && j + k <= n && 
             S[i + k] == S[j + k]) {
        ++ k;
      }
    }
    height[rk[i]] = k;
  }
}
void SuffixSort() {
  m = std :: max(n, 300);
  for (int i = 1; i <= n; ++ i) ++ cnt[rk[i] = S[i]];
  for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
  for (int i = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;
  for (int p, w = 1; w < n; w <<= 1) {
    p = 0;
    for (int i = n; i > n - w; -- i) id[++ p] = i;
    for (int i = 1; i <= n; ++ i) {
      if (sa[i] > w) id[++ p] = sa[i] - w;
    }
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) ++ cnt[(rkid[i] = rk[id[i]])];
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rkid[i]] --] = id[i];
    std ::swap(rk, oldrk);
    m = 0;
    for (int i = 1; i <= n; ++ i) {
      m += (cmp(sa[i], sa[i - 1], w) ^ 1);
      rk[sa[i]] = m;
    }
  }
  GetHeight();
}
void Build(int now_, int L_, int R_) {
  size[now_] = sum[now_] = 0ll;
  tag[now_] = false;
  if (L_ == R_) return ;
  int mid = (L_ + R_) >> 1;
  Build(lson, L_, mid), Build(rson, mid + 1, R_);
}
void Pushdown(int now_) {
  tag[lson] = tag[rson] = true;
  size[lson] = size[rson] = 0;
  sum[lson] = sum[rson] = 0;
  tag[now_] = false;
}
void Pushup(int now_) {
  size[now_] = size[lson] + size[rson];
  sum[now_] = sum[lson] + sum[rson];
}
ll Delete(int now_, int L_, int R_, int ql_, int qr_) {
  if (ql_ <= L_ && R_ <= qr_) {
    ll ret = size[now_];
    tag[now_] = true;
    size[now_] = sum[now_] = 0;
    return ret;
  }
  if(tag[now_]) Pushdown(now_);
  int mid = (L_ + R_) >> 1;
  ll ret = 0ll;
  if (ql_ <= mid) ret += Delete(lson, L_, mid, ql_, qr_);
  if (qr_ > mid) ret += Delete(rson, mid + 1, R_, ql_, qr_);
  Pushup(now_);
  return ret;
}
void Insert(int now_, int L_, int R_, int pos_, ll num) {
  if (! num) return ;
  if (L_ == R_) {
    size[now_] += num;
    sum[now_] += 1ll * num * (L_ - 1ll); //注意减去偏移量。
    return ;
  }
  if (tag[now_]) Pushdown(now_);
  int mid = (L_ + R_) >> 1;
  if (pos_ <= mid) Insert(lson, L_, mid, pos_, num);
  else Insert(rson, mid + 1, R_, pos_, num);
  Pushup(now_);
}
//=============================================================
int main() {
  scanf("%s", S + 1); n1 = strlen(S + 1);
  S[n1 + 1] = 'z' + 1;
  scanf("%s", S + n1 + 2); n = strlen(S + 1);
  SuffixSort();

  //正序枚举所有后缀，计算字典序 >sa_i 的 前半段的后缀的贡献。
  //当枚举到一个 后半段的后缀，仅用于更新 min(lcp)。
  //枚举到一个 前半段的后缀，用于更新 min(lcp)，且需新插入一个后缀。
  //由于 lcp 可能为 0，线段树维护的区间加了偏移量 1。
  for (int i = 2; i <= n; ++ i) {
    //计算 lcp > height(i) 的 前半段后缀的数量，并将他们删除。
    ll num = Delete(1, 1, n + 1, height[i] + 1 + 1, n + 1); 
    Insert(1, 1, n + 1, height[i] + 1, num + (sa[i - 1] <= n1)); //插入被删除的后缀 与 新后缀。注意边界。
    if (sa[i] > n1 + 1) ans += sum[1]; //若枚举到一个 后半段后缀，计算贡献。 注意边界。
  }
  Build(1, 1, n + 1); //清空线段树
  //倒序枚举所有后缀，计算字典序 >sa_i 的 前半段的后缀的贡献。
  for (int i = n; i >= 2; -- i) {
    ll num = Delete(1, 1, n + 1, height[i] + 2, n + 1);
    Insert(1, 1, n + 1, height[i] + 1, num + (sa[i] <= n1)); //注意边界
    if (sa[i - 1] > n1 + 1) ans += sum[1]; //注意边界
  }
  printf("%lld", ans);
  return 0;
}

「AHOI2013」 差异

给定一长度为 (n) 的字符串 (S)，令 (T_i) 表示从第 (i) 个字符开始的后缀，求：

[sum_{1le i<jle n}{operatorname{len}(T_i) +operatorname{len}(T_j) - 2 imes operatorname{lcp} (T_i,T_j)} ]

(operatorname{len}(a)) 表示字符串 (a) 的长度，(operatorname{lcp}(a,b)) 表示字符串 (a,b) 的最长公共前缀。
(1le nle 5 imes 10^5)。
1S，512MB。

化下式子：

[egin{aligned} &sum_{1le i<jle n}{operatorname{len}(T_i) +operatorname{len}(T_j) - 2 imes operatorname{lcp} (T_i,T_j)}\ =& sum_{1le i<jle n}{(n-i+1) +(n-j+1) - 2 imes operatorname{lcp} (T_i,T_j)}\ =& dfrac{(n-1) imes n imes (n+1)}{2} - 2sum_{1le i<jle n}operatorname{lcp} (T_i,T_j) end{aligned}]

考虑如何快速求后一半，即所有 (operatorname{lcp}) 之和。
各后缀与 (sa) 数组的元素，一一对应，则有下列等价关系：

[sum_{1le i<jle n}operatorname{lcp} (T_i,T_j) = sum_{1le i<jle n}operatorname{lcp}(sa_i, sa_j) ]

类似上题的套路，考虑枚举 (sa_j)，用权值线段树维护 (sa_i (i<j)) 的不同长度的 (operatorname{lcp}(sa_i, sa_j)) 的数量。

有一引理：

[forall 1le i < jle n,\, operatorname{lcp}(sa_i,sa_j) = minlimits_{k=i+1}^j{operatorname{height_k}} ]

模拟引理，当 (j+1) 时将权值线段树中所有 (>operatorname{height}_{j+1}) 的元素删除，并添加相同个数个 元素 (operatorname{height}_{j+1})。
添加一个 (operatorname{height}_{j+1})，代表新增的 (sa_j) 的贡献。
贡献求和即可，总复杂度 (O(nlog n))。

线段树太傻逼了，考虑单调数据结构。发现有下列等价关系：

[sum_{1le i<jle n}operatorname{lcp}(sa_i, sa_j) = sum_{1le i<jle n}min_{k=i+1}^{j}{operatorname{height}_k} ]

即求 (operatorname{height}) 中所有子区间的区间最小值之和。
这是个经典问题，单调栈维护 (operatorname{height}) 作为最小值的区间的 最大左右端点，答案即 (sumlimits_{i=2}^{n}(i-l_i) imes (r_i-i) imes operatorname{height}_i)。
总复杂度 (O(nlog n))。注意特判区间长度不能为 1。

线段树代码可以参考：这里。

//知识点：SA，单调性
/*
By:Luckyblock
*/
#include <cstdio>
#include <ctype.h>
#include <cstring>
#include <iostream>
#include <algorithm>
#define ll long long
const int kMaxn = 5e5 + 10;
//=============================================================
char S[kMaxn];
int n, m, sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1], height[kMaxn];
int cnt[kMaxn], id[kMaxn], rkid[kMaxn];
int top, st[kMaxn], l[kMaxn], r[kMaxn];
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
void GetMax(int &fir, int sec) {
  if (sec > fir) fir = sec;
}
void GetMin(int &fir, int sec) {
  if (sec < fir) fir = sec;
}
int cmp(int x, int y, int w) {
  return oldrk[x] == oldrk[y] && 
         oldrk[x + w] == oldrk[y + w];
}
void GetHeight() {
  for (int i = 1, k = 0; i <= n; ++ i) {
    if (rk[i] == 1) k = 0;
    else {
      if (k > 0) k --;
      int j = sa[rk[i] - 1];
      while (i + k <= n && j + k <= n && 
             S[i + k] == S[j + k]) {
        ++ k;
      }
    }
    height[rk[i]] = k;
  }
}
void SuffixSort() {
  n = strlen(S + 1);
  m = 1010;
  for (int i = 1; i <= n; ++ i) cnt[rk[i] = S[i]] ++;
  for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
  for (int i = n; i; -- i) sa[cnt[rk[i]] --] = i;
  
  for (int p, w = 1; w < n; w <<= 1) {
    p = 0;
    for (int i = n; i > n - w; -- i) id[++ p] = i;
    for (int i = 1; i <= n; ++ i) {
      if (sa[i] > w) id[++ p] = sa[i] - w;
    }
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) cnt[rkid[i] = rk[id[i]]] ++;
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i; -- i) sa[cnt[rkid[i]] --] = id[i]; 
    
    std :: swap(rk, oldrk);
    m = 0;
    for (int i = 1; i <= n; ++ i) {
      m += (cmp(sa[i], sa[i - 1], w) ^ 1);
      rk[sa[i]] = m;
    }
  }
  GetHeight();
}
//=============================================================
int main() {
  scanf("%s", S + 1);
  SuffixSort();
  ll ans = 1ll * ((n - 1ll) * n / 2ll) * (n + 1ll) ;
  st[(top = 1)] = 1;
	for (int i = 2; i <= n; ++ i) {
		while (top && height[st[top]] > height[i]) {
		  r[st[top]] = i;
		  top --;
    }
		l[i] = st[top];
		st[++ top] = i;
	} 
  while (top) r[st[top --]] = n + 1;
  for (int i = 2; i <= n; ++ i) {
    ans -= 2ll * (i - l[i]) * (r[i] - i) * height[i]; 
  }
  printf("%lld", ans);
  return 0;
}

「TJOI / HEOI2016」字符串

给定一长度为 (n) 的字符串，(m) 个询问。
每次询问给定参数 (a,b,c,d)，求子串 (S[a:b]) 的所有子串，与子串 (S[c:d]) 的最长公共前缀的最大值。
(1le n,mle 10^5, 1le ale ble n, 1le cle dle n)。
2S，256MB。

对于每一个询问，答案满足单调性，考虑二分答案。
设 (l) 为当前二分到的最长的，子串 (S[a:b]) 的 所有子串，与子串 (S[c:d]) 的最长公共前缀。
分析可知，若 (l) 合法，那么会存在至少一个后缀 (x) 满足：

  • 开头在 ([a:b-l+1]) 中。
  • (operatorname{lcp}(x, c)ge l)。

对于第二个限制，考虑 (operatorname{lcp}) 的单调性。(operatorname{lcp}(x,c)) 是一个在 (c) 处取极值的单峰函数。
则满足条件的 (x) 的取值，一定是 (sa) 数组上连续的一段。
套路地对 (operatorname{height}) 建立 st 表，即可 (O(1)) 查询 (operatorname{lcp})，于是可以通过二分排名快速得到后缀 (x) 排名的取值范围，将限制二也转化为了区间限制形式。

限制一限制了后缀的区间，限制二限制了 (rk) 的区间。查询这样的后缀的存在性，变成了一个静态二维数点问题。
对 (sa) 数组建立主席树维护即可，总复杂度 (O(nlog^2 n))。

//知识点：SA，二分答案，主席树
/*
By:Luckyblock
*/
#include <algorithm>
#include <cstdio>
#include <cstring>
#include <ctype.h>
#define ll long long
const int kMaxn = 1e5 + 10;
//=============================================================
char S[kMaxn];
int n, m, ans, cnt[kMaxn], id[kMaxn], rkid[kMaxn];
int sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1], height[kMaxn];
int MaxHeight[kMaxn][20], Log2[kMaxn];;
int lson[kMaxn << 5], rson[kMaxn << 5], size[kMaxn << 5];
int node_num, root[kMaxn];
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
void GetMax(int &fir, int sec) {
  if (sec > fir) fir = sec;
}
bool cmp(int x, int y, int w) { //判断两个子串是否相等。
  return oldrk[x] == oldrk[y] && 
         oldrk[x + w] == oldrk[y + w]; 
}
void GetHeight() {
  for (int i = 1, k = 0; i <= n; ++ i) {
    if (rk[i] == 1) k = 0;
    else {
      if (k > 0) k --;
      int j = sa[rk[i] - 1];
      while (i + k <= n && j + k <= n && 
             S[i + k] == S[j + k]) {
        ++ k;
      }
    }
    height[rk[i]] = k;
  }
}
int QueryLcp(int l_, int r_) {
  int k = Log2[r_ - l_ + 1];
  return std :: min(MaxHeight[l_][k], MaxHeight[r_ - (1 << k) + 1][k]);
}
void MakeSt() {
  for (int i = 2; i <= n; ++ i) MaxHeight[i][0] = height[i];
  for (int i = 2; i <= n; ++ i) {
    Log2[i] = Log2[i >> 1] + 1;
  }
  for (int j = 1; j < 20; ++ j) {
    for (int i = 1; i + (1 << j) - 1 <= n; ++ i) {
      MaxHeight[i][j] = std :: min(MaxHeight[i][j - 1], 
                                   MaxHeight[i + (1 << (j - 1))][j - 1]);
    }
  }
}
void SuffixSort() {
  int M = 300;
  for (int i = 1; i <= n; ++ i) ++ cnt[rk[i] = S[i]];
  for (int i = 1; i <= M; ++ i) cnt[i] += cnt[i - 1];
  for (int i = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;
  for (int p, w = 1; w < n; w <<= 1) {
    p = 0;
    for (int i = n; i > n - w; -- i) id[++ p] = i;
    for (int i = 1; i <= n; ++ i) {
      if (sa[i] > w) id[++ p] = sa[i] - w;
    }
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) ++ cnt[(rkid[i] = rk[id[i]])];
    for (int i = 1; i <= M; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rkid[i]] --] = id[i];
    std ::swap(rk, oldrk);
    M = 0;
    for (int i = 1; i <= n; ++ i) {
      M += (cmp(sa[i], sa[i - 1], w) ^ 1);
      rk[sa[i]] = M;
    }
  }
  GetHeight();
  MakeSt();
}
void Insert(int pre_, int &now_, int L_, int R_, int val_) {
  now_ = ++ node_num;
  size[now_] = size[pre_] + 1;
  lson[now_] = lson[pre_], rson[now_] = rson[pre_];
  if(L_ >= R_) return ;
  int mid = (L_ + R_) >> 1;
  if(val_ > mid) Insert(rson[pre_], rson[now_], mid + 1, R_, val_);
  else Insert(lson[pre_], lson[now_], L_, mid, val_);
}
int Query(int u_, int v_, int L_, int R_, int ql_, int qr_) {
  if (ql_ <= L_ && R_ <= qr_) return size[v_] - size[u_];
  int mid = (L_ + R_) >> 1, ret = 0;
  if (ql_ <= mid) ret += Query(lson[u_], lson[v_], L_, mid, ql_, qr_);
  if (qr_ > mid) ret += Query(rson[u_], rson[v_], mid + 1, R_, ql_, qr_);
  return ret;
}
bool Judge(int len_, int a_, int b_, int c_) {
  int l = 1, r = rk[c_], up, down;
  while (l < r) {
    int mid = (l + r) >> 1;
    if (QueryLcp(mid + 1, rk[c_]) < len_) l = mid + 1;
    else r = mid;
  }
  up = r, l = rk[c_], r = n;
  while (l < r) {
    int mid = (l + r + 1) >> 1;
    if (QueryLcp(rk[c_] + 1, mid) < len_) r = mid - 1;
    else l = mid;
  }
  down = r;
  return Query(root[up - 1], root[down], 1, n, a_, b_ - len_ + 1) > 0;
}
int Solve(int a_, int b_, int c_, int d_) {
  int l = 0, r = std :: min(b_ - a_ + 1, d_ - c_ + 1);
  while (l < r) {
    int len = (l + r + 1) >> 1;
    if (Judge(len, a_, b_, c_)) l = len;
    else r = len - 1;
  }
  return r;
}
//=============================================================
int main() {
  n = read(), m = read();
  scanf("%s", S + 1);
  SuffixSort();
  for (int i = 1; i <= n; ++ i) Insert(root[i - 1], root[i], 1, n, sa[i]);
  for (int i = 1; i <= m; ++ i) {
    int a = read(), b = read(), c = read(), d = read();
    printf("%d
", Solve(a, b, c, d));
  }
  return 0;
}

「NOI2015」品酒大会

「そして谁もいなくなるか？」

给定一字符串 (S)，位置 (i) 的属性值为 (a_i)。
定义位置 (p,q) 为「 (r) 相似」，当且仅当 (S[p:p+r-1] = S[q:q+r-1]) 成立。
特别地，对于任意 (1le p,qle n, p ot ={q})，(p,q) 都是「 (0) 相似」的。
求：选出两个 「 (0sim n-1) 相似」 的位置的 方案数，及选出两个 「 (0sim n-1) 相似」的位置属性值 乘积的最大值。
1S，256MB。

显然，若两个位置 (i,j) 是「(r) 相似」的，那么它们也是「(0sim (r-1)) 相似」的。它们会对「(0sim r) 相似」的答案做出贡献，于是考虑倒序枚举「(r) 相似」的位置计算。

发现「(r) 相似」问题的实质即 (operatorname{lcp}) 问题，先对原串跑 SA，求得 (operatorname{height}) 数组。
考虑「(r) 相似」的定义，则对于任意两位置 (i,j)，他们最大是 「(operatorname{lcp}(S[i:n],S[j:n])) 相似」的。

根据引理：

[forall 1le i < jle n,\, operatorname{lcp}(sa_i,sa_j) = min_{k=i+1}^j{operatorname{height_k}} ]

考虑按照 (operatorname{height}) 将后缀排序后的后缀进行划分。
若 (operatorname{height}_ige r)，将 (sa_{i-1}) 与 (sa_i) 划入一个集合，否则划入不同的集合。这样划分后，对于所有大小 (ge 2) 的集合，集合中所有后缀两两的 (operatorname{lcp}ge r)，它们都会对 「(r) 相似」的答案做出贡献。
将所有集合的贡献合并，即为 「(r) 相似」的答案。

定义上述划分方式为 「(r) 划分」，考虑如何在此基础上获得 「(r-1) 划分」。
显然，只需将 (operatorname{height}_i = r-1) 的后缀 (sa_{i-1}) 和 (sa_i) 所在集合合并即可。
考虑将 (operatorname{height}) 降序排序，用并查集维护集合，按上述过程依次进行合并，即可依次得到「(nsim 1) 相似」的答案。

考虑如何维护集合的贡献。

先考虑第一问，对于「(r) 划分」中一个大小 (ge 2) 的集合，集合中任意两后缀的 (operatorname{lcp} ge r)，该集合的贡献即为 ((operatorname{size}-1) imes operatorname{size})。
合并时直接将 (operatorname{size}) 累加即可。

再考虑第二问，由于可能存在 (a_i<0)，考虑维护集合的最大值，次大值，最小值，次小值。为保证答案合法，四个值中的任意两个都不能来自 同一个位置。
集合的贡献为 (max{max_1 imes max_2, min_1 imes min_2})。
合并时注意四个值的大小关系，代码中使用了 multiset 维护。

复杂度瓶颈是倍增，总复杂度 (O(nlog n)) 级别，使用炫酷 DC3/SA-IS 魔术可做到 (O(n)) 级别的复杂度。

//知识点：SA，并查集
/*
By:Luckyblock
*/
#include <algorithm>
#include <cstdio>
#include <cstring>
#include <ctype.h>
#include <set>
#define ll long long
const int kMaxn = 3e5 + 10;
//=============================================================
char S[kMaxn];
int n, m, a[kMaxn];
int cnt[kMaxn], id[kMaxn], rkid[kMaxn];
int sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1], height[kMaxn];
ll ans1[kMaxn], ans2[kMaxn];
int fa[kMaxn], size[kMaxn];
std :: multiset <int> s[kMaxn];
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
void GetMax(ll &fir, ll sec) {
  if (sec > fir) fir = sec;
}
bool cmp(int x, int y, int w) { //判断两个子串是否相等。
  return oldrk[x] == oldrk[y] && 
         oldrk[x + w] == oldrk[y + w]; 
}
void GetHeight() {
  for (int i = 1, k = 0; i <= n; ++ i) {
    if (rk[i] == 1) k = 0;
    else {
      if (k > 0) k --;
      int j = sa[rk[i] - 1];
      while (i + k <= n && j + k <= n && 
             S[i + k] == S[j + k]) {
        ++ k;
      }
    }
    height[rk[i]] = k;
  }
}
void SuffixSort() {
  m = 300;
  for (int i = 1; i <= n; ++ i) ++ cnt[rk[i] = S[i]];
  for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
  for (int i = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;
  for (int p, w = 1; w < n; w <<= 1) {
    p = 0;
    for (int i = n; i > n - w; -- i) id[++ p] = i;
    for (int i = 1; i <= n; ++ i) {
      if (sa[i] > w) id[++ p] = sa[i] - w;
    }
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) ++ cnt[(rkid[i] = rk[id[i]])];
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rkid[i]] --] = id[i];
    std ::swap(rk, oldrk);
    m = 0;
    for (int i = 1; i <= n; ++ i) {
      m += (cmp(sa[i], sa[i - 1], w) ^ 1);
      rk[sa[i]] = m;
    }
  }
  GetHeight();
}
bool Compare(int fir, int sec) {
  return height[fir] > height[sec];
}
int Find(int x) {
  return fa[x] == x ? x : fa[x] = Find(fa[x]);
}
void Union(int x, int y, int z) {
  x = Find(x), y = Find(y);
  if (size[x] < size[y]) std :: swap(x, y);
  ans1[z] += 1ll * size[x] * size[y];
  for (std :: set <int> :: iterator it = s[y].begin(); it != s[y].end(); ++ it) {
    s[x].insert(* it);
  }
  int t[5];
  t[1] = *s[x].begin(), t[2] = *(++ s[x].begin());
  t[3] = *(-- s[x].end()), t[4] = *(-- (-- s[x].end()));
  GetMax(ans2[z], std :: max(1ll * t[1] * t[2], 1ll * t[3] * t[4]));
  fa[y] = x;
  size[x] += size[y];
  if (s[x].size() > 5) {
    s[x].clear();
    for (int i = 1; i <= 4; ++ i) s[x].insert(t[i]);
  }
}
//=============================================================
int main() {
  n = read();
  scanf("%s", S + 1);
  for (int i = 1; i <= n; ++ i) a[i] = read();
  SuffixSort();
  memset(ans2, - 63, sizeof (ans2));
  for (int i = 2; i <= n; ++ i) id[i] = i;
  for (int i = 1; i <= n; ++ i) fa[i] = i;
  for (int i = 1; i <= n; ++ i) size[i] = 1;
  std :: sort(id + 2, id + n + 1, Compare);
  for (int i = 1; i <= n; ++ i) s[i].insert(a[i]);
  for (int i = 2; i <= n; ++ i) {
    Union(sa[id[i] - 1], sa[id[i]], height[id[i]]);
  }
  for (int i = n - 1; i >= 0; -- i) {
    ans1[i] += ans1[i + 1];
    GetMax(ans2[i], ans2[i + 1]);
  }
  for (int i = 0; i < n; ++ i) {
    printf("%lld %lld
", ans1[i], ans1[i] ? ans2[i] : 0); 
  }
  return 0;
}

「十二省联考 2019」字符串问题

神笔出题人居然卡清空/jk
调一上午发现把昨天的 TLE 代码的邻接表换成 vector 就过了/jk
草草草草草

可怜金发小女孩

(T) 组数据，每次给定一字符串 (S)。
在 (S) 中存在 (n_a) 个 A 类子串 ((la_i, ra_i))，(n_b) 个 B 类子串 ((lb_i,rb_i))。且存在 (m) 组支配关系，支配关系 ((x,y)) 表示第 (x) 个 A 类串支配第 (y) 个 B 类串。
要求构造一个目标串 (T)，满足：

  • (T) 由若干 A 类串拼接而成。
  • 对于分割中所有相邻的串，后一个串存在一个被前一个串支配的前缀。

求该目标串的最大长度，若目标串可以无限长输出 (-1)。
(1le Tle 100)，(n_a,n_b,|S|,mle 2 imes 10^5)。
6S，1G。

首先建立图论模型，从每个 A 类子串向其支配的 B 串连边，从每个 B 串向以它为前缀的 A 串连边。A 串节点的权值为其长度，B 串节点权值为 0。
在图上 DP 求得最长路即为答案，若图中存在环则无解。
第一类边有 (m) 条，但第二类边数可以达到 (n_an_b) 级别，考虑优化建图。

对于某 A 串 ((la_i, ra_i))，它以 B 串 ((lb_j, rb_j)) 作为一个前缀的充要条件是 (operatorname{lcp}(S[la_i:n],S[lb_j:n]) ge rb_j-lb_j+1) 且 (ra_i - la_i + 1ge rb_j-lb_j+1)。
对于限制一，考虑求得 (S) 的 SA，对 (operatorname{height}) 建立 ST 表，可在 (sa) 上二分求得满足 (operatorname{lcp}ge rb_j-lb_j+1) 的区间的左右边界，满足条件的 A 串一定都在这段区间内。第二类边转化为区间连边问题。
此时不考虑限制二直接线段树优化建图，可以拿到 80 分的好成绩。

限制二实际上限定了 B 连边的对象的长度。
考虑将所有 A，B 串按长度递减排序，按长度递减枚举 A 串并依次加入可持久化线段树。
对于每一个 B 串，先找到使得 A 串长度大于其长度的最晚的历史版本，此时线段树上的所有 A 串长度都大于其长度，再向这棵线段树上的节点连边。

时间复杂度 (O((|S| + n_a + n_b)log n))，空间复杂度 (O(m + (|S| + n_a + n_b)log n))，不需要刻意卡常就能稳过。

看见上面轻描淡写的是不是觉得这题太傻逼了？以下是菜鸡 Lb 在代码实现上的小问题。

边数在极限数据下可以达到 (10^7) 级别，不注意空间大小和清空时的实现会被卡到 60。这个时空限制显然就是给选手乱搞的，数组往大了开就行。

在线段树优化建图中，实点会在建树操作中就与虚点连好边。本题的实点是代表 A，B 串的节点，在本题的可持久化线段树优化中，实点与虚点的连边发生在动态开点的插入过程中。
在新建节点时，需要将该节点连向上一个版本中对应位置的节点。
对于 A 串 ((la_i, ra_i))，它应该被插入到线段树中 (rk_{la_i}) 的位置，即叶节点 (rk_{la_i}) 与该实点相连。

(operatorname{height}_1) 没有意义。

注意二分时的初始值。

long long

函数传参顺序 是通过栈传递的，因此是从右向左的。

//知识点：SA，可持久化线段树，优化建图，DAGDP
/*
By:Luckyblock
*/
#include <algorithm>
#include <cctype>
#include <cstdio>
#include <cstring>
#include <queue>
#define LL long long
const int kN = 2e5 + 10;
//=============================================================
struct Str {
  int l, r, lth, id;
} subs[kN << 1];
int node_num, n, na, nb, m, into[kN <<5];
int e_num, head[kN << 5], v[50 * kN], ne[50 * kN];
LL val[kN << 5], f[kN << 5];
char s[kN];
//=============================================================
inline int read() {
  int f = 1, w = 0;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) {
    w = (w << 3) + (w << 1) + (ch ^ '0');
  }
  return f * w;
}
void Chkmax(LL &fir_, LL sec_) {
  if (sec_ > fir_) fir_ = sec_;
}
void Chkmin(int &fir_, int sec_) {
  if (sec_ < fir_) fir_ = sec_;
}
bool cmp(Str fir_, Str sec_) {
  if (fir_.lth != sec_.lth) return fir_.lth > sec_.lth;
  return fir_.id < sec_.id;
}
void Add(int u_, int v_) {
  v[++ e_num] = v_, ne[e_num] = head[u_], head[u_] = e_num;
  into[v_] ++;
}
namespace ST {
  int Minn[kN][21], Log2[kN];
  void MakeST(int *a_) {
    for (int i = 1; i <= n; ++ i) Minn[i][0] = a_[i];
    for (int i = 2; i <= n; ++ i) Log2[i] = Log2[i >> 1] + 1;
    for (int j = 1; j <= 20; ++ j) {
      for (int i = 1; i + (1 << j) - 1 <= n; ++ i) { //
        Minn[i][j] = std::min(Minn[i][j - 1], Minn[i + (1 << (j - 1))][j - 1]);
      }
    }
  }
  int Query(int l_, int r_) {
    int k = Log2[r_ - l_ + 1];
    return std::min(Minn[l_][k], Minn[r_ - (1 << k) + 1][k]);
  }
}
namespace SA {
  int sa[kN], rk[kN << 1];
  int oldrk[kN << 1], cnt[kN], id[kN];
  int height[kN];
  void SuffixSort() {
    int rknum = std::max(n, 300);
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) cnt[rk[i] = s[i]] ++;
    for (int i = 1; i <= rknum; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;
    
    for (int w = 1, p; w < n; w <<= 1) {
      p = 0;
      for (int i = n; i > n - w; -- i) id[++ p] = i;
      for (int i = 1; i <= n; ++ i) {
        if (sa[i] > w) id[++ p] = sa[i] - w;
      }
      
      memset(cnt, 0, sizeof (cnt));
      for (int i = 1; i <= n; ++ i) ++ cnt[rk[id[i]]];
      for (int i = 1; i <= rknum; ++ i) cnt[i] += cnt[i - 1];
      for (int i = n; i >= 1; -- i) sa[cnt[rk[id[i]]] --] = id[i];
      
      for (int i = 1; i <= n; ++ i) oldrk[i] = rk[i];
      rknum = 0;
      for (int i = 1; i <= n; ++ i) {
        rknum += (oldrk[sa[i]] == oldrk[sa[i - 1]] && 
                  oldrk[sa[i] + w] == oldrk[sa[i - 1] + w]) ^ 1;
        rk[sa[i]] = rknum;
      }
    }
  }
  void GetHeight() {
    for (int i = 1, k = 0; i <= n; ++ i) {
      if (rk[i] == 1) {
        k = 0;
      } else {
        if (k) -- k;
        int j = sa[rk[i] - 1];
        while (i + k <= n && j + k <= n && 
               s[i + k] == s[j + k]) {
          ++ k;
        }
      }
      height[rk[i]] = k;
    }
  }
  int Lcp(int x_, int y_) {
    if (x_ > y_) std::swap(x_, y_);
    return ST::Query(x_ + 1, y_);
  }
  void Init() {
    SuffixSort();
    GetHeight();
    ST::MakeST(SA::height);
  }
}
namespace Hjt {
  #define ls lson[now_]
  #define rs rson[now_]
  #define mid ((L_+R_)>>1)
  int root[kN], lson[kN << 5], rson[kN << 5];
  void Insert(int &now_, int pre_, int L_, int R_, int pos_, int id_) {
    now_ = ++ node_num;
    ls = lson[pre_], rs = rson[pre_];
    if (pre_) Add(now_, pre_);
    if (L_ == R_) {
      Add(now_, id_);
      return ;
    }
    if (pos_ <= mid) {
      Insert(ls, lson[pre_], L_, mid, pos_, id_);
      Add(now_, ls);
    } else {
      Insert(rs, rson[pre_], mid + 1, R_, pos_, id_);
      Add(now_, rs);
    }
  }
  void AddEdge(int now_, int L_, int R_, int l_, int r_, int id_) {
    if (! now_) return ;
    if (l_ <= L_ && R_ <= r_) {
      Add(id_, now_);
      return ;
    }
    if (l_ <= mid) AddEdge(ls, L_, mid, l_, r_, id_);
    if (r_ > mid) AddEdge(rs, mid + 1, R_, l_, r_, id_);
  }
  #undef ls
  #undef rs
  #undef mid
}
void Init() {
  e_num = 0;
  for (int i = 0; i <= node_num; ++ i) {
    head[i] = val[i] = into[i] = f[i] = 0;
  }
  
  scanf("%s", s + 1);
  n = strlen(s + 1);
  SA::Init();
  na = read();
  for (int i = 1; i <= na; ++ i) {
    int l_ = read(), r_ = read();
    subs[i] = (Str) {l_, r_, r_ - l_ + 1, i};
    val[i] = subs[i].lth;
  }
  nb = read();
  for (int i = 1; i <= nb; ++ i) {
    int l_ = read(), r_ = read();
    subs[na + i] = (Str) {l_, r_, r_ - l_ + 1, na + i};
  }
  m = read();
  for (int i = 1; i <= m; ++ i) {
    int u_ = read(), v_ = read();
    Add(u_, v_ + na); //Add(read(), read()+na) 会倒着读
  }
  node_num = na + nb;
}
bool Check(int x_, int y_, int lth_) {
  return SA::Lcp(x_, y_) >= lth_;
}
void AddEdgeB(int id_, int now_) {
  int pos = SA::rk[subs[id_].l], l_ = pos, r_ = pos; //l_，r_ 初始值
  for (int l = 1, r = pos - 1; l <= r; ) {
    int mid = (l + r) >> 1;
    if (Check(mid, pos, subs[id_].lth)) {
      r = mid - 1;
      l_ = mid;
    } else {
      l = mid + 1;
    }
  }
  for (int l = pos + 1, r = n; l <= r; ) {
    int mid = (l + r) >> 1;
    if (Check(pos, mid, subs[id_].lth)) {
      l = mid + 1;
      r_ = mid;
    } else {
      r = mid - 1; 
    }
  }
  Hjt::AddEdge(Hjt::root[now_], 1, n, l_, r_, subs[id_].id);
}
void Build() {
  node_num = na + nb;
  std::sort(subs + 1, subs + na + nb + 1, cmp);
  for (int now = 0, i = 1; i <= na + nb; ++ i) {
    if (subs[i].id > na) {
      AddEdgeB(i, now);
      continue;
    }
    ++ now;
    Hjt::Insert(Hjt::root[now], Hjt::root[now - 1], 1, n, SA::rk[subs[i].l], 
                subs[i].id);
  }
}
void TopSort() {
  std::queue <int> q;
  for (int i = 1; i <= node_num; ++ i) {
    if (!into[i]) {
      f[i] = val[i];
      q.push(i);
    }
  }
  while (! q.empty()) {
    int u_ = q.front(); q.pop();
    for (int i = head[u_]; i; i = ne[i]) {
      int v_ = v[i];
      Chkmax(f[v_], f[u_] + val[v_]);
      -- into[v_];
      if (!into[v_]) q.push(v_);
    }
  }
  LL ans = 0;
  for (int i = 1; i <= node_num; ++ i) {
    Chkmax(ans, f[i]);
    if (into[i]) {
      printf("-1
");
      return ;
    }
  }
  printf("%lld
", ans);
}
//=============================================================
int main() {
  int T = read();
  while (T --) {
    Init();
    Build();
    TopSort();
  }
  return 0;
}

与后缀树的联系

后缀树相关介绍可以参考：「笔记」后缀树。

SA 可以看作由后缀树的所有叶结点按字典序排列得到的数组。
SA 的 (operatorname{height}) 数组可以看做后缀树中相邻两关键点的 (operatorname{lca})，从这个角度上可以直观地理解 SA 的许多性质。

同时，这也证明了后缀树的适用范围一定比 SA 大。
奈何 SA 配合 (operatorname{height}) 可以实现大部分后缀树的功能，又好写，时间上与后缀树又几乎相同，空间又小，所以我喜欢 SA（

做题的时候可以先直观地从后缀树的角度出发思考，再用 SA 进行实现。
可以避免抽象的思考过程，比较适合我这样的无脑选手（

后缀数组结合虚树的思想可以用于构建后缀树，详见这篇博客：利用后缀数组构造后缀树_AZUI

写在最后

参考资料：

OI-wiki SA
后缀数组详解 - 自为风月马前卒
「后缀排序SA」学习笔记 - Rainy7
后缀数组-许智磊
后缀数组学习笔记 _ Menci's Blog
利用后缀数组构造后缀树_AZUI
P5284 [十二省联考2019]字符串问题 - xht37 的洛谷博客

作者@Luckyblock，转载请声明出处。
【推广】 免费学中医，健康全家人
原文地址：https://www.cnblogs.com/luckyblock/p/14217034.html