作者自己胡的×要是重复造轮子了就当搬运吧,感觉挺适合往 OI 里搬。
本文将介绍一种可以在 $O(\log n)$ 的时间内支持单点修改区间查询,并能够在 $O(1)$ 时间内完成简单末尾追加的数据结构,与它的一点树上扩展。
注:本文中提及的“线段树”均指朴素递归线段树。序列版本不管是码长还是常数都打不过高级线段树。
- 斜二进制
斜二进制是一种奇怪的进制。它从低到高第 $i\ge1$ 位位权为 $2^i-1$。用斜二进制表示一个数时,需要满足至多一位为二且低于此位者均为零,其余位不超过一。定义一个斜二进制数的最低有效位是它最低位值非零的位。
我们将使用一种归纳构造的方式来生成每个自然数的斜二进制分解。0 的斜二进制是 0。我们假设已经对于 $n$ 生成了它的斜二进制分解,那么 $n+1$ 的斜二进制将在此基础上进行分类:
- $n$ 的斜二进制分解最低有效位为二。
不妨这是第 $i$ 位,那么乘上二的位值即为 $2\times(2^i-1)=2^{i+1}-2$,再加上新多出来的一就是 $2^{i+1}-1$。于是我们清除这两位并向前一位进一即可。
- $n$ 的斜二进制分解最低有效位不为二。
在第一位上加一即可。
不难验证以此法生成的斜二进制分解符合它应有的性质。
这么讲可能比较抽象,我们来举点例子。
0 的斜二进制是 0。
1 的斜二进制是 1。因为 0 的最低有效位不为二。
2 的斜二进制是 2。因为 1 的最低有效位不为二。
3 的斜二进制是 10。因为 2 的最低有效位为二,向前进一。
4 的斜二进制是 11。因为 10 的最低有效位不为二。
5 的斜二进制是 12。因为 11 的最低有效位不为二。
6 的斜二进制是 20。因为 12 的最低有效位为二,向前进一。
7 的斜二进制是 100。因为 20 的最低有效位为二,向前进一。
8 的斜二进制是 101。因为 100 的最低有效位不为二。
9 的斜二进制是 102。因为 101 的最低有效位不为二。
10 的斜二进制是 110。因为 102 的最低有效位为二,向前进一。
11 的斜二进制是 111。因为 110 的最低有效位不为二。
12 的斜二进制是 112。因为 111 的最低有效位不为二。
13 的斜二进制是 120。因为 112 的最低有效位为二,向前进一。
总之大概就是这样。
- 后树
一棵后树维护一个长度为 $2^n-1$ 的序列,其中 $n\ge 1$。
后树的根管辖整个区间,根的左子树是区间的前 $2^{n-1}-1$ 个元素构建的后树,右子树是区间的接下来 $2^{n-1}-1$ 个元素构建的后树。注意最后一个元素不属于任何一棵子树,它是特别的——也正因如此,我们可以将一棵后树的根编号为它最后一个元素的下标。
不难发现,如果假设整个序列从一开始编号并把后树的每个节点的斜二进制写出,那么根管辖的斜二进制区间将是 $(0,100\cdots0]$,左子树的管辖区间是 $ (0,10\cdots0]$,右子树的管辖区间是 $(10\cdots0,20\cdots0]$。我们不难证明,一个点 $i$ 为根的子树的管辖区间的左端点 $j=i-\textit{skew_lowbit}(i)$,其中 $\textit{skew_lowbit}(i)$ 表示 $i$ 在斜二进制下的最低有效位。
- 数据结构
我们考虑将 $n$ 进行斜二进制分解并用一系列后树来维护整个序列。类似这样。
我们对于每个下标 $i$ 维护四个东西:
$a_i$:序列第 $i$ 个元素的值。
$\textit{tr}_i$:序列上 $(i-\textit{skew_lowbit}(i),i]$ 下标的元素合并后的结果(本质上是在后树上的子树“和”)。
$\textit{lb}_i$:$i-\textit{skew_lowbit}(i)$。
$\textit{tf}_i$:点 $i$ 在后树上的父亲(不存在则置零)。
(这几个东西在缩写之前分别是:array
,tree
,left_bound
,tree_father
。)
我们以下将对几个常见操作分别说明。以和为例。
- 节点更新
在除了 $\textit{tr}_i$ 未知以外所有必要信息均已知的情况下计算 $\textit{tr}_i$,俗称 pushup
。由于并不知道 $i$ 会有零个还是两个儿子,我们的手段相当暴力:遍历所有子树进行累和。
void pushup(int x) {
tr[x] = a[x];
for (int y = x - 1; tf[y] == x; y = lb[y]) tr[x] += tr[y];
}
有群 u 指出应该实现得更精细一点所以也写了一版没那么暴力的,效率差距不大:
void pushup(int x) {
if (int p = x - 1; p == lb[x])
tr[x] = a[x];
else
tr[x] = a[x] + tr[p] + tr[lb[p]];
}
- 末尾追加
$a$ 可以直接把追加的值薅过去,$\textit{tr}$ 可以调用 pushup
处理,那么我们实际上要做的就是把 $\textit{tf}$ 和 $\textit{lb}$ 维护好。
注意到只有新加点的 $\textit{lb}$ 和它两个儿子(若有)的 $tf$ 可能需要变化,我们仿照斜二进制分解的归纳构造操作就行。对于最低有效位为二的情况,我们可以简单地找到新加点的两个儿子;对于最低有效位不为二的情况,新加点没有儿子。
获取一个数 $i$ 的斜二进制最低有效位可以直接 $i-\textit{lb}_i$ 简单处理。
void append(int x) {
int p = n, q = lb[n++], r = lb[q];
if (a[n] = x, !q || p - q != q - r)
lb[n] = n - 1;
else
lb[tf[p] = tf[q] = n] = r;
pushup(x);
}
- 建立
嗯做就行。本质是挨个 append
。
void build() {
for (int i = 1; i <= n; ++i) {
int p = i - 1, q = lb[p], r = lb[q];
if (!q || p - q != q - r)
lb[i] = i - 1;
else
lb[tf[p] = tf[q] = i] = r;
pushup(i);
}
}
- 单点修改。
注意到树的结构没有变化,$a$ 变化的只有被修改的位置。于是我们主要考虑 $\textit{tr}$ 的变化。注意到 $tr_i$ 受影响仅当 $i$ 是被修改位置在后树上的祖先。于是我们沿着 $\textit{tf}$ 一路往上遍历祖先就行。
下面这个例子是单点加。
void update(int x, int d) {
for (a[x] += d; x; x = tf[x]) pushup(x);
}
- 区间查询。
从待查区间 $[l,r]$ 的右端点 $r$ 开始向 $l$ 跳,有 $\textit{tr}$ 跳 $\textit{tr}$ 没 $\textit{tr}$ 跳 $a$ 就行。可能讲得比较抽象,下面挂个代码。标准库函数 exchange(u, v)
的含义等同于是执行 u = v
并返回执行前的 u
。
int query(int l, int r) {
int ans = 0;
while (r >= l)
if (int x = lb[r]; x < l - 1)
ans += a[r--];
else
ans += tr[exchange(r, x)];
return ans;
}
接下来我们证一下这个东西的复杂度。我们定义一个位置 $i$ 的高度 $h_i$ 是以它为根的后树的大小加一以二为底的对数。不难发现每个点是不超过 $O(\log n)$ 的正整数。特别地,我们定义零的高度是全局最高高度加一。
于是我们有了这么一条事情:
- $h_{\textit{lb}_{\textit{lb}_i}}\ge h_i+1$。
写个斜二进制,这是显然的。
这条告诉我们,在 if
的第一个分支被执行到之前,$h_r$ 至多两步一升,复杂度是 $O(\log n)$ 的。
而在此之后我们又会有这么第二条事情:
- 若 $r$ 被执行第一个分支时 $h_r=x$,那么此后 $h_r
因为 $r$ 之前最后一个满足 $h_i\ge x$ 的位置 $i$ 是 $\textit{lb}_r\lt x-1$,它不可能进入循环体。
而根据上面的第一条,走第一个分支之后,至多再走一步第二个分支,一定会再回到第一个分支(否则连走两步之后 $h_r$ 加一,与第二条矛盾)。于是我们得知,在 if
的第一个分支被执行到之后,$h_r$ 至多两步一减,复杂度是 $O(\log n)$ 的。
综上,区间查询的复杂度是 $O(\log n)$ 的。
挂个【模板】树状数组 1 的完整代码吧:
#include <bits/stdc++.h>
using namespace std;
constexpr int N = 5e5 + 9;
int n, m, lb[N], tf[N], a[N], tr[N];
void pushup(int x) {
tr[x] = a[x];
for (int y = x - 1; tf[y] == x; y = lb[y]) tr[x] += tr[y];
}
void build() {
for (int i = 1; i <= n; ++i) {
int p = i - 1, q = lb[p], r = lb[q];
if (!q || p - q != q - lb[q])
lb[i] = i - 1;
else
lb[tf[p] = tf[q] = i] = r;
pushup(i);
}
}
void update(int x, int d) {
for (a[x] += d; x; x = tf[x]) pushup(x);
}
int query(int l, int r) {
int ans = 0;
while (r >= l)
if (int x = lb[r]; x < l - 1)
ans += a[r--];
else
ans += tr[exchange(r, x)];
return ans;
}
signed main() {
cin.tie(nullptr)->sync_with_stdio(false);
cin >> n >> m;
for (int i = 1; i <= n; ++i) cin >> a[i];
build();
for (int op; m; --m)
if (cin >> op; op == 1) {
int x, d;
cin >> x >> d;
update(x, d);
} else {
int l, r;
cin >> l >> r;
cout << query(l, r) << '\n';
}
return cout << flush, 0;
}
不难发现,该算法在简洁方面对比线段树有较大优势,功能强度方面介于树状数组与线段树之间。常数应该不小。
它的作用?总之大概能提供一个新的选择吧。写起来比线段树快,比树状数组逻辑清晰(不需要前缀和差分的转)。
事实上这个东西还能上树。不难发现每个点拖一个区间这个事情可以直接挂在每个点上往祖先去拖。不难发现这样依然可以在 $O(\log n)$ 的时间内刻画出树上的一条路径。另一方面,序列上需要倍增的场景其实并不多见,即使有也许多为双向,但树上倍增一般都是往祖先,这使得它可以更好地适应形如 “找到第一个满足某条件”的祖先之类的问题。
这种数据结构的预处理时空是 $O(n)$,查询是 $O(\log n)$,同时代码相当简短,在大部分情况下可以作为树上倍增的上位替代。功能性除了支持动态加叶子是优势以外一般不如树剖,但时空常数和码量均更优(笔者个人实现的时间是树剖的三分之一),同时逻辑也更清晰。严格劣 ST 树。
初始化非常简单:
void dfs(int x) {
int p = fa[x], q = lb[p], r = lb[q];
lb[x] = d[p] - d[q] != d[q] - d[r] ? p : r;
for (d[x] = d[p] + 1; int y : es[x])
if (y != fa[x]) fa[y] = x, dfs(y);
}
这里给出一份查 LCA 的代码,不难发现略加修改即可查询静态链信息:
int lca(int u, int v) {
if (d[u] < d[v]) swap(u, v);
while (d[u] > d[v])
if (d[lb[u]] < d[v])
u = fa[u];
else
u = lb[u];
while (u != v)
if (lb[u] == lb[v])
u = fa[u], v = fa[v];
else
u = lb[u], v = lb[v];
return u;
}
同时可以看到其实逻辑和倍增是很类似的,易于配套理解。
如果要维护路径信息的话,以[JLOI2015] 城池攻占 为例,初始化稍微改一下:
void dfs(int x) {
int p = fa[x], q = lb[p], r = lb[q];
if (!q || d[p] - d[q] != d[q] - d[r])
lb[x] = p, tr[x] = vl[x];
else
lb[x] = r, tr[x] = vl[x] + tr[p] + tr[q];
for (d[x] = d[p] + 1; int y : es[x]) dfs(y);
}
(这题是输入每个点的爹所以向下搜不需要判。)
查询就是每次跳 $\textit{lb}$,跳不了就跳爹:
int calc(int x, int h) {
while (x)
if (func f = tr[x]; h >= f.xl && h <= f.xr)
h = f(h), x = lb[x];
else if (f = vl[x]; h >= f.xl && h <= f.xr)
h = f(h), x = fa[x];
else
break;
return x;
}
Upd on 2024.3.1:来活了,注意到 DAG 支配树需要支持的操作是动态加叶子和求 LCA,专业对口。目前是洛谷最优解。
#include <bits/stdc++.h>
using namespace std;
constexpr int S = 1 << 21, N = 65535;
char buf[S], *p1, *p2, obuf[S], *O = obuf;
inline char gc() {
if (p1 == p2) {
p2 = (p1 = buf) + cin.read(buf, S).gcount();
if (p1 == p2) return EOF;
}
return *p1++;
}
inline int rd() {
char ch;
while (!isdigit(ch = gc()))
;
int x = ch & 0xf;
while (isdigit(ch = gc())) x = x * 10 + (ch & 0xf);
return x;
}
inline void prtln(int x) {
char s[7];
int t = 0;
if (!x)
*O++ = '0';
else {
while (x) s[t++] = x % 10 | '0', x /= 10;
while (t) *O++ = s[--t];
}
*O++ = '\n';
}
basic_string<int> es[N];
int n, fa[N], d[N], lb[N], dg[N], sz[N];
inline void lca(int& u, int v) {
if (!~u) return u = v, void();
if (d[u] < d[v]) swap(u, v);
while (d[u] > d[v]) d[lb[u]] < d[v] ? (u = fa[u]) : (u = lb[u]);
while (u != v)
lb[u] == lb[v] ? (u = fa[u], v = fa[v]) : (u = lb[u], v = lb[v]);
}
inline void addf(int x) {
int p = fa[x], q = lb[p], r = lb[q];
d[x] = d[p] + 1, lb[x] = d[p] - d[q] != d[q] - d[r] ? p : r;
}
inline void build() {
memset(fa + 1, -1, n * sizeof(int));
int q[N], l = 0, r = 0, u;
for (int i = 1; i <= n; ++i)
if (!dg[i]) fa[q[r++] = i] = 0, d[i] = 1;
while (l < r)
for (addf(u = q[l++]); int v : es[u])
if (lca(fa[v], u), !--dg[v]) q[r++] = v;
while (r--)
if (int f = fa[u = q[r]]) sz[f] += sz[u] + 1;
}
signed main() {
cin.tie(nullptr)->sync_with_stdio(false);
n = rd();
for (int i = 1; i <= n; ++i)
while (int j = rd()) es[j].push_back(i), ++dg[i];
build(), for_each_n(sz + 1, n, prtln);
return cout.write(obuf, O - obuf).flush(), 0;
}
Upd on 2024.4.10:整了个例题 https://www.luogu.com.cn/problem/U421630。
同时这个东西还可以用来写全局平衡二叉树。但这个时候我们会意识到一个很难受的事情:根的位置被固定了,无法再保证每个点子树大小和父亲的关系。但是这个事情我们魔怔一下就能解决:把本来的左右儿子之间插一个中儿子就行,等于是一个点有仨儿子,其中中间那个的子树大小为一。那么我们的性质就又没问题了。写出来会比朴素 GBST 稍微少一点特判。
例题:[ZJOI2008] 树的统计。写的线性建树。
#include <bits/stdc++.h>
#define int long long
using namespace std;
constexpr int N = 3e4 + 9;
int n, m, fa[N], top[N], sn[N], sz[N], d[N], lb[N], tf[N], stk[N], sum[N], tp;
vector<int> es[N];
struct info {
int sm, mx;
info() : sm(0), mx(INT_MIN) {}
info(int x) : sm(x), mx(x) {}
info(int sm, int mx) : sm(sm), mx(mx) {}
info& operator+=(const info& to) {
return sm += to.sm, mx = max(mx, to.mx), *this;
}
info operator+(const info& to) const { return info(*this) += to; }
} a[N], tr[N];
int dfs1(int x) {
d[x] = d[fa[x]] + 1;
for (int y : es[x]) {
erase(es[y], fa[y] = x);
if (sz[x] += dfs1(y), sz[y] > sz[sn[x]]) sn[x] = y;
}
return ++sz[x];
}
inline void pushup(int x) {
tr[x] = a[x];
for (int y = fa[x]; tf[y] == x; y = lb[y]) tr[x] += tr[y];
}
int build(int l, int r) {
int x = stk[r--];
if (lb[x] = stk[l - 1], l <= r) {
int i = upper_bound(sum + l, sum + r + 1, (sum[l - 1] + sum[r]) >> 1) - sum;
int y = stk[i];
tr[y] = a[y], lb[y] = stk[i - 1], tf[y] = x;
l < i && (tf[build(l, i - 1)] = x), i < r && (tf[build(i + 1, r)] = x);
}
return pushup(x), x;
}
void dfs2(int x, int t) {
top[stk[++tp] = x] = t, sum[tp] = sum[tp - 1] + sz[x];
if (int z = sn[x]) {
for (sum[tp] -= sz[z], dfs2(z, t); int y : es[x])
if (y != z) dfs2(y, y);
} else
build(tp - (d[x] - d[t]), tp);
--tp;
}
inline int lca(int u, int v) {
while (top[u] != top[v]) {
if (d[top[u]] < d[top[v]]) swap(u, v);
u = fa[top[u]];
}
return d[u] < d[v] ? u : v;
}
inline void update(int x, int t) {
for (a[x] = t; x; x = tf[x]) pushup(x);
}
inline info qlink(int u, int z) {
info ans;
while (d[u] >= z)
if (int x = lb[u]; d[x] < z - 1)
ans += a[exchange(u, fa[u])];
else
ans += tr[exchange(u, x)];
return ans;
}
inline info query(int u, int v) {
int z = d[lca(u, v)];
return qlink(u, z) + qlink(v, z + 1);
}
signed main() {
cin.tie(nullptr)->sync_with_stdio(false);
cin >> n;
for (int i = 1, u, v; i < n; ++i)
cin >> u >> v, es[u].push_back(v), es[v].push_back(u);
for (int i = 1, x; i <= n; ++i) cin >> x, a[i] = x;
dfs1(1), dfs2(1, 1);
for (cin >> m; m; --m) {
string s;
int x, y;
if (cin >> s >> x >> y, s == "CHANGE")
update(x, y);
else if (auto [sm, mx] = query(x, y); s == "QMAX")
cout << mx << '\n';
else
cout << sm << '\n';
}
return cout << flush, 0;
}
Upd on 2024.3.9:
感谢 @ClHg2 的指出,直接把中儿子塞到左儿子里复杂度也是对的,因为一个点的爷爷权和一定至少是自己的两倍所以树高依然是 $\log$ 的。常数会更小一点。
改完的代码:
#include <bits/stdc++.h>
#define int long long
using namespace std;
constexpr int N = 3e4 + 9;
int n, m, fa[N], top[N], sn[N], sz[N], d[N], lb[N], tf[N], stk[N], sum[N], tp;
vector<int> es[N];
struct info {
int sm, mx;
info() : sm(0), mx(INT_MIN) {}
info(int x) : sm(x), mx(x) {}
info(int sm, int mx) : sm(sm), mx(mx) {}
info& operator+=(const info& to) {
return sm += to.sm, mx = max(mx, to.mx), *this;
}
info operator+(const info& to) const { return info(*this) += to; }
} a[N], tr[N];
int dfs1(int x) {
d[x] = d[fa[x]] + 1;
for (int y : es[x]) {
erase(es[y], fa[y] = x);
if (sz[x] += dfs1(y), sz[y] > sz[sn[x]]) sn[x] = y;
}
return ++sz[x];
}
inline void pushup(int x) {
tr[x] = a[x];
for (int y = fa[x]; tf[y] == x; y = lb[y]) tr[x] += tr[y];
}
int build(int l, int r) {
int x = stk[r--];
if (lb[x] = stk[l - 1], l <= r) {
int i = upper_bound(sum + l, sum + r + 1, (sum[l - 1] + sum[r]) >> 1) - sum;
tf[build(l, i)] = x, i < r && (tf[build(i + 1, r)] = x);
}
return pushup(x), x;
}
void dfs2(int x, int t) {
top[stk[++tp] = x] = t, sum[tp] = sum[tp - 1] + sz[x];
if (int z = sn[x]) {
for (sum[tp] -= sz[z], dfs2(z, t); int y : es[x])
if (y != z) dfs2(y, y);
} else
build(tp - (d[x] - d[t]), tp);
--tp;
}
inline int lca(int u, int v) {
while (top[u] != top[v]) {
if (d[top[u]] < d[top[v]]) swap(u, v);
u = fa[top[u]];
}
return d[u] < d[v] ? u : v;
}
inline void update(int x, int t) {
for (a[x] = t; x; x = tf[x]) pushup(x);
}
inline info qlink(int u, int z) {
info ans;
while (d[u] >= z)
if (int x = lb[u]; d[x] < z - 1)
ans += a[exchange(u, fa[u])];
else
ans += tr[exchange(u, x)];
return ans;
}
inline info query(int u, int v) {
int z = d[lca(u, v)];
return qlink(u, z) + qlink(v, z + 1);
}
signed main() {
cin.tie(nullptr)->sync_with_stdio(false);
cin >> n;
for (int i = 1, u, v; i < n; ++i)
cin >> u >> v, es[u].push_back(v), es[v].push_back(u);
for (int i = 1, x; i <= n; ++i) cin >> x, a[i] = x;
dfs1(1), dfs2(1, 1);
for (cin >> m; m; --m) {
string s;
int x, y;
if (cin >> s >> x >> y, s == "CHANGE")
update(x, y);
else if (auto [sm, mx] = query(x, y); s == "QMAX")
cout << mx << '\n';
else
cout << sm << '\n';
}
return cout << flush, 0;
}
具体有啥更多的应用可能还有待进一步发掘。感觉是个非常适合引入 OI 的东西,尤其是平替树上倍增的部分(?
以上。