算法与数据结构 13 - 字符串算法+ 数论大杂烩

0. 部分运算的规定

若无特殊说明，所有变量均为整数，整数除法自动向下取整。

本文约定 \(\log x=\log_2x\)。

1.1 border 与 KMP 算法

除特殊说明外，各下标从 \(0\) 开始。

对于字符串 \(s\) 和 \(0\le r<|s|\)，若 \(s\) 长度为 \(r\) 的前缀与长度为 \(r\) 的后缀相等，则称 \(s\) 长度为 \(r\) 的前缀是 \(s\) 的 border。例如，abcab 中有 border ab。

为了在 \(O(n)\) 的时间复杂度内计算 \(s\) 每个前缀的 border 长度，考虑设 \(p_i\) 表示 \(s\) 长度为 \(i\) 的前缀的 border 长度。假设在计算 \(p_i\) 时 \(p_{0\cdots i-1}\) 都计算过了，因此有 \(s[0\cdots p_{i}-1]=s[i-p_{i}+1\cdots i]\)。如果有 \(s_i=s_{p_{i-1}}\)，则 \(s[0\cdots p_i]=s[i-p_i\cdots i]\)，所以 \(p_i=p_{i-1}+1\)；否则，则考虑 \(s[p_{p_{i-1}-1}]\overset{\text{?}}{=}s[i]\)，因为 \(s[0\cdots p_{p_i-1}-1]=s[p_i-p_{p_i-1}\cdots p_i-1]=s[i-p_{p_i-1}+1\cdots i-1]\)，以此类推。

vector<int> prefix_function(string s) {
    int n = (int)s.length();
    vector<int> pi(n);
    for (int i = 1; i < n; i++) {
        int j = pi[i - 1];
        while (j > 0 && s[i] != s[j]) j = pi[j - 1];
        if (s[i] == s[j]) j++;
        pi[i] = j;
    }
    return pi;
}

由此，可以得到一个 \(O(|s|+|t|)\) 地计算 \(s\) 在 \(t\) 中出现的所有位置的算法：构造字符串 \(s\)+#+\(t\)，跑一遍上述算法，找到所有 \(p_i=|s|\) 的 \(i>|s|\) 即可。另外，这个做法的空间复杂度为 \(O(|s|)\)，因为在计算 \(t\) 部分的 \(p_i\) 时只需要 \(p_{i-1}\) 和 \(p_{0\cdots |s|-1}\) 即可。这个算法被称为 KMP，上文提到的 \(p_i\) 数组也被称为 \(nxt\) 数组或 \(fail\) 数组。如果在下标从 \(1\) 开始的 \(s\) 上求 \(nxt\)，部分细节需要调整。下面的代码是按下标从 \(1\) 开始写的，其中 \(a\) 是匹配串，\(b\) 是模式串。

p[1] = 0;
for (int i = 2, j = 0; i <= m; i++) {
	while (j && b[j + 1] != b[i]) j = p[j];
	if (b[j + 1] == b[i]) j++;
	p[i] = j;
}
for (int i = 1, j = 0; i <= n; i++) {
	while (j && b[j + 1] != a[i]) j = p[j];
	if (b[j + 1] == a[i]) j++;
	if (j == m) j = p[j], cout << i - m + 1 << '\n';
}

例题：P2375 [NOI2014] 动物园

设 \(num_i\) 为 \(s\) 长为 \(i\) 中长度不大于 \(\frac i2\) 的不同公共前后缀的数量，求 \(\Pi (num_i+1)\mod (10^9+7)\)。

根据上面的 KMP 过程，容易想到：设 \(s\) 的下标从 \(1\) 开始，则 \(s[1\cdots nxt[i]]\)、\(s[1\cdots nxt[nxt[i]]]\)、\(s[1\cdots nxt[nxt[nxt[i]]]]\)……都是 \(s[1\cdots i]\) 的公共前后缀。因此有个 \(O(n^2)\) 的思路就是求出 \(nxt_i\) 之后暴力跳转计数。

考虑优化，发现可以在 KMP 过程中先求出不考虑长度问题的公共前后缀数量 \(g_i\)，方法是用类似 DP 的 \(g_i=g_{nxt[i]}+1\) 递推，于是有 \(num_i=g_{nxt[nxt[...nxt[i]]]}\)，仍然会超时。考虑在求出 \(num_{i-1}\) 时保留 \(j=nxt[nxt[...nxt[i]]]\le frac{i-1}2\)，这样就不需要从头开始枚举了，时间复杂度平摊 \(O(n)\)。

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define gc getchar
template<typename T> void read(T &x) {
    char c = gc(); int f = 1; x = 0;
    while (!isdigit(c)) f = (c == '-' ? -f : f), c = gc();
    while (isdigit(c)) x = (x << 3) + (x << 1) + c - '0', c = gc();
    x *= f;
}
template<typename T, typename ... Args> void read(T &x, Args& ... args) {read(x); read(args...);}
#define rep(i, s, t) for(int i=(s);i<=(t);i++)
#define per(i, s, t) for(int i=(s);i>=(t);i--)
#define ll long long
#define ull unsigned long long
#define pii pair<int,int>
#define pll pair<long long, long long>
#define lll __int128
template<typename T> void chkmax(T &x, const T &y) {x = max(x, y);}
template<typename T> void chkmin(T &x, const T &y) {x = min(x, y);}
const int N = 1e6 + 5;
const int mod = 1e9 + 7;
int n, nxt[N], g[N];
char s[N];
void pd() {
    scanf("%s", s + 1);
    n = strlen(s + 1);
    int j = 0; nxt[1] = 0; g[1] = 1;
    rep(i, 2, n) {
        while (j && s[j + 1] != s[i]) j = nxt[j];
        if (s[j + 1] == s[i]) j++;
        nxt[i] = j; g[i] = g[j] + 1;
    }
    ll ans = 1; j = 0;
    rep(i, 2, n) {
        while (j && s[j + 1] != s[i]) j = nxt[j];
        if (s[j + 1] == s[i]) j++;
        while (j * 2 > i) j = nxt[j];
        ans = ans * (g[j] + 1) % mod;
    }
    printf("%lld\n", ans);
}
int main() {
    int _;
    read(_);
    while (_--) pd();
}

1.2 Z 函数（扩展 KMP）

定义 \(z[i]\) 为 \(s\) 和 \(s[i\cdots n-1]\) 的最长公共前缀（LCP）长度，此之谓 \(s\) 的 Z 函数或扩展 KMP。特别地，\(z[0]=0\)。

同样的，在考虑求 \(z[i]\) 时，假设已经求完了 \(z[0\cdots i-1]\)。称 \(s[i\cdots i+z[i]-1]\) 为 \(i\) 的匹配段（Z-box）。计算时维护右端点最靠右的匹配段，设为 \([l\cdots r]\)。有 \(s[0\cdots r-l]=s[l\cdots r]\)。考虑如果 \(i\le r\)，则有 \(s[i\cdots r]=s[i-l\cdots r-l]\)，因此 \(z[i]\ge \min(z[i-l],r-i+1)\)。如果 \(z[i-l]<r-i+1\)，则说明 \(z[i]=z[i-l]\)；反之，则可以从 \(r-i+1\) 暴力扩展 \(z[i]\)，然后用 \([i\cdots i+z[i]-1]\) 更新 \([l,r]\)。如果 \(i>r\) 则直接暴力扩展，再更新 \([l\cdots r]\) 即可。

vector<int> z_function(string s) {
    int n = (int)s.length();
    vector<int> z(n);
    for (int i = 1, l = 0, r = 0; i < n; ++i) {
        if (i <= r && z[i - l] < r - i + 1) {
            z[i] = z[i - l];
        } else {
            z[i] = max(0, r - i + 1);
            while (i + z[i] < n && s[z[i]] == s[i + z[i]]) ++z[i];
        }
        if (i + z[i] - 1 > r) l = i, r = i + z[i] - 1;
    }
    return z;
}

例题：P13790 「CZOI-R6」Border
过几天再做。

1.3 AC 自动机

AC 自动机是一种建立在 Trie 上的、用于解决多模式匹配任务的算法。其核心思想是在 Trie 上建立失配指针（Fail），当匹配失败时沿着 Fail 移动到当前状态的最长后缀状态以继续匹配。例如，his 的失配指针可以指向 is，如果这是另一个模式串的前缀。

构造失配指针的思路也比较好理解：假设现在 Trie 已经建好，当前在节点 \(u\)，父节点是 \(p\)，\(p\) 指向 \(u\) 的字符是 \(c\)。如果 \(fail(p)\) 有边指向 \(c\)，则让 \(u\) 的 Fail 指向 \(fail(p)\) 那条 \(c\) 边指向的节点，相当于在 \(fail(p)\) 的后面添加一个字符 \(c\)；反之，则检查 \(fail(fail(p))\) 有没有 \(c\) 边，重复判断，直到找到或到达根节点。如果始终找不到，则让 \(u\) 的 Fail 指向根节点，表示不存在这样的最长后缀。

当然，真正建立 AC 自动机时我们不可能像这样一直递归，否则遇到 aaa...aaab 这样的串就完犊子了。因此，需要使用一些方法对这个“跳转”过程进行压缩，在 \(O(1)\) 的时间复杂度内解决问题。注意到反复跳转发生在一直找不到 \(fail(fail(...fail(u)))\rightarrow ^c\) 的边，因此考虑利用不存在的 \(\rightarrow ^c\)：在为 \(u\) 建立 Fail 指针时，假如发现该节点没有指向 \(c\) 的边，但 \(c\) 出现在字符集里，则考虑让 \(u\rightarrow^c\) 指向 \(fail(u)\rightarrow^c\)，这样就对跳转完成了压缩，后面为 \(p\rightarrow^c u\) 添加 Fail 时直接让 \(u\) 的 Fail 指向 \(fail(p)\rightarrow^c\) 即可。

例题：P14363 [CSP-S 2025] 谐音替换

把 \(s_1,s_2\) 转化为 ABC、ADC，其中 AC 分别是最长公共前后缀，于是把它处理为 A#BD#C，其中 # 表示特殊字符，\(t_1,t_2\) 照此处理，则替换合法当且仅当 \(s\) 串是 \(t\) 串的子串。于是对 \(s\) 串建立 AC 自动机即可。代码稍后补。

1.4 Manacher

Manacher 用于解决求 \(s\) 中第 \(i\) 个字符的最大回文半径类问题。考虑维护最靠右的已求出的回文子串 \(s[l\cdots r]\)，现在要求 \(i\) 为中心的回文半径长度 \(d[i]\)。假如 \(i>r\)，则暴力找 \(d[i]\)，然后更新 \(s[l\cdots r]\)；否则，因为有 \(s[i\cdots r]=s[l\cdots l+r-i]\)，置 \(d[i]\) \(d[l+r-i]\)，若 \(i+d[i]-1\ge r\) 则暴力扩展，否则就这样。

2.1 中国剩余定理（CRT）

中国剩余定理的作用是求出 \(x\)，使得对于给出的 \(a_i,n_i\)，满足 \(x\equiv a_i(\mod n_i)\)。保证 \(n_i\) 两两互质。

过程：

计算所有模数的乘积 \(n\)；
对于第 \(i\) 组限制，
1. 计算 \(m_i=\frac{n}{n_i}\)，
2. 计算 \(m_i\) 在模 \(n_i\) 意义下的逆元 \(m^{-1}_i\)，
3. 计算 \(c_i=m_im^{-1}_i\)；
\(x\equiv \sum a_ic_i\mod n\)。

这样做是对的，证明如下：

当 \(x\) 对 \(n_i\) 取模时，\(c_i=m_im^{-1}_i\equiv1\mod n_i\)，\(c_j=m_j\times\text{something}=n_i\times\text{something}\equiv0\mod n_i(j\not=i)\)，故 \(x\equiv 0+0+0+\cdots+0+a_i\equiv a_i\mod n_i\)。

但当 \(n_i\) 之间不互质时，那可就麻烦了，因为上面那种 \(c_j\equiv[i=j]\) 的奇妙性质顿时不复存在了。这时需要扩展中国剩余定理（exCRT）。

exCRT 的思想是把这些形如 \(x\equiv a_i\mod n_i\) 的式子两两合并，最终只剩一个，或者报告无解。因此，只要处理好两个式子合并就能合并所有式子。

假设有两个式子：\(x\equiv a_1\mod n_1,x\equiv a_2\mod n_2\)。把式子变形为 \(k_1\times n_1+a_1=x=k_2\times n_2+a_2\)，因此有 \(a_1-a_2=k_2\times n_2-k_1\times n_1\)。根据裴蜀定理，如果 \(\gcd(n_1,n_2)|a_1-a_2\) 则必然有解。

裴蜀定理：设 \(a,b\) 是不全为 \(0\) 的整数，则对于任意整数 \(x,y\) 都有 \(\gcd(a,b)|ax+by\)，而且存在整数 \(x,y\) 使得 \(gcd(a,b)=ax+by\)。

那么，再求出 \(k_1,k_2\) 后，就可以得到 \(x\equiv k_1\times n_1+a_1\equiv k_2\times n_2+a_2\mod \operatorname{lcm}(n_1,n_2)\)。所以怎么求 \(k_1,k_2\) 呢？答案是扩展欧几里得算法（exGCD）。

int exgcd(int a, int b, int &x, int &y) {
    if (!b) {
        x = 1;
        y = 0;
        return a;
    }
    int d = Exgcd(b, a % b, x, y);
    int t = x;
    x = y;
    y = t - (a / b) * y;
    return d;
}

2.2 线性筛与积性函数

若想知道 \(\le n\) 的质数有多少个，则需要使用筛法。如果对每个数跑一边质因数分解的话那就太慢了，所以考虑使用线性筛。

线性筛又叫欧拉筛，是一种 \(O(n)\) 的筛法。对于求质数，一种朴素的想法是对每个访问到的数 \(x\)，枚举 \(2x,3x,4x,\cdots\)，然后把这些数全都标记为合数。可以证明这样的时间复杂度是 \(O(n\log\log n)\)，与我们想的有一些差距，这是因为同一个合数被标记的多次。其实，只需要枚举 \(x\) 的质数倍数，然后在枚举质数时顺便判断这个质数是不是 \(x\) 的因数即可，因为如果这个质数是 \(x\) 的倍数则 \(x\) 的更高质数倍一定被这个质数标记过，因此直接退出循环即可。

for (int i = 2; i <= n; i++) {
        if (prime[i]) {
            list[++list[0]] = i;
        }
        for (int j = 1; j <= list[0]; j++) {
            if (i * list[j] > n) break;
            prime[i * list[j]] = false;
            if (i % list[j] == 0) break;
        }
    }

可以发现，对于 \(i\times list[j]\)，\(list[j]\) 是它的最小质因子。这对后面的积性函数有作用。

什么是积性函数呢？就是对于互质的 \(a,b\)，有 \(f(ab)=f(a)f(b)\) 的函数。假如对于一个积性函数 \(f\) 和任意一个质数 \(p\)，能快速（通常 \(O(1)\)）地求出 \(f(p^k)\)（\(k\) 是任意正整数），则可以通过线性筛在 \(O(n)\) 内求出 \(f(i)\)。

注意到，当进行线性筛时，若 \(i\) 不是 \(list[j]\) 的倍数，则 \(f(i\times list[j])=f(i)\times f(list[j])\)；反之，则可以先从 \(i\) 中去处所有 \(list[j]\)，得到 \(a=\frac{i}{list[j]^k}\) 和 \(b=list[j]^k\)，然后有 \(f(i\times list[j])=f(ab)=f(a)f(b)\)。但这么做是 \(O(n\log n)\) 的~~（尽管快得飞起）~~，因此考虑优化。考虑在线性筛的过程中记录下每个数最大的「最小质因子的次幂」因数，记为 \(low[i]\)，然后就有 \(f(i\times list[j])=f(\frac{i\times list[j]}{low[i\times list[j]]})\times f(low[i\times list[j]))\)。

sig[1] = 1; low[1] = 1;
for (int i = 2; i <= n; i++) {
    if (prime[i]) {
        list[++list[0]] = i;
        low[i] = i;
        sig[i] = ... ;                                        // 求 f(p) 的值
    }
    for (int j = 1; j <= list[0]; j++) {
        if (i * list[j] > n) break;
        prime[i * list[j]] = false;
        if (i % list[j] == 0) {
            low[i * list[j]] = low[i] * list[j];
            int b = i * list[j] / low[i * list[j]], a = low[i * list[j]];
            if (b == 1) sig[i * list[j]] = ... ;              // 求 f(p^k) 的值
            else sig[i * list[j]] = sig[a] * sig[b];          // 将 f(i*list[j]) 转化为 f(a*b)=f(a)*f(b)
            break;
        } else {
            low[i * list[j]] = list[j];
            sig[i * list[j]] = sig[i] * sig[list[j]];         // 按一般情况处理
        }
    }
}

常见的积性函数有不少，如：

欧拉函数（不大于 \(n\) 的正整数中与 \(n\) 互质的数的个数）；
莫比乌斯函数（\(\mu(x)=1\text{ when }x=1,0 \text{ when } x=a\times a\times b,(-1)^k \text{ when } x=a_1\times a_2\times \cdots a_k \text{ and they are all prime}\)）。

2.3 容斥原理

容斥原理的目的是用一堆集合的交集求一堆集合的并集，有算式

\[|\bigcup_{i=1}^nS_i|=\sum_{k=1}^n(-1)^{k-1}\sum_{a_i<a_{i+1}}|\bigcap_{i=1}^kS_{a_i}| \]

容斥原理一般有两个应用：1. 正难则反；2. 条件计数。

2.4 卢卡斯定理

卢卡斯定理是：对于质数 \(p\)，有 \(\binom{n}{k}=\binom{n/p}{k/p}\binom{n\%p}{k\%p}\%p\)，其中 \(/\) 表示整数除法，\(\%\) 表示取模运算。卢卡斯定理通常要求模数为质数且在 \(10^6\) 附近，这是因为预处理阶乘逆元是 \(O(p)\) 的。一种常见的实现：

long long C(long long n, long long k) {
    if (n < k) return 0;
    return fac[n] * invf[k] % p * invf[n - k] % p;
}
long long Lucas(long long n, long long k) {
  if (k == 0) return 1;
  return (C(n % p, k % p) * Lucas(n / p, k / p)) % p;
}

2.5 卡特兰数与斯特林数

2.6 矩阵运算

2.7 高斯消元