k-D Tree学习笔记

某人奶了一口省选会考，于是认真学了学。

形态

k-D Tree 是一棵二叉搜索树，每个结点都对应 $k$ 维空间中的一个点，每个子树中的点都在且都仅在一个 $k$ 维超长方体内。

构建

假设已知 $k$ 维空间中的 $n$ 个互不相同的点的坐标，要将其构建为一棵 k-D Tree。

1. 当前超长方体中只有一个点时直接返回该点。
2. 选择一个维度与该维度上的一个切割点，将当前超长方体按这个维度及这个切割点分为两个超长方体。
3. 将选择的切割点作为该子树的根，递归处理分出的两个超长方体。

容易发现，对于步骤 2，每次切割点选为该维度上的中位数，树高可以在 $O(\log n)$。

如何寻找中位数？可以使用类似快速排序的方法，并且每次只递归包含中位数的一侧，对于 $n$ 个元素的序列期望时间复杂度即为 $O(n)$。（可用algorithm库中的nth_element()实现）

每次选什么维度来割？理论上应当每次选方差最大的维度割常数应当最小，然而算法竞赛中较常用且方便的写法是直接轮流按每一维割。

插入与删除

点集可变，可以利用替罪羊树不平衡重构的方法把树高稳定在 $O(\log n)$。

然而实验证明定期整树重构通常来说常数更小。

查询超长方体

若该结点完全被包含于超长方体内，或完全与其无交可直接返回。

否则应递归两边。

可以证明复杂度最坏为 $O(n^{1-\frac{1}{k}})$。

习题

[Luogu 4148]简单题

初值全为 0 的 $n\times n$ 二维矩阵上， $q$ 个操作：

1. 单点修改。
2. 矩形查询和。

强制在线，内存限制20M，答案及过程量均int范围内。$n\le 5\times 10^5, q\le 2\times 10^5$。


空间卡树套树，强制在线卡 CDQ 分治，于是只能 KDT。

不要不平衡重构，定期重构常数更小。

复制

#include <bits/stdc++.h>
namespace my_std {
using namespace std;
#define reg register
#define Rint register int
#define FOR(i, a, b) for (register int i = (a), ed_##i = (b); i <= ed_##i; ++i)
#define ROF(i, a, b) for (register int i = (a), ed_##i = (b); i >= ed_##i; --i)
#define FORit(templ, arr, i, a, b) \
for (register templ *i = (arr) + (a), *ed_##i = (arr) + (b) + 1; i != ed_##i; ++i)
#define ROFit(templ, arr, i, a, b) \
for (register templ *i = (arr) + (a), *ed_##i = (arr) + (b)-1; i != ed_##i; --i)
#define GO(x, p, e, i, v) for (register int i = p[x].head, v; i; i = e[i].link)
#define MEM(x, v) memset(x, v, sizeof(x))
#define fir first
#define sec second
#define pq priority_queue
#define PB push_back
#define MP make_pair
typedef long long LL;
typedef double DB;
typedef pair<int, int> PII;
#define Templ(T) template <typename T>
inline int read() {
    reg int ans = 0, f = 1;
    reg char c = getchar();
    while (!isdigit(c)) f ^= (c == '-'), c = getchar();
    for (; isdigit(c); c = getchar()) ans = (ans << 1) + (ans << 3) + (c ^ 48);
    return f ? ans : -ans;
}
Templ(_Tp) inline int chkmin(_Tp &x, _Tp y) { return x > y ? x = y, 1 : 0; }
Templ(_Tp) inline int chkmax(_Tp &x, _Tp y) { return x < y ? x = y, 1 : 0; }
#define using_mod
const int mod = 998244353;
#ifdef using_mod
inline void inc(int &x, const int &y) { x += y; if (x >= mod) x -= mod; }
inline void dec(int &x, const int &y) { x -= y; if (x < 0) x += mod; }
inline int ksm(int x, LL y) {
    int res = 1;
    for (; y; y >>= 1, x = 1ll * x * x % mod)
        if (y & 1) res = 1ll * res * x % mod;
    return res;
}
inline int qmo(const int &x) { return x + ((x >> 31) ? mod : 0); }
#endif
Templ(_Tp) inline _Tp gcd(_Tp x, _Tp y) { return y ? gcd(y, x % y) : x; }
#define FILE(s) freopen(s ".in", "r", stdin), freopen(s ".out", "w", stdout)
#define PBTXDY
} // namespace my_std
using namespace my_std;

const int N = 200010;
struct Node{
    int e[2][2], b[2], c[2], val, sum, siz;
};

namespace KDT{
Node a[N];
static const DB alpha = 0.75;
int ldr[N], ccnt, root;
static int d_d;
inline bool c_d(const Node &l, const Node &r){
    return l.b[d_d] != r.b[d_d] ? l.b[d_d]
        < r.b[d_d] : l.b[d_d ^ 1] < r.b[d_d ^ 1];
}
#define lc a[t].c[0]
#define rc a[t].c[1]
inline void Print(const int &t){
    printf("Node[%d]: --- <%d, %d>\n", t, lc, rc);
    printf("[%d, %d, %d], [%d, %d, %d], %d, %d, %d\n",
        a[t].e[0][0], a[t].b[0], a[t].e[0][1],
        a[t].e[1][0], a[t].b[1], a[t].e[1][1],
        a[t].val, a[t].sum, a[t].siz);
}
inline void PushUp(const int &t){
    a[t].siz = a[lc].siz + a[rc].siz + 1;
    a[t].sum = a[lc].sum + a[rc].sum + a[t].val;
    if(lc){
        chkmin(a[t].e[0][0], a[lc].e[0][0]);
        chkmax(a[t].e[0][1], a[lc].e[0][1]);
        chkmin(a[t].e[1][0], a[lc].e[1][0]);
        chkmax(a[t].e[1][1], a[lc].e[1][1]);
    }
    if(rc){
        chkmin(a[t].e[0][0], a[rc].e[0][0]);
        chkmax(a[t].e[0][1], a[rc].e[0][1]);
        chkmin(a[t].e[1][0], a[rc].e[1][0]);
        chkmax(a[t].e[1][1], a[rc].e[1][1]);
    }
}
inline bool CanRbu(const int &t){
    return alpha * a[t].siz <= (DB)max(a[lc].siz, a[rc].siz);
}
inline void Build(const int &l, const int &r, const int &d){
    if(l > r) return;
    const int t = (l + r) >> 1;
    d_d = d, nth_element(a + l, a + t, a + r + 1, c_d);
    a[t].e[0][0] = a[t].e[0][1] = a[t].b[0];
    a[t].e[1][0] = a[t].e[1][1] = a[t].b[1];
    a[t].siz = 1, a[t].sum = a[t].val;
    Build(l, t - 1, d ^ 1), Build(t + 1, r, d ^ 1);
    lc = l > t - 1 ? 0 : (l + t - 1) >> 1;
    rc = t + 1 > r ? 0 : (t + r + 1) >> 1;
    return PushUp(t);
}
//inline void Flatten(int &ldc, const int &t){
//    if(!t) return;
//    Flatten(ldc, lc);
//    ldr[ldc++] = t;
//    Flatten(ldc, rc);
//}
//inline int Build(const int &l, const int &r, const int &d){
//    if(l >= r) return 0;
//    const int mid = (l + r) >> 1;
//    d_d = d, nth_element(ldr + l, ldr + mid, ldr + r, c_d);
//    const int t = ldr[mid];
//    lc = Build(l, mid, d ^ 1), rc = Build(mid + 1, r, d ^ 1);
//    return PushUp(t), t;
//}
//inline void Rebuild(int &t, const int &d){
////    printf("Rebuild: %d\n", t);
//    int ldc = 0;
//    Flatten(ldc, t);
//    t = Build(0, ldc, d);
//}
inline void Insert(int &t, const Node &v, const int &d){
    if(!t){
        t = ++ccnt;
        a[t] = v;
        a[t].e[0][0] = a[t].e[0][1] = a[t].b[0];
        a[t].e[1][0] = a[t].e[1][1] = a[t].b[1];
        return;
    }
    if(v.b[0] == a[t].b[0] && v.b[1] == a[t].b[1]){
        a[t].val += v.val;
        a[t].sum += v.sum;
        return;
    }
    if(v.b[d] <= a[t].b[d]) Insert(lc, v, d ^ 1);
    else Insert(rc, v, d ^ 1);
    PushUp(t);
//    if(CanRbu(t)) Rebuild(t, d);
}
inline int Query(const int &t, const Node &v){
    if(!t) return 0;
//    Print(t);
    if(v.e[0][0] <= a[t].e[0][0] && a[t].e[0][1] <= v.e[0][1] &&
        v.e[1][0] <= a[t].e[1][0] && a[t].e[1][1] <= v.e[1][1]){
        return a[t].sum;
    }
    if(v.e[0][0] > a[t].e[0][1] || a[t].e[0][0] > v.e[0][1] ||
        v.e[1][0] > a[t].e[1][1] || a[t].e[1][0] > v.e[1][1]){
        return 0;
    }
    int res = 0;
    if(v.e[0][0] <= a[t].b[0] && a[t].b[0] <= v.e[0][1] &&
        v.e[1][0] <= a[t].b[1] && a[t].b[1] <= v.e[1][1]){
        res = a[t].val;
    }
    return res + Query(lc, v) + Query(rc, v);
}
#undef lc
#undef rc
}

int n;

int main() {
    n = read();
    Rint opt, lst = 0, tot = 0;
    reg const int blo = 6666;
    reg Node v;
    v.siz = 1, v.c[0] = v.c[1] = 0;
    while(1){
        opt = read();
        if(opt == 3) break;
        if(opt == 1){
            ++ tot;
            v.b[0] = read() ^ lst, v.b[1] = read() ^ lst;
            v.val = v.sum = read() ^ lst;
            KDT::Insert(KDT::root, v, 0);
            if(tot == blo){
                tot = 0;
                KDT::Build(1, KDT::ccnt, 0);
                KDT::root = (1 + KDT::ccnt) >> 1;
            }
        }
        else{
            v.e[0][0] = read() ^ lst, v.e[1][0] = read() ^ lst,
            v.e[0][1] = read() ^ lst, v.e[1][1] = read() ^ lst;
//            printf("Qn: [%d, %d], [%d, %d]\n",
//                v.e[0][0], v.e[0][1], v.e[1][0], v.e[1][1]);
            printf("%d\n", lst = KDT::Query(KDT::root, v));
        }
    }
    return 0;
}

参考资料

OI Wiki