OI模板梳理

明天自招考试!
花点时间把一些老博客上的一些板子转移过来,顺便也复习一下。

1 图论

最大流

此处采用多路增广的 Dinic 算法,并加入了当前弧优化。注意,该算法的复杂度为 $O(n^2m)$,对于存在一个边的容量均是 1 的割边集的图,该算法的复杂度为 $O(n\sqrt{m})$。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
namespace DINIC {
class fs{
public: int to, nxt, cap, fl;
void init(int t, int n, int c, int f) {to = t, nxt = n, cap = c, fl = f;}
} E[CN * 10]; int hd[CN], cur[CN], ecnt;
void clr() {memset(hd, 0, sizeof(hd)), ecnt = 1;}
void adde(int x, int y, int z) {
E[++ecnt].init(y, hd[x], z, 0), hd[x] = ecnt;
E[++ecnt].init(x, hd[y], 0, 0), hd[y] = ecnt;
}
int dep[CN]; queue<int> Q;
bool bd(int s, int t){
memset(dep, 0, sizeof(dep)), dep[s] = 1, Q.push(s);
while(!Q.empty()){
int u = Q.front(); Q.pop();
for(int k = hd[u]; k; k = E[k].nxt){
fs e = E[k];
if(e.cap > e.fl && !dep[e.to]) dep[e.to] = dep[u] + 1, Q.push(e.to);
}
}
return dep[t];
}
int aug(int u, int t, int rst){
if(u == t) return rst;
int usd = 0;
for(int k = cur[u]; k; k = E[k].nxt){
fs e = E[k]; cur[u] = k;
if(e.cap > e.fl && dep[e.to] == dep[u] + 1){
int add = aug(e.to, t, min(rst - usd, e.cap - e.fl));
if(add){
E[k].fl += add, E[k ^ 1].fl -= add, usd += add;
if(usd == rst) return usd;
}
}
}
return usd;
}
int mf(int s, int t){
int res = 0;
while(bd(s, t)) memcpy(cur, hd, sizeof(hd)), res += aug(s, t, INF);
return res;
}
}

最小费用流

此处采用在稠密图上效率更高的 ZKW 费用流,该算法的复杂度为 $O(n^2m)$。关于最短路增广实现的费用流,请参考站内相关博客。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
namespace ZKW{
class fs {
public: int to, nxt, cap, fl, w;
void init(int t, int n, int c, int f, int ww) {
to = t, nxt = n, cap = c, fl = f, w = ww;
}
} E[CN << 2];
int hd[CN], ecnt = 1;
void adde(int x, int y, int z, int c){
E[++ecnt].init(y, hd[x], z, 0, c), hd[x] = ecnt;
E[++ecnt].init(x, hd[y], 0, 0, -c), hd[y] = ecnt;
}
int cst, mf, d[CN]; bool vis[CN];
bool bd(int t){
if(vis[t]) return 1;
int del = INF;
for(int i = 2; i <= ecnt; i++){
int u = E[i ^ 1].to, v = E[i].to;
if(vis[u] && !vis[v] && E[i].cap > E[i].fl)
del = min(del, d[u] + E[i].w - d[v]);
}
if(del == INF) return 0;
for(int i = 1; i <= n; i++) if(vis[i]) d[i] -= del;
return 1;
}
int aug(int u, int t, int rst){
vis[u] = 1;
if(u == t) return rst;
for(int k = hd[u]; k; k = E[k].nxt){
fs e = E[k];
if(e.cap > e.fl && !vis[e.to] && d[u] + e.w == d[e.to]){
int add = aug(e.to, t, min(rst, e.cap - e.fl));
if(add) return E[k].fl += add, E[k ^ 1].fl -= add, add;
}
}
return 0;
}
void mcmf(int s, int t){
memset(d, 0, sizeof(d)), cst = mf = 0;
do{
memset(vis, 0, sizeof(vis));
int add = aug(s, t, INF);
mf += add, cst -= add * d[s];
}
while(bd(t));
printf("%d %d\n", mf, cst);
}
}

其它

以下内容请在本站搜索相应文章。
  • 树链剖分 \ LCA
  • 树的重心 \ 树的直径
  • 强连通分量(缩点) \ 双连通分量
  • 倍增LCA
  • 二分图匹配(匈牙利)
  • k 短路问题(A* 寻路,请参阅这篇文章

2 数论

素数筛

Eratosthenes 筛法,此筛法非严格的线性。有关线性筛,请见另一篇博文

1
2
3
4
5
6
7
8
9
10
11
const int CN = 1e7 + 7;
bool isp[CN]; int prime[CN];
void GetPrime(int n){
memset(isp, 1, sizeof(isp));
for(int i = 2; i <= n; i++){
if(isp[i]){
prime[++prime[0]] = i;
for(int k = 2; i * k <= n; k++) isp[i * k] = 0;
}
}
}

分解质因数

试除法,复杂度$O(\sqrt{x})$。
注:一个数的质因数数量是$\log$级别。

1
2
3
4
5
6
7
8
9
10
11
#define LL long long
const int CN = 101;
int cnt = 0; LL p[CN], a[CN];
void Div(LL x){
for(LL i = 2; i * i <= x; i++){
if(x % i) continue;
p[++cnt] = i;
while(!(x % i)) {a[cnt]++; x /= i;}
}
if(x > 1) p[++cnt] = x,a[cnt] = 1;
}

高斯消元

1
2
3
4
5
6
7
8
9
10
11
12
13
14
#define DB double
bool equ(int a[][CN]){
for(int i = 1; i <= n; i++){
int p = i;
while(!a[p][i] && p <= n) p++;
if(!a[p][i]) return puts("No Solution"), 0;
swap(a[p], a[i]);
for(int j = 1; j <= n; j++){
if(i == j) continue;
DB t = a[j][i] / a[i][i];
for(int k = i; k <= n + 1; k++) a[j][k] -= t * a[i][k];
}
}
}

扩展欧几里得算法

扩展欧几里得算法用于求解形如 $ax+by=c$ 的二元不定方程的整数解。根据裴蜀定理,该不定方程有解当且仅当 $(a,b)|c$。容易发现,若 $(x,y)$ 是解,则 $(x+kb, y-ka), k\in \mathbb{Z}$ 也是解。

1
2
3
4
5
6
7
8
9
int gcd(int a, int b) {return b ? gcd(b, a % b) : a;}
void exgcd(int a, int &x, int b, int &y){
if(!b) return (void)(x = 1, y = 0);
exgcd(b, x, a % b, y); int t = x; x = y, y = t - (a / b) * y;
}
bool solve(int a, int &x, int b, int &y, int c){
int g = gcd(a, b); if(c % g) return false;
return exgcd(a, x, b, y), c /= g, x *= c, y *= c, true;
}

类欧几里得算法

类欧几里得算法仅适用于处理斜率和截距非负的线段。当斜率为负时,需要通过对称变换使得斜率为正;当截距为负时,需要通过平移坐标轴使截距为正。

求解:

$$f(n,A,B,C)=\sum_{i=0}^n\left\lfloor \frac{Ai+B}{C} \right\rfloor$$

其中满足 $A,B,C \ge 0$。它的几何意义是:在第一象限内,一条斜率和截距非负的线段下方的整点的数量。

实数版:

1
2
3
4
5
6
7
8
#define LL long long
LL s2(LL n) {return n * (n + 1) / 2;}
LL f(LL n, LL a, LL b, LL c){
if(!a) return (b / c) * (n + 1);
if(a >= c || b >= c) return s2(n) * (a / c) + (n + 1) * (b / c) + f(n, a % c, b % c, c);
LL m = (a * n + b) / c;
return n * m - f(m - 1, c, c - b - 1, a);
}

取模版:

1
2
3
4
5
6
7
8
9
#define LL long long
const int P = 998244353;
int s2(LL n) {return (n * (n + 1) / 2) % P;}
int f(LL n, LL a, LL b, LL c){
if(!a) return (b / c) * (n + 1) % P;
if(a >= c || b >= c) return (s2(n) * (a / c) % P + (n + 1) * (b / c) % P + f(n, a % c, b % c, c)) % P;
LL m = (a * n + b) / c;
return (n * m % P - f(m - 1, c, c - b - 1, a) + P) % P;
}

BSGS

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
int mul(int x, int y, int P) {return (long long)x * y % P;}
int qp(int a, int b, int P){
int r = 1;
for(; b; a = mul(a, a, P), b >>= 1) if(b & 1) r = mul(r, a, P);
return r;
}
int B; map<int, int> vis;
void bd(int a, int P){
B = ceil(sqrt(P)), vis.clear(); int t = qp(a, B, P);
for(int i = 1, p = t; i <= B; i++, p = mul(p, t, P))
if(!vis.count(p)) vis[p] = i * B;
}
int qu(int a, int b, int P){
int ans = 2e9;
for(int i = 0, p = 1; i < B; i++, p = mul(p, a, P)){
int cur = mul(b, p, P);
if(vis.count(cur)) ans = min(ans, vis[cur] - i);
}
return ans < int(2e9) ? ans : -1;
}

Cipolla’s Algorithm

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
class COMP {public: int x, y;}; int w;
COMP mk(int a, int b) {COMP o; o.x = a, o.y = b; return o;}
COMP mul(COMP a, COMP b, int p){
COMP r;
r.x = (1ll * a.x * b.x % p + 1ll * a.y * b.y % p * w % p) % p;
r.y = (1ll * a.x * b.y % p + 1ll * a.y * b.x % p) % p;
return r;
}
int qp(int a, int b, int p){
int r = 1;
for(; b; b >>= 1, a = 1ll * a * a % p) if(b & 1) r = 1ll * r * a % p;
return r;
}
COMP qp(COMP a, int b, int p){
COMP r = mk(1, 0);
for(; b; b >>= 1, a = mul(a, a, p)) if(b & 1) r = mul(r, a, p);
return r;
}
bool ck(int n, int p) {return qp(n, (p - 1) / 2, p) == 1;}
int sqrt(int n, int p){
n %= p;
if(p == 2) return n; if(!n) return 0;
if(!ck(n, p)) return -1;
int a = rand() % p;
while(!a || ck((1ll * a * a % p - n + p) % p, p)) a = rand() % p;
COMP x = mk(a, 1); w = (1ll * a * a % p - n + p) % p;
return qp(x, (p + 1) / 2, p).x;
}

其它

以下内容请在本站搜索相应文章。
  • 矩阵快速幂(矩阵加速递推)
  • lucas(组合数取模)
  • gcd(最大公约数) \ exgcd(关于其求解不定方程的模板,请参阅CRT&exCRT
  • 逆元
  • CRT & exCRT(同余方程组)

3 多项式

有关多项式的总结请参见多项式学习笔记

FFT / 快速傅里叶变换

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
class COMP{
public: LDB x, y;
COMP operator + (const COMP &o) const{
COMP r = *this; r.x += o.x, r.y += o.y;
return r;
}
COMP operator - (const COMP &o) const{
COMP r = *this; r.x -= o.x, r.y -= o.y;
return r;
}
COMP operator * (const COMP &o) const{
COMP r; r.x = x * o.x - y * o.y, r.y = x * o.y + y * o.x;
return r;
}
} ;
COMP mk(LDB a, LDB b) {COMP o; o.x = a, o.y = b; return o;}
int rev[CN << 2];
void cg(COMP t[], int n){
for(int i = 0; i < n; i++) rev[i] = (rev[i >> 1] >> 1) + (i & 1) * (n >> 1);
for(int i = 0; i < n; i++) if(i < rev[i]) swap(t[i], t[rev[i]]);
}
void fft(COMP t[], int n, int tp){
cg(t, n);
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1; COMP gn = mk(cos(2 * PI / (LDB)w), sin(2 * tp * PI / (LDB)w));
for(int i = 0; i < n; i += w){
COMP g = mk(1, 0);
for(int j = i; j < i + l; j++){
COMP u = t[j], v = t[j + l] * g;
t[j] = u + v, t[j + l] = u - v, g = g * gn;
}
}
}
if(tp ^ 1){
for(int i = 0; i < n; i++) t[i].x /= (LDB)n, t[i].y /= (LDB)n;
}
}

NTT / 快速数论变换

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
int qp(LL a, LL b){
LL r = 1;
for(; b; b >>= 1, a = a * a % P) if(b & 1) r = r * a % P;
return r;
}
void cg(int t[], int n){
for(int i = 0; i < n; i++) rev[i] = (rev[i >> 1] >> 1) + (i & 1) * (n >> 1);
for(int i = 0; i < n; i++) if(i < rev[i]) swap(t[i], t[rev[i]]);
}
void ntt(int t[], int n, int tp){
cg(t, n); int G = 3, iG = qp(G, P - 2);
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1; LL gn = qp(tp ? G : iG, (P - 1) / w), g = 1;
for(int i = 0; i < n; g = 1, i += w) for(int j = i; j < i + l; j++){
LL u = t[j], v = g * t[j + l];
t[j] = (u + v) % P, t[j + l] = (u - v) % P, g = g * gn % P;
}
}
for(int i = 0; i < n; i++) t[i] = (t[i] + P) % P;
if(!tp) {LL iv = qp(n, P - 2); for(int i = 0; i < n; i++) t[i] = iv * t[i] % P;}
}

FWT / 快速沃尔什变换

按位异或卷积

1
2
3
4
5
6
7
8
9
10
void fwt(int t[], int n, int tp){
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1;
for(int i = 0; i < n; i += w) for(int j = i; j < i + l; j++){
int u = t[j], v = t[j + l];
t[j] = add(u, v), t[j + l] = add(u, P - v);
if(!tp) t[j] = 1ll * t[j] * i2 % P, t[j + l] = 1ll * t[j + l] * i2 % P;
}
}
}

按位或卷积

1
2
3
4
5
6
7
8
9
void fwt(int t[], int n, int tp){
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1;
for(int i = 0; i < n; i += w) for(int j = i; j < i + l; j++){
if(tp) t[j + l] = add(t[j + l], t[j]);
else t[j + l] = add(t[j + l], P - t[j]);
}
}
}

按位与卷积

1
2
3
4
5
6
7
8
9
void fwt(int t[], int n, int tp){
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1;
for(int i = 0; i < n; i += w) for(int j = i; j < i + l; j++){
if(tp) t[j] = add(t[j], t[j + l]);
else t[j] = add(t[j], P - t[j + l]);
}
}
}

求逆 / 指对运算

相当于多项式全家桶了。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
int rev[CN << 2];
void cg(int t[], int n){
for(int i = 0; i < n; i++) rev[i] = (rev[i >> 1] >> 1) + (i & 1) * (n >> 1);
for(int i = 0; i < n; i++) if(i < rev[i]) swap(t[i], t[rev[i]]);
}
void ntt(int t[], int n, int tp){
cg(t, n); int og = 3, ig = invx(og);
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1, gn = qp(tp ? og : ig, (P - 1) / w);
for(int i = 0; i < n; i += w){
int g = 1;
for(int j = i; j < i + l; j++){
int u = t[j], v = 1ll * t[j + l] * g % P;
t[j] = add(u, v), t[j + l] = add(u, P - v), g = 1ll * g * gn % P;
}
}
}
if(!tp){
ig = invx(n);
for(int i = 0; i < n; i++) t[i] = 1ll * t[i] * ig % P;
}
}
int c[CN << 2];
void inv(int a[], int b[], int n){
for(int i = 0; i < (n << 1); i++) b[i] = c[i] = 0;
b[0] = invx(a[0]);
for(int w = 2; w < (n << 1); w <<= 1){
for(int i = 0; i < w; i++) c[i] = a[i];
ntt(b, w << 1, 1), ntt(c, w << 1, 1);
for(int i = 0; i < (w << 1); i++) b[i] = 1ll * b[i] * add(2, P - (1ll * c[i] * b[i] % P)) % P;
ntt(b, w << 1, 0);
for(int i = w; i < (w << 1); i++) b[i] = 0;
}
for(int i = n; i < (n << 1); i++) b[i] = 0;
}
int ia[CN << 2];
void ln(int a[], int b[], int n){
int N = 1; while(N < (n << 1)) N <<= 1;
for(int i = 0; i < N; i++) ia[i] = b[i] = 0;
inv(a, ia, n);
for(int i = 0; i < N; i++) c[i] = 0;
for(int i = 0; i < n - 1; i++) c[i] = 1ll * (i + 1) * a[i + 1] % P;
ntt(ia, N, 1), ntt(c, N, 1);
for(int i = 0; i < N; i++) c[i] = 1ll * c[i] * ia[i] % P;
ntt(c, N, 0);
for(int i = 1; i < n; i++) b[i] = 1ll * c[i - 1] * invx(i) % P;
b[0] = 0; for(int i = n; i < (n << 1); i++) b[i] = 0;
}
int lnb[CN << 2];
void exp(int a[], int b[], int n){
for(int i = 0; i < (n << 1); i++) b[i] = lnb[i] = 0; b[0] = 1;
for(int w = 2; w < (n << 1); w <<= 1){
ln(b, lnb, w);
for(int i = 0; i < w; i++) lnb[i] = add(a[i], P - lnb[i]);
lnb[0] = add(lnb[0], 1);
ntt(b, w << 1, 1), ntt(lnb, w << 1, 1);
for(int i = 0; i < (w << 1); i++) b[i] = 1ll * b[i] * lnb[i] % P;
ntt(b, w << 1, 0);
for(int i = w; i < (w << 1); i++) b[i] = 0;
}
for(int i = n; i < (n << 1); i++) b[i] = 0;
}

任意模数卷积

拆系数fft,细节参考毛啸(myy)的论文。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
#define LDB long double
class COMP{
public: LDB x, y;
COMP operator + (const COMP &o) const{
COMP r = *this; r.x += o.x, r.y += o.y;
return r;
}
COMP operator - (const COMP &o) const{
COMP r = *this; r.x -= o.x, r.y -= o.y;
return r;
}
COMP operator * (const COMP &o) const{
COMP r; r.x = x * o.x - y * o.y, r.y = x * o.y + y * o.x;
return r;
}
} ;
COMP mk(LDB a, LDB b) {COMP o; o.x = a, o.y = b; return o;}
COMP conj(COMP o) {o.y = -o.y; return o;}
int rev[CN << 2];
void cg(COMP t[], int n){
for(int i = 0; i < n; i++) rev[i] = (rev[i >> 1] >> 1) + (i & 1) * (n >> 1);
for(int i = 0; i < n; i++) if(i < rev[i]) swap(t[i], t[rev[i]]);
}
void fft(COMP t[], int n, int tp){
cg(t, n);
for(int w = 2; w <= n; w <<= 1){
int l = w >> 1; COMP gn = mk(cos(2 * PI / (LDB)w), sin(2 * tp * PI / (LDB)w));
for(int i = 0; i < n; i += w){
COMP g = mk(1, 0);
for(int j = i; j < i + l; j++){
COMP u = t[j], v = t[j + l] * g;
t[j] = u + v, t[j + l] = u - v, g = g * gn;
}
}
}
if(tp ^ 1){
for(int i = 0; i < n; i++) t[i].x /= (LDB)n, t[i].y /= (LDB)n;
}
}
COMP p[CN << 2], q[CN << 2], x[CN << 2], y[CN << 2], z[CN << 2], w[CN << 2];
void conv(int a[], int b[], int n){ // a = a * b
int B = (1 << 15) - 1, N = 1; while(N < (n << 1)) N <<= 1;
for(int i = 0; i < n; i++) p[i] = mk(a[i] >> 15, a[i] & B); // k1 r1
for(int i = 0; i < n; i++) q[i] = mk(b[i] >> 15, b[i] & B); // k2 r2
fft(p, N, 1), fft(q, N, 1);
for(int i = 0; i < N; i++){
int j = (N - 1) & (N - i);
COMP k1, r1, k2, r2;
k1 = (p[i] + conj(p[j])) * mk(0.5, 0);
r1 = (p[i] - conj(p[j])) * mk(0, -0.5);
k2 = (q[i] + conj(q[j])) * mk(0.5, 0);
r2 = (q[i] - conj(q[j])) * mk(0, -0.5);
x[i] = k1 * k2, y[i] = r1 * r2, z[i] = k1 * r2, w[i] = k2 * r1;
}
for(int i = 0; i < N; i++) p[i] = x[i] + y[i] * mk(0, 1);
for(int i = 0; i < N; i++) q[i] = z[i] + w[i] * mk(0, 1);
fft(p, N, -1), fft(q, N, -1);
for(int i = 0; i < N; i++){
LL X = (LL)(0.5 + p[i].x), Y = (LL)(0.5 + p[i].y), Z = (LL)(0.5 + q[i].x), W = (LL)(0.5 + q[i].y);
X = (X % P + P) % P, Y = (Y % P + P) % P, Z = (Z % P + P) % P, W = add((W % P + P) % P, Z);
a[i] = add((X << 30) % P, add((W << 15) % P, Y));
}
}

4 字符串

KMP

有关 KMP 的总结请参见KMP学习笔记

1
2
3
4
5
6
7
8
9
10
int k = 0; nxt[1] = 0, nxt[0] = -1;
for(int i = 2; i <= m; i++){
while(k ^ -1 && t[k + 1] != t[i]) k = nxt[k];
nxt[i] = (k += 1);
}
k = 0;
for(int i = 1; i <= n; i++){
while(k ^ -1 && t[k + 1] != s[i]) k = nxt[k];
if((k += 1) == m) printf("%d", i - m + 1), puts("");
}

ACAM / AC 自动机

有关AC自动机的总结请参见KMP学习笔记

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
const int CN = 1e6 + 6;
class ACAM {
public: int son[CN][26], fail[CN], e[CN], idx; queue<int> Q;
void ins(char *s){ // 建立 tire 结构
int u = 0;
for(int i = 0; s[i]; i++){
if(!son[u][ s[i] - 'a' ]) son[u][ s[i] - 'a' ] = ++idx;
u = son[u][ s[i] - 'a' ];
}
e[u]++;
}
void bd(){ // 建立 fail 指针
for(int i = 0; i < 26; i++) if(son[0][i]) Q.push( son[0][i] );
while(!Q.empty()){
int u = Q.front(); Q.pop();
for(int i = 0; i < 26; i++)
if(son[u][i]) fail[ son[u][i] ] = son[ fail[u] ][i], Q.push(son[u][i]);
else son[u][i] = son[ fail[u] ][i];
}
}
} D;

SA / 后缀数组

后缀数组通过将后缀按字典序排序来获得一些优美的性质。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
int rk[CN << 1], sa[CN], ht[CN];
/* rk[] 将 s[i:] 映射到排序后的排名 rk[i],有类似于离散化的作用 */
namespace SA{
int prk[CN << 1], id[CN], px[CN], cnt[CN];
void sort(char a[], int n){
int m = max(n, 300);
for(int i = 1; i <= n; i++) rk[i] = a[i];
for(int i = 1; i <= n; i++) cnt[ rk[i] ]++;
for(int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for(int i = n; i; i--) sa[ cnt[ rk[i] ]-- ] = i; // id[i] = i
for(int w = 1; w < n; w <<= 1){
memset(cnt, 0, sizeof(cnt));
for(int i = 1; i <= n; i++) id[i] = sa[i];
for(int i = 1; i <= n; i++) cnt[ px[i] = rk[id[i] + w] ]++;
for(int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for(int i = n; i; i--) sa[ cnt[ px[i] ]-- ] = id[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[ px[i] = rk[id[i]] ]++;
for(int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for(int i = n; i; i--) sa[ cnt[ px[i] ]-- ] = id[i];
memcpy(prk, rk, sizeof(rk)), m = 0;
for(int i = 1; i <= n; i++)
if(prk[sa[i]] == prk[sa[i - 1]] && prk[sa[i] + w] == prk[sa[i - 1] + w])
rk[ sa[i] ] = m;
else rk[ sa[i] ] = ++m;
if(m == n) break;
}
for(int p = 0, i = 1; i <= n; i++){
if(p) p--;
while(a[i + p] == a[sa[rk[i] - 1] + p]) p++;
ht[ rk[i] ] = p;
}
}
}

SAM / 后缀自动机

后缀自动机可以将字符串的每一个子串双射在有向单词无环图(DAWG)上,从而获得一系列优美的性质。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
const int CN = 1e6 + 6;
class SAM{
public: int len[CN << 1], nxt[CN << 1], last, sz; int son[CN << 1][26];
SAM() {memset(son, 0, sizeof(son)), len[0] = 0, nxt[0] = -1, sz = 1, last = 0;}
void extend(int c){
int u = sz++, p = last;
last = u, len[u] = len[p] + 1;
while(p != -1 && !son[p][c]) son[p][c] = u, p = nxt[p];
if(p == -1) return (void)(nxt[u] = 0);
int d = son[p][c];
if(len[d] == len[p] + 1) return (void)(nxt[u] = d);
int v = sz++;
len[v] = len[p] + 1, nxt[v] = nxt[d], nxt[u] = nxt[d] = v;
memcpy(son[v], son[d], sizeof(son[d]));
while(p != -1 && son[p][c] == d) son[p][c] = v, p = nxt[p];
}
};

广义后缀自动机

先建立 Tire 树的结构,再在 Tire 树上建立 Parent 树和 DAWG 即可。注意在 SAM 的结构建立完成后,原有的 Tire 结构会被破坏。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
class PAIR {public: int x, y;};
PAIR mp(int x, int y) {PAIR o; o.x = x, o.y = y; return o;}
class SAM {
public: int nxt[CN << 1], son[CN << 1][26], len[CN << 1], idx; queue<PAIR> Q;
SAM() {nxt[0] = -1;}
int et(int p, int c){ // 在 SAM 结构上的 p 节点后扩展出字符 c
int u = son[p][c]; if(len[u]) return u;
len[u] = len[p] + 1, p = nxt[p];
while(p ^ -1 && !son[p][c]) son[p][c] = u, p = nxt[p];
if(p == -1) return u;
int d = son[p][c];
if(len[d] == len[p] + 1) return nxt[u] = d, u;
int v = ++idx;
len[v] = len[p] + 1, nxt[v] = nxt[d], nxt[d] = nxt[u] = v;
for(int i = 0; i < 26; i++) if(len[son[d][i]]) son[v][i] = son[d][i];
while(p ^ -1 && son[p][c] == d) son[p][c] = v, p = nxt[p];
return u;
}
void ins(char ch[]){
int u = 0;
for(int i = 1; ch[i]; i++){
if(!son[u][ch[i] - 'a']) son[u][ch[i] - 'a'] = ++idx;
u = son[u][ch[i] - 'a'];
}
}
void bd(){
for(int i = 0; i < 26; i++) if(son[0][i]) Q.push(mp(0, i));
while(!Q.empty()){
int u, x = Q.front().x, y = Q.front().y; Q.pop();
u = et(x, y);
for(int i = 0; i < 26; i++) if(son[u][i]) Q.push(mp(u, i));
}
}
} ;

Manacher

众所周知,Manacher 是一种优雅的暴力。

1
2
3
4
5
6
7
8
9
10
cin >> (s + 1); n = strlen(s + 1); 
c[0] = '$', c[1] = '#';
for(int i = 1; i <= n; i++) c[i << 1] = s[i], c[i << 1 | 1] = '#';
n = n << 1 | 1, c[n + 1] = '.';
int k, i0 = 1; r[1] = 1;
for(int i = 2; i <= n; i++){
k = min(i0 + r[i0] - i, r[(i0 << 1) - i]);
while(c[i + k] == c[i - k]) k++;
r[i] = k; if(i + r[i] > i0 + r[i0]) i0 = i;
}

PAM / 回文自动机

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
int n; char ch[CN];
class PAM{
public: int len[CN], nxt[CN], son[CN][26], idx, lst;
PAM() {len[0] = nxt[0] = -1, idx = lst = 1;}
void et(int r){
int p = lst;
while(ch[r - len[p] - 1] ^ ch[r]) p = nxt[p];
if(!son[p][ch[r] - 'a']){
son[p][ch[r] - 'a'] = ++idx, len[idx] = len[p] + 2;
int pp = nxt[p];
while(pp ^ -1 && ch[r - len[pp] - 1] ^ ch[r]) pp = nxt[pp];
nxt[idx] = pp ^ -1 ? son[pp][ch[r] - 'a'] : 1;
}
lst = son[p][ch[r] - 'a'];
}
} ;

Lyndon Factorize & Runs

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
const int P = 191011109, B = 39; int pB[CN], ha[CN];
int get(int l, int r) {return add(ha[r], P - (1ll * ha[l - 1] * pB[r - l + 1] % P));}
int LCP(int l1, int r1, int l2, int r2){
if(ch[l1] != ch[l2]) return 0;
int l = 1, r = min(r1 - l1 + 1, r2 - l2 + 1);
while(l < r){
int mid = (l + r + 1) >> 1;
if(get(l1, l1 + mid - 1) == get(l2, l2 + mid - 1)) l = mid;
else r = mid - 1;
}
return l;
}
int LCS(int l1, int r1, int l2, int r2){
if(ch[r1] != ch[r2]) return 0;
int l = 1, r = min(r1 - l1 + 1, r2 - l2 + 1);
while(l < r){
int mid = (l + r + 1) >> 1;
if(get(r1 - mid + 1, r1) == get(r2 - mid + 1, r2)) l = mid;
else r = mid - 1;
}
return l;
}
bool le(int l1, int r1, int l2, int r2){
if(l1 == l2) return r1 < r2;
int l = LCP(l1, r1, l2, r2);
if(l1 + l > r1 || l2 + l > r2) return r1 - l1 < r2 - l2;
return ch[l1 + l] < ch[l2 + l];
}
class PAIR{
public: int l, r, p;
bool operator < (const PAIR &o) const{
return l ^ o.l ? l < o.l : (r ^ o.r ? p < o.p : r < o.r);
}
bool operator == (const PAIR &o) const{
return l == o.l && r == o.r;
}
} ;
PAIR mp(int a, int b, int c) {PAIR o; o.l = a, o.r = b, o.p = c; return o;}
vector<PAIR> runs; int stk[CN], ed[CN], top;
void lyndon(){
top = 0;
for(int i = n; i; i--){
stk[++top] = i;
while(top > 1 && le(i, stk[top], stk[top] + 1, stk[top - 1])) top--;
ed[i] = stk[top];
}
for(int i = 1; i <= n; i++){
int j = ed[i], lcs = LCS(1, i - 1, 1, j), lcp = LCP(i, n, j + 1, n), l, r;
l = i - lcs, r = j + lcp;
if((r - l + 1) / (j - i + 1) > 1) runs.pb(mp(l, r, j - i + 1));
}
}
void bd(){
lyndon();
for(int i = 1; i <= n; i++) ch[i] = 'a' + 'z' - ch[i];
lyndon(), sort(runs.begin(), runs.end());
int len = unique(runs.begin(), runs.end()) - runs.begin();
}
作者

ce-amtic

发布于

2019-03-29

更新于

2021-03-21

许可协议

CC BY-NC-SA 4.0

评论

Your browser is out-of-date!

Update your browser to view this website correctly.&npsb;Update my browser now

×