__gnu_cxx::rope的真面目——浅谈可持久化并查集测试

最近又遇到了可以用可持久化并查集解决的问题,尝试用__gnu_cxx::rope水过,结果被卡成暴力分。鉴于已经有前车之鉴——一次模拟赛也是这样,但是BZOJ上的可持久化并查集却能水过,于是我打算做一些实验。

可持久化数组

由于可持久化并查集的基础是可持久化数组,因此我设计了一道测试用的模板题。

为了简单起见,给定数组 ,进行mm次修改操作,然后qq个询问,求某次修改后某个位置上的数。由于仅用于测试速度,所以数据随机。0a[i]1090\le a[i]\le 10^9

1
2
3
4
5
6
7
8
9
10
11
12
13
5
1 2 3 4 5
4
2 6
3 2
2 3
2 1
5
1 3
4 3
3 5
4 2
2 2
1
2
3
4
5
3
2
5
1
6
id n,m,qn,m,q
1 10410^4
2 3×1043\times 10^4
3 10510^5
4 3×1053\times 10^5
5 10610^6
gen.cpp
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
#include <fstream>
#include <random>
#include <ctime>
using namespace std;
ofstream fout("rope.in");
const int n = 3e5, m = 3e5, q = 3e5, val = 1e9;
minstd_rand gen(time(NULL));
int random(int l, int r)
{
return uniform_int_distribution<>(l, r)(gen);
}
int main()
{
fout << n << endl;
for (int i = 1; i <= n; i++)
fout << random(0, val) << ' ';
fout << endl
<< m << endl;
for (int i = 1; i <= m; i++)
fout << random(1, n) << ' ' << random(0, val) << endl;
fout << q << endl;
for (int i = 1; i <= q; i++)
fout << random(1, m) << ' ' << random(1, n) << endl;
return 0;
}

测试在Manjaro下进行,不计输入输出时间,估计折算净内存使用(有较大的误差,仅作参考)。即用CCR-Plus的monitor来获取内存使用(iorun不支持),因此设置限制为“10 s, 2 GB”。在使用前,需要把monitor代码中的“进程被阻塞”判定删掉。

__gnu_cxx::rope

定义在<ext/rope>中。可以看出,这和pb_ds一样是一个GNU库,因此只有libstdc++中有。据说也是由SGI开发的,SGI更为著名的是STL。这个模板类原本是为了快速字符串编辑,像是不少古老的文本编辑器题目可以用rope水过。但也可以把存储类型改成int,再加上特殊的存储方式,使复制的复杂度低于线性,于是就成了“可持久化数组”。

使用方法可以查阅man __gnu_cxx::rope,这可能是用NOI Linux的优势,可以随便查函数。为了实现快速复制,一般用指针管理更快,调用复制构造就可以快速复制了。rope的operator []at()都不能作为左值,修改一个数要用replace()方法。

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
#include <fstream>
#include <ext/rope>
#include <ctime>
#include <iostream>
using namespace std;
using __gnu_cxx::rope;
ifstream fin("rope.in");
ofstream fout("rope.out");
const int N = 1000005, M = 1000005, Q = 1000005;
int a[N], ans[Q];
pair<int, int> mod[M], que[Q];
rope<int> *root[M];
int main()
{
int n;
fin >> n;
for (int i = 1; i <= n; i++)
fin >> a[i];
int m;
fin >> m;
for (int i = 1; i <= m; i++)
fin >> mod[i].first >> mod[i].second;
int q;
fin >> q;
for (int i = 1; i <= q; i++)
fin >> que[i].first >> que[i].second;
clock_t t = clock();
root[0] = new rope<int>(a, a + n + 1);//比一个一个加入更高效
for (int i = 1; i <= m; i++)
{
root[i] = new rope<int>(*root[i - 1]);
root[i]->replace(mod[i].first, mod[i].second);
}
for (int i = 1; i <= q; i++)
ans[i] = root[que[i].first]->at(que[i].second);
cerr << clock() - t << endl;
for (int i = 1; i <= q; i++)
fout << ans[i] << endl;
return 0;
}
id O0 time O2 time mem
1 81 ms 43 ms 30.6 MB
2 302 ms 157 ms 107.37 MB
3 1184 ms 639 ms 425.05 MB
4 4313 ms 2227 ms 1425.61 MB
5 N/A N/A N/A

librope

宣传的比SGI rope更高效的实现,参见https://gist.github.com/josephg/3474848#file-results。不过仔细观察就可以发现,librope是用跳表实现的,而不是树,而且测试的项目是大量插入和少量删除操作。

不幸的是,librope是用C写的,不支持模板,只支持字符。更讨厌的是只能用C风格字符串初始化和修改,我不得不把int转为4个字符,而且还得保证中间没有0,一种简单的方法是用255进制。

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
#include <fstream>
#include "rope.h"
#include <ctime>
#include <iostream>
#include <cstdint>
using namespace std;
ifstream fin("rope.in");
ofstream fout("rope.out");
const int N = 1000005, M = 1000005, Q = 1000005;
int a[N], ans[Q];
pair<int, int> mod[M], que[Q];
rope *root[M];
int main()
{
int n;
fin >> n;
for (int i = 1; i <= n; i++)
fin >> a[i];
int m;
fin >> m;
for (int i = 1; i <= m; i++)
fin >> mod[i].first >> mod[i].second;
int q;
fin >> q;
for (int i = 1; i <= q; i++)
fin >> que[i].first >> que[i].second;
clock_t t = clock();
uint8_t *as = new uint8_t[n * 4 + 1];
for (int i = 1; i <= n; i++)
for (int j = 0; j < 4; j++)
{
as[(i - 1) * 4 + j] = a[i] % 255 + 1;
a[i] /= 255;
}
as[n * 4] = 0;
root[0] = rope_new_with_utf8(as);
for (int i = 1; i <= m; i++)
{
root[i] = rope_copy(root[i - 1]);
uint8_t as[5];
for (int j = 0; j < 4; j++)
{
as[j] = mod[i].second % 255 + 1;
mod[i].second /= 255;
}
as[4] = 0;
rope_del(root[i], (mod[i].first - 1) * 4, 4);
rope_insert(root[i], (mod[i].first - 1) * 4, as);
}
for (int i = 1; i <= q; i++)
{
rope *tmp = rope_copy(root[que[i].first]);
rope_del(tmp, 0, (que[i].second - 1) * 4);
rope_del(tmp, 4, (n - que[i].second) * 4);
uint8_t *as = rope_create_cstr(tmp);
for (int j = 3; j >= 0; j--)
ans[i] = ans[i] * 255 + as[j] - 1;
}
cerr << clock() - t << endl;
for (int i = 1; i <= q; i++)
fout << ans[i] << endl;
return 0;
}

注意要先用gcc编译librope,再用g++编译这份代码。测试结果悲惨的证明,根本不支持快速复制,看来上当了!

id O0 time O2 time mem
1 652 ms 574 ms 680.71 MB
2 N/A N/A N/A

主席树

最为简单的可持久化数据结构,用作可持久化数组时,非叶节点不用保存值。其实可以把叶节点的值放在左儿子(或右儿子)来节省不少空间。空间复杂度和倍增一样达到O(nlogn)\mathcal O(n\log n),每次操作都是O(logn)\mathcal O(\log n)

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
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
#include <fstream>
#include <ctime>
#include <iostream>
using namespace std;
ifstream fin("rope.in");
ofstream fout("rope.out");
const int N = 1000005, M = 1000005, Q = 1000005;
int a[N], ans[Q], cnt;
pair<int, int> mod[M], que[Q];
typedef int ptr;//容易区分是下标还是值
struct node
{
int ls, rs;//这个没办法,相当于手动union
} tree[M * 25];
ptr root[M];
void build(ptr &id, int l, int r)
{
id = ++cnt;
if (l == r)
tree[id].ls = a[l];
else
{
int mid = (l + r) / 2;
build(tree[id].ls, l, mid);
build(tree[id].rs, mid + 1, r);
}
}
void modify(ptr &id, ptr ref, int l, int r, int x, int val)
{
id = ++cnt;
if (l == r)
tree[id].ls = val;
else
{
int mid = (l + r) / 2;
if (x <= mid)
{
modify(tree[id].ls, tree[ref].ls, l, mid, x, val);
tree[id].rs = tree[ref].rs;//不修改的直接复制
}
else
{
modify(tree[id].rs, tree[ref].rs, mid + 1, r, x, val);
tree[id].ls = tree[ref].ls;
}
}
}
int query(ptr id, int l, int r, int x)
{
if (l == r)
return tree[id].ls;
int mid = (l + r) / 2;
if (x <= mid)
return query(tree[id].ls, l, mid, x);
return query(tree[id].rs, mid + 1, r, x);
}
int main()
{
int n;
fin >> n;
for (int i = 1; i <= n; i++)
fin >> a[i];
int m;
fin >> m;
for (int i = 1; i <= m; i++)
fin >> mod[i].first >> mod[i].second;
int q;
fin >> q;
for (int i = 1; i <= q; i++)
fin >> que[i].first >> que[i].second;
clock_t t = clock();
build(root[0], 1, n);
for (int i = 1; i <= m; i++)
{
root[i] = root[i - 1];
modify(root[i], root[i], 1, n, mod[i].first, mod[i].second);
}
for (int i = 1; i <= q; i++)
ans[i] = query(root[que[i].first], 1, n, que[i].second);
cerr << clock() - t << endl;
//cerr << cnt << endl;
for (int i = 1; i <= q; i++)
fout << ans[i] << endl;
return 0;
}
id O0 time O2 time
1 5 ms 7 ms
2 18 ms 18 ms
3 92 ms 56 ms
4 334 ms 215 ms
5 1308 ms 871 ms

由于开的是静态数组,内存估计恒定为195 MB。

分块

显然可以做,但是好像复杂度不平衡。

留坑待填

Random Access List

参见https://www.cs.oberlin.edu/~jwalker/ra-list/,一种没听说过的数据结构,还可以做我以前出的数列

留坑待填

统计

可持久化并查集

按秩合并

一种常见的伎俩,即合并时把少的合并到多的,常称为启发式合并。可以保证并查集的树高为O(logn)\mathcal O(\log n)级别。若还是mm次操作qq次询问,复杂度O((m+q)log2n)\mathcal O((m+q)\log^2 n)

路径压缩?

据说不能路径压缩,但是其实好像也没有坏处,但会慢不少。在按秩合并的基础上,复杂度O(αmlogn+qlog2n)\mathcal O(\alpha m\log n+q\log^2 n),因为询问的高度最坏可以达到O(logn)\mathcal O(\log n)

例题

不要认为这个没什么用,其实有时候可以暴力过一些题。

模板

如果有空,可以用上面的方法去试一下:洛谷P3919、洛谷P3402,当然还有著名的BZOJ 3673、3674。

NOIp2013 货车运输

大意是给定一个无向图,多次询问两点之间的路径,最大化路径上最小边。标准做法是建最大生成树,然后树上倍增。最近才做了这题,然后就开始研究可持久化并查集。

但是我一开始没有想到只有最大生成树上的边有效,于是开始折腾。首先注意到可以二分答案,在并查集上连上大于等于答案的边,就可以判断连通性。复杂度O(qmlogz)\mathcal O(qm\log z)。然后发现不用二分,离线处理,直接从大到小加边,每次判断所有询问。复杂度O(qm)\mathcal O(qm),这样就是暴力的60分。

还能怎么改进呢?如果用可持久化并查集,就能保存每次加边后的并查集,询问只要二分判断连通即可。然后我就用了SGI rope,结果还是只有60分,剩下的MLE。

我不肯就此罢休,于是改成了主席树。不幸的是还是60分,剩下的RE,恰好是满数据。主席树的节点应该够了呀。搞来大数据才发现其实是写错了,但只有询问足够大才能遇到这种错误。还好这种错误我在写的时候就注意到了,由于我用了路径压缩,如果询问直接调用getf()就可能会修改前面的版本,而造成后面的版本不兼容。因此要么先复制一份,要么写一个“只读”版本的getf()。加上按秩合并优化常数,最后终于过了。感觉完全是在作死,比赛不太可能写出来。

复杂度O(αmlogn+qlog3n)\mathcal O(\alpha m\log n+q\log^3 n),好像很慢的样子,其实后面的3次log不会满的,大约只有2次。但在当年可能还是会被卡常,毕竟标算只有1次。另外后来发现还是不要路径压缩更好,看来传言还是有道理的。

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
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
#include <iostream>
#include <algorithm>
using namespace std;
const int N = 10005, M = 50005, NODES = 5e6;
struct edge
{
int u, v, w;
bool operator<(const edge &rhs) const
{
return w > rhs.w;
}
} e[M];
typedef int ptr;
struct node
{
int ls, rs;
} tree[NODES];
ptr root[M];
int n, m, cnt, sz[N];
void build(ptr &id, int l, int r)
{
id = ++cnt;
if (l == r)
tree[id].ls = l;
else
{
int mid = (l + r) / 2;
build(tree[id].ls, l, mid);
build(tree[id].rs, mid + 1, r);
}
}
void modify(ptr &id, ptr rt, int l, int r, int x, int val)
{
id = ++cnt;
if (l == r)
tree[id].ls = val;
else
{
int mid = (l + r) / 2;
if (x <= mid)
{
modify(tree[id].ls, tree[rt].ls, l, mid, x, val);
tree[id].rs = tree[rt].rs;
}
else
{
modify(tree[id].rs, tree[rt].rs, mid + 1, r, x, val);
tree[id].ls = tree[rt].ls;
}
}
}
int query(ptr id, int l, int r, int x)
{
if (l == r)
return tree[id].ls;
int mid = (l + r) / 2;
if (x <= mid)
return query(tree[id].ls, l, mid, x);
return query(tree[id].rs, mid + 1, r, x);
}
int getf(ptr &id, int x)
{
int fat = query(id, 1, n, x);
if (fat == x)
return x;
int rt = getf(id, fat);
modify(id, id, 1, n, x, rt);
return rt;
}
int getf_readonly(ptr id, int x)
{
int fat = query(id, 1, n, x);
if (fat == x)
return x;
return getf_readonly(id, fat);
}
int main()
{
ios::sync_with_stdio(false);
cin >> n >> m;
for (int i = 1; i <= m; i++)
cin >> e[i].u >> e[i].v >> e[i].w;
sort(e + 1, e + m + 1);
build(root[0], 1, n);
fill(sz + 1, sz + n + 1, 1);
for (int i = 1; i <= m; i++)
{
root[i] = root[i - 1];
int ru = getf(root[i], e[i].u), rv = getf(root[i], e[i].v);//改成getf_readonly更快
if (ru != rv)
{
if (sz[ru] > sz[rv])
{
modify(root[i], root[i], 1, n, rv, ru);
sz[ru] += sz[rv];
}
else
{
modify(root[i], root[i], 1, n, ru, rv);
sz[rv] += sz[ru];
}
}
}
int q;
cin >> q;
while (q--)
{
int u, v;
cin >> u >> v;
int l = 1, r = m;
while (l < r)
{
int mid = (l + r) / 2;
if (getf_readonly(root[mid], u) == getf_readonly(root[mid], v))
r = mid;
else
l = mid + 1;
}
if (getf_readonly(root[l], u) != getf_readonly(root[l], v))
cout << -1 << endl;
else
cout << e[l].w << endl;
}
return 0;
}

SillyHook's history

随便取了个名字,我们这届的2017-07-17模拟赛T3(去年暑假集训),也有公开的https://wenku.baidu.com/view/16dc09b8250c844769eae009581b6bd97e19bc74.html

当年直接用了SGI rope,结果卡成30分,本来可以写写其他部分分的。