用于解决 \(k\) 维区间问题的数据结构。
P4148 简单题
这是一种 2-D 树,因为只有两维。
K-D 树用的是类似平衡树的思路,也就是根节点本身也代表一个节点信息,故选择一个维度,用 nth_element 函数找到这个维度居中的点,把它作为根节点即可。
对于一个正方形,横着切一刀竖着切一刀时间复杂度显然优,但是要是对于一个平放的木棍,横着切完全没必要,所以每次按照方差较大的维度进行分治。
插入点之后,不平衡了怎么办?直接用替罪羊树思路暴力重建。具体地,设定一个 \(\alpha=0.75\),若某个点的儿子的子树大小,已经超过了这个点的子树大小的 \(\alpha\) 倍,就把这整个子树的点记录下来,直接重建。
查询的时候,如果是像例题中这种求和,就只能暴力求了,但是如果是求距离最小值 / 最大值,显然可以用各种爆搜用的贪心剪枝,时间复杂度就是 \(O(?)\) 了。
超长代码预警。
namespace KD {
#define mid (l + r >> 1)const double alpha = 0.7;struct TREE {int ls, rs, x, y, v;int x1, y1, x2, y2, s, sz;int d;} t[N];int p[N], np, root, idx;void push_up(int rt) {t[rt].s = t[t[rt].ls].s + t[t[rt].rs].s + t[rt].v;t[rt].sz = t[t[rt].ls].sz + t[t[rt].rs].sz + 1;t[rt].x1 = t[rt].x2 = t[rt].x;t[rt].y1 = t[rt].y2 = t[rt].y;if (t[rt].ls) {t[rt].x1 = min(t[rt].x1, t[t[rt].ls].x1);t[rt].y1 = min(t[rt].y1, t[t[rt].ls].y1);t[rt].x2 = max(t[rt].x2, t[t[rt].ls].x2);t[rt].y2 = max(t[rt].y2, t[t[rt].ls].y2);}if (t[rt].rs) {t[rt].x1 = min(t[rt].x1, t[t[rt].rs].x1);t[rt].y1 = min(t[rt].y1, t[t[rt].rs].y1);t[rt].x2 = max(t[rt].x2, t[t[rt].rs].x2);t[rt].y2 = max(t[rt].y2, t[t[rt].rs].y2);}}void getp(int rt) {if (!rt)return;p[++np] = rt;getp(t[rt].ls), getp(t[rt].rs);}bool cmpx(int x, int y) {return t[x].x < t[y].x;}bool cmpy(int x, int y) {return t[x].y < t[y].y;}void build(int &rt, int l, int r) {if (l > r) {rt = 0;return;}double avx = 0, avy = 0, sx = 0, sy = 0;for (int i = l; i <= r; i++)avx += t[p[i]].x, avy += t[p[i]].y;avx = 1.0 * avx / (r - l + 1), avy = 1.0 * avy / (r - l + 1);for (int i = l; i <= r; i++) {sx += (t[p[i]].x - avx) * (t[p[i]].x - avx);sy += (t[p[i]].y - avy) * (t[p[i]].y - avy);}if (sx > sy) {nth_element(p + l, p + mid, p + r + 1, cmpx);rt = p[mid];t[rt].d = 1;} else {nth_element(p + l, p + mid, p + r + 1, cmpy);rt = p[mid];t[rt].d = 2;}build(t[rt].ls, l, mid - 1);build(t[rt].rs, mid + 1, r);push_up(rt);}void rebuild(int &rt) {if (max(t[t[rt].ls].sz, t[t[rt].rs].sz) <= t[rt].sz * alpha)return;np = 0;getp(rt);build(rt, 1, np);}void update(int &rt, int v) {if (!rt) {rt = v;push_up(rt);return;}if (t[rt].d == 1)t[v].x <= t[rt].x ? update(t[rt].ls, v) : update(t[rt].rs, v);elset[v].y <= t[rt].y ? update(t[rt].ls, v) : update(t[rt].rs, v);push_up(rt);rebuild(rt);}int query(int rt, int x1, int y1, int x2, int y2) {if (!rt or t[rt].x1 > x2 or t[rt].y1 > y2 or t[rt].x2 < x1 or t[rt].y2 < y1)return 0;if (x1 <= t[rt].x1 and t[rt].x2 <= x2 and y1 <= t[rt].y1 and t[rt].y2 <= y2)return t[rt].s;int ret = 0;if (x1 <= t[rt].x and t[rt].x <= x2 and y1 <= t[rt].y and t[rt].y <= y2)ret = t[rt].v;return ret + query(t[rt].ls, x1, y1, x2, y2) + query(t[rt].rs, x1, y1, x2, y2);}
}
using namespace KD;