最近共同祖先：FC&B算法（Lowest Common Ancestor : Farach-Colton & Bender Algorithm）

假设 $G$ 是一棵树。对于每一个形式为 $(u,v)$ 的查询，我们想找到节点 $u$ 和 $v$ 的最近共同祖先，也就是说，我们想找到一个节点 $w$ ，它位于从 $u$ 和 $v$ 到根节点的路径上，如果有多个节点，我们选择离根节点最远的那个。换句话说，所需的节点 $w$ 是 $u$ 和 $v$ 的最近祖先。特别情况下，如果 $u$ 是 $v$ 的祖先，那么 $u$ 就是它们的最近共同祖先。

本文将描述的算法是由 Farach-Colton 和 Bender 共同开发的，是一种渐进式的最优算法。

算法

我们使用经典的 LCA 问题归约为 RMQ 问题。用 DFS 遍历树的所有节点，并保留一个包含所有访问过的节点和这些节点的高度的数组。两个节点 $u$ 和 $v$ 的 LCA 是在遍历中出现的 $u$ 和 $v$ 之间的节点，它的高度最小。

在下图（表）中，你可以看到一种情况下的欧拉遍历和访问过的节点以及它们的高度。

讯享网
在 Lowest Common Ancestor （可能需要科学上网）一文中读到更多关于这种归约的信息。在上述这篇文章中，一个范围的最小值是通过平方分解求得的 $O\left(\sqrt{N}\right)$ 或者用线段树求得的 $O\left(\log N\right)$ 。

而在本文中，我们将探讨如何能在 $O\left(1\right)$ 的时间内解决给定的范围最小值查询，同时仍然只用 $O\left(N\right)$ 的时间复杂度进行预处理。

需要注意的是，归约的 RMQ 问题是非常具体的：数组中任何两个相邻的元素正好相差 $1$（因为数组的元素只不过是按遍历顺序访问的节点的高度，我们要么去找一个后代，在这种情况下，下一个元素大 $1$ ；要么回到祖代，在这种情况下，下一个元素小 $1$ ）。Farach-Colton 和 Bender 算法描述的正是这种专门的 RMQ 问题的解决方案。

让我们用 $A$ 来表示我们要对其进行范围最小查询的数组。而 $N$ 将是 $A$ 的大小。

有一个简单的数据结构，需要进行 $O\left(N\log N\right)$ 的预处理，我们可以用它来解决 RMQ 问题，其中每个查询为 $O(1)$ 的 稀疏表 。我们创建一个表 $T$ ，其中每个元素 $T[i][j]$ 等于 $A$ 在区间 $[i,i+2^j-1]$ 中的最小值。显然 $0\le j\le \lceil \log N\rceil$ ，因此稀疏表的大小将是 $O\left(N\log N\right)$ ，你可以在 $O(N\log N)$ 内轻松建立表，其中 $T[i][j]=min(T[i][j-1],T[i+2^{j-1}][j-1])$ 。

我们如何用这种数据结构在 $O(1)$ 内回答一条查询 RMQ ？假设收到的查询是 $[l,r]$ ，那么答案是 $min(T[l][sz],T[r-2^{sz}+1][sz])$ ，其中 $sz$ 是最大的指数，使 $2^{sz}$ 不大于范围长度 $r-l+1$ 。事实上，我们可以把范围 $[l,r]$ 覆盖在两个长度为 $2^{sz}$ 的片段上，一端从 $l$ 开始，另一端在 $r$ 结束。这些片段是重叠的，但这并不影响我们的计算。为了真正实现每个查询的时间复杂度为 $O(1)$ ，我们需要知道从 $1$ 到 $N$ 的所有可能的长度的 $sz$ 值，这是很容易预先计算的。

现在我们想将预处理的复杂度降低到 $O\left(N\right)$ 。

我们把数组 $A$ 分成大小为 $K=\frac{{\log }_{2}N}{2}$ 的块。对于每个块，我们计算出最小的元素，并将其存储在数组 $B$ 中， $B$ 的大小为 $\frac{N}{K}$ 。我们从数组 $B$ 构建一个稀疏表。它的规模和时间的复杂度为： $$\frac{N}{K}\log \left(\frac{N}{K}\right)=\frac{2N}{logN}\log \left(\frac{2N}{\log N}\right)=\frac{2N}{\log N}\left(1+\log \left(\frac{N}{\log N}\right)\right)\le \frac{2N}{\log N}+2N=O\left(N\right)$$

现在我们只需要学习如何快速解决每个区块内的最小范围查询。事实上，如果收到的范围最小查询是 $[l,r]$ ，并且 $l$ 和 $r$ 在不同的区块中，那么答案就是以下三个值的最小值：从 $l$ 开始的 $l$ 区块的后缀的最小值，从 $r$ 结束的 $r$ 区块的前缀的最小值，以及这些区块之间的最小值。这之间的区块的最小值可以用稀疏表在 $O\left(1\right)$ 内回答。因此，我们只剩下区块内的最小值查询。

这里我们将利用数组的属性。请记住，数组中只是树中的高度值将总是相差 $1$ 。如果我们去掉一个区块的第一个元素，并将其与区块中的其他每一项相减，那么每个区块都可以由一个长度为 $K-1$ 的序列来识别，该序列由数字 $1$ 和 $0$ 组成，由于这些区块非常小，只有几个不同的序列可以出现。可能的序列数量为： $$2^{K-1}=2^{\frac{\log N}{2}-1}=\frac{\sqrt{2^{\log N}}}{2}=\frac{\sqrt{N}}{2}$$

因此，不同区块的数量是 $O\left(\sqrt{N}\right)$ ，因此我们可以在 $O\left(\sqrt{N}{K}^2\right)=O\left(\sqrt{N}{\log }^{2}N\right)=O\left(N\right)$ 时间内预先计算出所有不同区块内的范围最小值查询结果。为了实现这一点，我们可以用一个长度为 $K-1$ 的比特掩码来描述一个区块的特征（适合于标准int），并将最小值的索引存储在一个大小为 $O\left(\sqrt{N}{\log }^{2}N\right)$ 的数组 $block[mask][l][r]$ 中。

因此，我们学会了如何在每个区块内预先计算范围内的最小查询，以及在一个区块范围内的范围内的最小查询，所有这些都在 $O\left(N\right)$ 的时间中。通过这些预计算，我们最多使用四个预先计算的值：包含 $l$ 的区块的最小值，包含 $r$ 的区块的最小值，以及它们之间区块重叠段的两个最小值，就可以在 $O\left(1\right)$ 内解决每个查询。

实现

int n; vector<vector<int>> adj; int block_size, block_cnt; vector<int> first_visit; vector<int> euler_tour; vector<int> height; vector<int> log_2; vector<vector<int>> st; vector<vector<vector<int>>> blocks; vector<int> block_mask; void dfs(int v, int p, int h) { 
    first_visit[v] = euler_tour.size(); euler_tour.push_back(v); height[v] = h; for (int u : adj[v]) { 
    if (u == p) continue; dfs(u, v, h + 1); euler_tour.push_back(v); } } int min_by_h(int i, int j) { 
    return height[euler_tour[i]] < height[euler_tour[j]] ? i : j; } void precompute_lca(int root) { 
    // get euler tour & indices of first occurences first_visit.assign(n, -1); height.assign(n, 0); euler_tour.reserve(2 * n); dfs(root, -1, 0); // precompute all log values int m = euler_tour.size(); log_2.reserve(m + 1); log_2.push_back(-1); for (int i = 1; i <= m; i++) log_2.push_back(log_2[i / 2] + 1); block_size = max(1, log_2[m] / 2); block_cnt = (m + block_size - 1) / block_size; // precompute minimum of each block and build sparse table st.assign(block_cnt, vector<int>(log_2[block_cnt] + 1)); for (int i = 0, j = 0, b = 0; i < m; i++, j++) { 
    if (j == block_size) j = 0, b++; if (j == 0 || min_by_h(i, st[b][0]) == i) st[b][0] = i; } for (int l = 1; l <= log_2[block_cnt]; l++) { 
    for (int i = 0; i < block_cnt; i++) { 
    int ni = i + (1 << (l - 1)); if (ni >= block_cnt) st[i][l] = st[i][l-1]; else st[i][l] = min_by_h(st[i][l-1], st[ni][l-1]); } } // precompute mask for each block block_mask.assign(block_cnt, 0); for (int i = 0, j = 0, b = 0; i < m; i++, j++) { 
    if (j == block_size) j = 0, b++; if (j > 0 && (i >= m || min_by_h(i - 1, i) == i - 1)) block_mask[b] += 1 << (j - 1); } // precompute RMQ for each unique block int possibilities = 1 << (block_size - 1); blocks.resize(possibilities); for (int b = 0; b < block_cnt; b++) { 
    int mask = block_mask[b]; if (!blocks[mask].empty()) continue; blocks[mask].assign(block_size, vector<int>(block_size)); for (int l = 0; l < block_size; l++) { 
    blocks[mask][l][l] = l; for (int r = l + 1; r < block_size; r++) { 
    blocks[mask][l][r] = blocks[mask][l][r - 1]; if (b * block_size + r < m) blocks[mask][l][r] = min_by_h(b * block_size + blocks[mask][l][r], b * block_size + r) - b * block_size; } } } } int lca_in_block(int b, int l, int r) { 
    return blocks[block_mask[b]][l][r] + b * block_size; } int lca(int v, int u) { 
    int l = first_visit[v]; int r = first_visit[u]; if (l > r) swap(l, r); int bl = l / block_size; int br = r / block_size; if (bl == br) return euler_tour[lca_in_block(bl, l % block_size, r % block_size)]; int ans1 = lca_in_block(bl, l % block_size, block_size - 1); int ans2 = lca_in_block(br, 0, r % block_size); int ans = min_by_h(ans1, ans2); if (bl + 1 < br) { 
    int l = log_2[br - bl - 1]; int ans3 = st[bl+1][l]; int ans4 = st[br - (1 << l)][l]; ans = min_by_h(ans, min_by_h(ans3, ans4)); } return euler_tour[ans]; }