KD Tree 瞎写

这几天学了 \(\rm KDTree\),但我只会两种操作,就是 OI wiki 上的两种操作,然后写了点题,写点东西总结一下初学 \(\rm KDTree\) 的收获。

KD Tree 是个啥

\(\rm KDTree\),当中的 \(\rm D\) 表示 \(\rm Dimension\),即维度,所以顾名思义,\(\rm KDTree\) 是被用来解决 \(k\) 维问题的一棵树,更具体的讲,是一种存储 \(k\) 维点集合的数据结构。

\(\rm KDTree\) 的结构是 \(\rm BST\),即二叉搜索树,每个子树维护的是一个 \(k\) 维的长方体,而每个点维护的为其所在子树长方体内的一个点。可以发现,这样要求儿子的长方体包含在父亲的内,且一个点的子树内含有这个 \(k\) 维长方体内的所有点。这样 \(\rm KDTree\) 的结构已经呼之欲出了,我们只需要从当前的长方体内选一个维度和一个点,按照这个点的这个维度把长方体分成两个部分,把这两个部分作为这个点的左右子树,之后对这两个子树递归建树,直到长方体退化为单个点。

当然,很多情况下我们并不需要用到长方体内的所有点,这时候可以用动态开点减少树上结点的个数。

KD Tree 咋建啊

这里要把问题分两类讨论,一类是静态的点集,一类是动态的点集。

对于静态的点集,我们只需要在 \(\rm KDTree\) 结构的描述基础上,稍微修改一些部分就可以得到建树的方法。以长方体的切割为基本思路,每次选出一个点集内的点作为子树的根,然后选出一个维度切割,将切割出的点作为左右子树,如果一个长方体内没有点,则这个子树置为空即可。

当然,上面的过程说了跟没说一样,因为根本不知道选哪个点和选哪个维度啊。我们可以感性理解一下,肯定是要尽可能把点分得均匀,才能保证 \(\rm BST\) 的层数尽可能低,所以维度应该选择方差最大的那个,点应该选择当前维度上的中位数。这样,我们可以保证建出的 \(\rm BST\) 树高是 \(\mathcal{O}(\log n)\) 级别的。

对应到具体代码实现上,求方差没啥问题,但求中位数怎么办?如果用 std::sort,时间复杂度是 \(\mathcal{O}(n\log^2n)\) 的。但其实,如果您做过 P1138 第 k 小整数 就会知道,求中位数是可以做到 \(\mathcal{O}(n)\) 的,从而整个算法的时间复杂度可以降至 \(\mathcal{O}(n\log n)\)。当然,我们也不需要手写快排,只需调用 std::nth_element 即可,它可以在 \(\mathcal{O}(n)\) 的时间复杂度内把某一位数放置到它排序后该在的位置上,且它左边全是小于它的,右边全是大于它的。

// 这里函数的参数 l,r,即表示当前正在处理 [l,r] 内的点构成的长方体对应的子树。
int build(int l, int r)
{
    if (l > r) return 0;
    int mid = (l + r) >> 1;
    db av1 = 0, av2 = 0, s1 = 0, s2 = 0;
    for (int i = l; i <= r; ++i) av1 += p[i].x, av2 += p[i].y;
    av1 /= (r - l + 1); av2 /= (r - l + 1);
    for (int i = l; i <= r; ++i)
        s1 += (av1 - p[i].x) * (av1 - p[i].x), s2 += (av2 - p[i].y) * (av2 - p[i].y);
    if (s1 > s2) std::nth_element(p + l, p + mid, p + r + 1, [&](const ::node& x, const ::node& y) { return x.x < y.x; }), a[mid].D = 1;
    else std::nth_element(p + l, p + mid, p + r + 1, [&](const ::node& x, const ::node& y) { return x.y < y.y; }), a[mid].D = 2;
	// D 表示按照哪个维度分割的
    a[mid].ls = build(l, mid - 1); a[mid].rs = build(mid + 1, r);
    return pushup(mid), mid;
}

而对于动态的点集,插入操作类似平衡树的插入操作,根据记录的分割维度,不断往下递归,直到找到空节点。但这个空节点分割的维度就只能任意指定一个了,毕竟我们不知道它的子树会是啥样。但这样下去,肯定会出现不平衡的情况,这时,因为父亲儿子之间除了大小关系,还有包含关系,不能轻易旋转,这里我们通常采用的是替罪羊树的暴力重构思想。(这个方法似乎被 lxl 干掉了)即如果出现一个子树不平衡,那就把这个子树全部拍扁重构,具体来讲就是把这个子树内的点集当做一个静态点集,按照上面的方法重新建 \(\rm KDTree\),代替原来的子树。对于不平衡的定义,我们定义一个常数 \(\alpha\),当有一个儿子的子树的占整个树的比例超过 \(\alpha\) 时,我们就认为这个树不平衡。\(\alpha\) 我喜欢取 \(0.75\)。当然,这要求我们维护子树的大小 siz

对于删除操作,可以继续用替罪羊树的思路,惰性删除,删除后不再计算它对子树信息的贡献即可,类似单点修改。

这样做,类似替罪羊树,\(\rm KDTree\) 的树高依然是 \(\mathcal{O}(\log n)\) 的。

// 获得点集
void get(int x) { if (x) get(a[x].ls), st[++tp] = x, get(a[x].rs); }
// 根据获得的点集建 KDTree
int build(int l, int r)
{
    if (l > r) return 0;
    int mid = (l + r) >> 1;
    double av1 = 0, av2 = 0, s1 = 0, s2 = 0;
    for (int i = l; i <= r; ++i) av1 += p[st[i]].x, av2 += p[st[i]].y;
    av1 /= (r - l + 1); av2 /= (r - l + 1);
    for (int i = l; i <= r; ++i) 
        s1 += (av1 - p[st[i]].x) * (av1 - p[st[i]].x), s2 += (av2 - p[st[i]].y) * (av2 - p[st[i]].y);
    if (s1 > s2) std::nth_element(st + l, st + mid, st + r + 1, [&](const int& x, const int& y) { return p[x].x < p[y].x; }), a[st[mid]].D = 1;
    else std::nth_element(st + l, st + mid, st + r + 1, [&](const int& x, const int& y) { return p[x].y < p[y].y; }), a[st[mid]].D = 2;
    a[st[mid]].ls = build(l, mid - 1); a[st[mid]].rs = build(mid + 1, r);
    return pushup(st[mid]), st[mid];
}
void rebuild(int& x) { tp = 0; get(x); x = build(1, tp); }
bool bad(int x) { return r * a[x].siz <= (double)std::max(a[a[x].ls].siz, a[a[x].rs].siz); }
void insert(int& x, int v)
{
    if (!x) return pushup(x = v);
    if (a[x].D == 1) p[v].x <= p[x].x ? insert(a[x].ls, v) : insert(a[x].rs, v);
    else p[v].y <= p[x].y ? insert(a[x].ls, v) : insert(a[x].rs, v);
	// 如果加入点后,这个子树不平衡了,就重构
    pushup(x); if (bad(x)) rebuild(x);
}

KD Tree 咋用啊

首先是邻域查询,其实我也不知道这个名词是啥意思但是 OI wiki 就是这么写的那还是直接搬上来吧,,这个主要是查询点集内距离某点最近,最远,或第 \(k\) 近的点。

首先,这个东西只是在 \(\rm KDTree\) 上搜索的话,单次查询是 \(\mathcal{O}(n)\) 的,因为它只是一种存储方式,并不能优化你找所求点的过程,还是只能挨个枚举。

但因为 \(\rm KDTree\) 上单个节点能存储它的子树的信息,我们可以根据这个信息来决定是否进入这个子树,从而在随机的数据下,达到较好的效率。对于一个子树,我们在它的根节点,维护这个子树代表的长方体,具体来讲,就是维护每一维度下,子树内点的最小值和最大值。从而,我们就能知道查询点到子树内的点距离理论最小值是多少。如果这个理论最小值都大于当前的答案,则这个子树一定不会造成贡献,不进入即可。

不仅如此,当我们面对两个子树都可以进入的情况时,优先进入最可能更新答案的子树,之后再决定是否进入另一棵子树也能优化效率。对于最可能更新答案,我们通常用估价函数来比较。本题中的估价函数是查询点到子树内点的理论最小值,也就是刚刚和当前答案比较的那个东西。对于其他邻域查询的题,估价函数要根据题意变化,比如把最小值换成最大值之类的。

当然,即使有以上的种种优化,还是能构造数据把 \(\rm KDTree\) 卡成单次查询 \(\mathcal{O}(n)\) 的,这点可以见 P7883 平面最近点对(加强加强版)

过个加强版 还是没啥问题的。

注意,使用 \(\rm KDTree\) 解决问题时,都应该考虑一下,是否能用一些常用的剪枝来加速搜索过程,有时候一个剪枝的影响是非常大的。

然后是一些邻域查询的板子,基本都是把最小改为最大,或者稍微魔改魔改,就不详细说了。

然后就是进行矩阵查询。这个东西我们就能保证复杂度了,提前上结论,查询 \(k\) 维长方体的信息时间复杂度可以做到单次最坏 \(\mathcal{O}(n^{1-k})\)。这个东西相比 \(\rm cdq\) 分治优势在在线,相比树套树优势在空间,共有的优势是当 \(k\) 大的时候码量不会相应增加许多。

这个东西需要以下三种剪枝:

  • 查询矩阵完全覆盖的子树不再进入,直接计入子树上记录的整棵子树的贡献。
  • 查询矩阵与当前子树完全不相交的子树不再进入。
  • 当查询的内容是最值相关时,当前子树无法对答案造成贡献的子树不再进入,这需要引入一个估价函数。

第三点的估价函数,我们通常可以采用这个子树在极限情况下造成的贡献,如果这个都无法更新答案,就不再进入。

P4148 简单题 就是为 \(\rm KDTree\) 量身打造的一个单点修改矩阵查询的题。这里查询和就没有第三条优化了。直接 上代码

能发现,这个矩阵查询,可以处理 \(k\) 维偏序的问题。板子题就是单纯的矩阵查询没啥好说的,这里给出两道优化 \(\rm dp\) 的例子。

P3769 [CH弱省胡策R2]TATT,一个四维意义下的最长不下降子序列。\(\rm dp\) 很好想:

\[f_i=\max_{a_j\le a_i,b_j\le b_i,c_j\le c_i,d_j\le d_i}\{f_j\}+1 \]

然后我们把其中一维作为主关键字排序后,就能得到转移顺序了。转移时就是一个三维偏序,或者说,三维矩阵最值查询。这个用 \(\rm KDTree\) 即可做到 \(\mathcal{O}(n^{\frac{5}{3}})\) 的复杂度。这里要注意,优化 \(\rm dp\) 时,由于我们每个坐标都是已知的,相当于是静态点集,所以没必要每次都插入新的点。对于 \(\rm dp\) 值的改变带来的修改,可以用类似线段树单点修改的方式递归修改 \(\mathcal{O}(\log n)\) 个结点,从而减小常数。\(\tt code\)

L. Il Derby della Madonnina,我快 -50 了还没过的一道题,正解是似乎二维偏序我非要莽三维偏序。首先,有个比较好想的 \(\rm dp\)

\[f_i=\max_{|a_i-a_j|\le v(t_i-t_j),i>j}\{f_j\}+1 \]

其中我们需要新建一个 \(0\) 号位置表示初始在 \(0\) 的位置。

发现绝对值很恶心,考虑拆开,则原转移方程变为两个:

\[f_i=\begin{cases}\max\limits_{j<i,a_i-vt_i\le a_j-vt_j}\{f_j\}+1&a_i\ge a_j\\\max\limits_{j<i,a_j+vt_j\le a_i+vt_i}\{f_j\}+1&a_i<a_j\end{cases} \]

容易发现是三维偏序的转移形式。当然,已经蕴含了转移顺序是按照编号从大到小,所以 \(\rm KDTree\) 上只需要做二维的矩阵查询即可,时间复杂度最坏 \(\mathcal{O}(n\sqrt{n})\)扔个我想不通为啥 \(\tt TLE\ on\ test\ 7\) 的代码在这儿

KD Tree 还能干啥

除了上面的两种用法,别的用法我就不太知道了。但一道题,P4475 巧克力王国,属实把我震撼住了。这道题看不出来上述的任何特征,但它就是能用 \(\rm KDTree\) 水过去,只利用 \(ax+by\) 是关于 \(x,y\) 线性的这一性质。怎么说呢,其实只要我们能利用维护的子树信息判断该子树有无更新答案的可能性,都可以用 \(\rm KDTree\) 的搜索剪枝掉一些子树,从而获得最坏单次 \(\mathcal{O}(n)\),最好单次 \(\mathcal{O}(?)\) 的神奇算法,而且随机数据跑得飞快。

posted @ 2022-05-16 10:45  zhiyangfan  阅读(110)  评论(0编辑  收藏  举报