動態幾何 史蒂芙的泡泡 (解法 2)

在處理完數以百計的政事後，受盡折磨的史蒂芙，打算回家好好地休息。 拖著疲倦的身軀，再也無法再容納任何一點複雜計算。從王宮走回寢居的路上， 發現身邊所見的事物都不再圓滑，看起來就像是粗糙的幾何多邊形構成的一切。

打算享受著泡泡浴的史蒂芙，看著眼前的多邊形泡泡，失去原本應有的色澤，那透涼的心境更蒙上了一層灰影

「為什麼是我呢？」感嘆道

伸出手戳著眼前的泡泡，卻飄了過去

「區區的泡泡也跟我作對，嗚嗚」

將一個泡泡視為一個簡單多邊形 $A$ ，方便起見用一個序列 $a_0, a_1, ..., a_{n-1}$ 表示多邊形 $A$ 的每一個頂點，則會有 $n$ 個線段 $\overline{a_0 a_1}, \overline{a_1 a_2}, \cdots, \overline{a_{n-1} a_0}$ 。

Sample Input

5 14

0 0
10 0 
10 10
0 10
5 5

1 0 0
1 1 0
1 2 1
1 3 2
3 5.5 5
3 10 -1
3 10 5
1 4 1
3 5.5 5
3 10 -1
3 10 5
2 3
3 5 7.5
3 5 2.5

Sample Output

Solution

參閱 動態幾何 史蒂芙的泡泡 (解法 1)

相較於先前的解法 $\mathcal{O}(\log n) - \mathcal{O}(\ast \sqrt{n})$ ，相當不穩定的嵌套 KD-BRH 的解法，實際上如果單純針對這一題，可以拋開 region search 的操作，只有完全的詢問點是否在多邊形內部，則可以做到 $\mathcal{O}(\log^2 n)$ 。

如同一般的射線法，對詢問點拉出無限延長的線，找到與多邊形的邊相交個數。如果單純把每一條邊拿出來看，最暴力的複雜度為 $\mathcal{O}(n)$ ，現在要減少查閱的邊數，且透過 rank tree 在 $\mathcal{O}(\log n)$ 累計相交數量。

由於給定的座標不需要動態，則著手離散化 $X = x_0, x_1, \cdots, x_n$ ，線段樹中的每一個點 $u$ 維護一個 rank tree，把包含者個區段 $[x_l, x_r]$ 的線段通通放入，且保證放入的線段彼此之間不相交 (除了兩端點)。如此一來，當詢問一個點，需要探訪 $\mathcal{O}(\log n)$ 個節點，每個節點 $u$ 需要 $\mathcal{O}(\log n)$ 時間來計算相交數量，最後詢問的複雜度為 $\mathcal{O}(\log^2 n)$ 。同理，修改點的操作也會在 $\mathcal{O}(\log^2 n)$ 。

實作時，特別注意到與射線方向相同的線段不予處理，按照下面的寫法則是不處理垂直的線段，一來是射線法也不會去計算，二來是線段樹劃分的時候會造成一些邊界問題，由於我們對於點離散，父節點 $[x_l, x_r]$ ，左右子樹控制的分別為 $[x_l, x_\text{mid}]$ 、 $[x_\text{mid}, x_r]$ ，劃分時會共用中間點。

即使有了上述的概念來解題，我們仍需要維護 rope data structrue 來維護節點的相鄰關係，可以自己實作任何的 binary tree 來達到 $\mathcal{O}(\log n)$ ，這裡採用 splay tree 示範。

接下來介紹內建黑魔法 PBDS (policy-based data structure) 和 rope。很多人都提及到這些非正式的 STL 函式庫，只有在 gcc/g++ 裡面才有附錄，如果是 clang 等其他編譯器可能是沒有辦法的。所以上傳相關代碼要看 OJ 是否採用純正的 gcc/g++。

參考資料：

• C++ STL: Policy based data structures - Codeforces
• 比STL还STL？——平板电视
• C++ STL: SGI rope

PBDS 比較常用在 rank tree 和 heap，由於代碼量非常多，用內建防止 code length exceed limitation 的出現，也不妨是個好辦法。用 rank tree 的每一個詢問操作皆在 $\mathcal{O}(\log n)$ ，而 heap 選擇 thin heap 或 pairing heap，除了 pop 操作外，皆為 $\mathcal{O}(1)$ ，在對最短路徑問題上別有優勢。

而這一題不太適用 SGI rope，原因在於雖為 $\mathcal{O}(\log n)$ 操作，但它原本就為了仿作 string 而強迫變成可持久化的引數結構，導致每一個操作需要額外的開銷來減少內存使用。由於此題經常在單一元素上操作，SGI rope 對於單一元素效能不彰，以造成嚴重的逾時。

這裡仍提及和示範這些概念的資料結構，哪天正式編入標準函式庫，想必效能問題都已解決。

• KD BRH 0.2s
• PBDS + SPLAY 0.3s
• PBDS + SGI rope 6.0s

把 splay tree 換成 SGI rope 整整慢了二十多倍，測試資料也許還不夠大到適合使用。

PBDS + SPLAY

#include<bits/stdc++.h>
using namespace std;

#include<ext/pb_ds/assoc_container.hpp>
#include<ext/pb_ds/tree_policy.hpp>
using namespace __gnu_pbds;

struct Mesh{
static const int MAXN = 1e5 + 5;
    int pt[MAXN][2];
    vector<int> X;
    void read(int n){
        for (int i = 0; i < n; i++)
            scanf("%d %d", &pt[i][0], &pt[i][1]);
        X.clear(); X.reserve(n);
        for (int i = 0; i < n; i++)
            X.push_back(pt[i][0]);
        sort(X.begin(), X.end());
        X.erase(unique(X.begin(), X.end()), X.end());
    }
} mesh;

class SplayTree{
public:
    struct Node{
        Node *ch[2], *fa;
        int size; int data;
        Node() {
            ch[0] = ch[1] = fa = NULL;
            size = 1;
        }
        bool is_root(){
            return fa->ch[0] != this && fa->ch[1] != this;
        }
    };
    Node *root, *EMPTY;
    void pushdown(Node *u){}
    void pushup(Node *u){
        if (u->ch[0] != EMPTY)	pushdown(u->ch[0]);
        if (u->ch[1] != EMPTY)	pushdown(u->ch[1]);
        u->size = 1 + u->ch[0]->size + u->ch[1]->size;
    }
    void setch(Node *p, Node *u,int i){
        if (p != EMPTY)	p->ch[i] = u;
        if (u != EMPTY)	u->fa = p;
    }
    SplayTree() {
        EMPTY = new Node();
        EMPTY->fa = EMPTY->ch[0] = EMPTY->ch[1] = EMPTY;
        EMPTY->size = 0;
    }
    void init(){
        root = EMPTY;
    }
    Node*newNode(){
        Node *u = new Node();
        u->fa = u->ch[0] = u->ch[1] = EMPTY;
        return u;
    }
    void rotate(Node *x){
        Node *y;
        int d;
        y = x->fa, d = y->ch[1] == x ? 1 : 0;
        x->ch[d^1]->fa = y, y->ch[d] = x->ch[d^1];
        x->ch[d^1] = y;
        if (!y->is_root())
            y->fa->ch[y->fa->ch[1] == y] = x;
        x->fa = y->fa, y->fa = x;
        pushup(y);
    }
    void deal(Node *x){
        if (!x->is_root())	deal(x->fa);
        pushdown(x);
    }
    void splay(Node *x, Node *below){
        if (x == EMPTY)	return ;
        Node *y, *z;
        deal(x);
        while (!x->is_root() && x->fa != below) {
            y = x->fa, z = y->fa;
            if (!y->is_root() && y->fa != below) {
                if (y->ch[0] == x ^ z->ch[0] == y)
                    rotate(x);
                else
                    rotate(y);
            }
            rotate(x);
        }
        pushup(x);
        if (x->fa == EMPTY)	root = x;
    }
    Node*prevNode(Node *u){
        splay(u, EMPTY);
        return maxNode(u->ch[0]);
    }
    Node*nextNode(Node *u){
        splay(u, EMPTY);
        return minNode(u->ch[1]);
    }
    Node*minNode(Node *u){
        Node *p = u->fa;
        for (; pushdown(u), u->ch[0] != EMPTY; u = u->ch[0]);
        splay(u, p);
        return u;
    }
    Node*maxNode(Node *u){
        Node *p = u->fa;
        for (; pushdown(u), u->ch[1] != EMPTY; u = u->ch[1]);
        splay(u, p);
        return u;
    }
    Node*findPos(int pos){ // [0...]
        for (Node *u = root; u != EMPTY;) {
            pushdown(u);
            int t = u->ch[0]->size;
            if (t == pos) {
                splay(u, EMPTY);
                return u;
            }
            if (t > pos)
                u = u->ch[0];
            else
                u = u->ch[1], pos -= t + 1;
        }
        return EMPTY;
    }
    tuple<int, int, int> insert(int data, int pos) {	// make [pos] = data
        Node *p, *q = findPos(pos);
        Node *x = newNode(); x->data = data;
        if (q == EMPTY) {
            p = maxNode(root), splay(p, EMPTY);
            setch(x, p, 0);
            splay(x, EMPTY);
        } else {
            splay(q, EMPTY), p = q->ch[0];
            setch(x, p, 0), setch(x, q, 1);
            setch(q, EMPTY, 0);
            splay(q, EMPTY);
            p = prevNode(x);
        }
        if (p == EMPTY)	p = maxNode(root);
        if (q == EMPTY)	q = minNode(root);
        return make_tuple(p->data, data, q->data);
    }
    tuple<int, int, int> remove(int pos) {
        Node *x = findPos(pos), *p, *q;
        p = prevNode(x), q = nextNode(x);
        if (p != EMPTY && q != EMPTY) {
            setch(p, q, 1);
            p->fa = EMPTY, splay(q, EMPTY);
        } else if (p != EMPTY) {
            p->fa = EMPTY, root = p;
        } else {
            q->fa = EMPTY, root = q;
        }
        int del = x->data;
        free(x);
        if (p == EMPTY)	p = maxNode(root);
        if (q == EMPTY)	q = minNode(root);
        return make_tuple(p->data, del, q->data);
    }
    int size(){
        return root == EMPTY ? 0 : root->size;
    }
} mlist;

struct Pt{
    double x, y;
    Pt() {}
    Pt(int xy[]):Pt(xy[0], xy[1]) {}
    Pt(double x, double y):x(x), y(y) {}
    bool operator<(const Pt &o) const {
        if (x != o.x)	return x < o.x;
        return y < o.y;
    }
};

struct PtP{
    static double x;
    Pt p, q;
    PtP(Pt a, Pt b) {
        p = a, q = b;
        if (q < p)
            swap(p, q);
    }
    double interpolate(const Pt& p1, const Pt& p2, double& x)const {
        if (p1.x == p2.x) return min(p1.y, p2.y);
        return p1.y + (p2.y - p1.y) / (p2.x - p1.x) * (x - p1.x);
    }
    bool operator<(const PtP &o) const {
        return interpolate(p, q, x) < interpolate(o.p, o.q, x);
    }
};
double PtP::x = 1;

struct SegSeg{
    struct Node{
        Node *lson, *rson;
        tree<pair<PtP, int>, null_type, less<pair<PtP, int>>, rb_tree_tag, tree_order_statistics_node_update> segs;
        Node() {
            lson = rson = NULL;
        }
    };
    Node *root;
    int xn;
    Node*newNode(){
        return new Node();
    }
    void freeNode(Node *u){
        free(u);
    }
    void init(){
        root = NULL;
        xn = mesh.X.size();
    }
    void insert(tuple<int,int,int> e){
        int p, q, r; tie(p, q, r) = e;
        remove(0, xn-1, {PtP(Pt(mesh.pt[p]), Pt(mesh.pt[r])), r});
        insert(0, xn-1, {PtP(Pt(mesh.pt[p]), Pt(mesh.pt[q])), q});
        insert(0, xn-1, {PtP(Pt(mesh.pt[q]), Pt(mesh.pt[r])), r});
    }
    void remove(tuple<int,int,int> e){
        int p, q, r; tie(p, q, r) = e;
        remove(0, xn-1, {PtP(Pt(mesh.pt[p]), Pt(mesh.pt[q])), q});
        remove(0, xn-1, {PtP(Pt(mesh.pt[q]), Pt(mesh.pt[r])), r});
        insert(0, xn-1, {PtP(Pt(mesh.pt[p]), Pt(mesh.pt[r])), r});
    }	
    int inside(double x, double y){
        PtP::x = x;
        return count(root, 0, xn-1, x, y)&1;
    }
    int count(Node* u,int l, int r, double x, double y){
        if (u == NULL)
            return 0;
        int ret = 0;
        if ((mesh.X[l] > x) != (mesh.X[r] > x))
            ret += u->segs.order_of_key({PtP(Pt(x, y), Pt(x, y)), -1});
        int m = (l+r)/2;
        if (x <= mesh.X[m])
            ret += count(u->lson, l, m, x, y);
        if (x >= mesh.X[m])
            ret += count(u->rson, m, r, x, y);
        return ret;
    }
    void insert(int l, int r, pair<PtP, int> s){
        if (s.first.p.x != s.first.q.x)
            insert(root, l, r, s);
    }
    void remove(int l, int r, pair<PtP, int> s){
        if (s.first.p.x != s.first.q.x)
            remove(root, l, r, s);
    }
    void insert(Node* &u,int l, int r, pair<PtP, int> s){
        if (u == NULL)
            u = newNode();
        if (s.first.p.x <= mesh.X[l] && mesh.X[r] <= s.first.q.x) {
            PtP::x = (mesh.X[l] + mesh.X[r])/2.0;
            u->segs.insert(s);
            return;
        }
        int m = (l+r)/2;
        if (s.first.q.x <= mesh.X[m])			insert(u->lson, l, m, s);
        else if (s.first.p.x >= mesh.X[m])	insert(u->rson, m, r, s);
        else	insert(u->lson, l, m, s), insert(u->rson, m, r, s);
    }
    void remove(Node* u,int l, int r, pair<PtP, int> s){
        if (u == NULL)
            return;
        if (s.first.p.x <= mesh.X[l] && mesh.X[r] <= s.first.q.x) {
            PtP::x = (mesh.X[l] + mesh.X[r])/2.0;
            u->segs.erase(s);
            return;
        }
        int m = (l+r)/2;
        if (s.first.q.x <= mesh.X[m])			remove(u->lson, l, m, s);
        else if (s.first.p.x >= mesh.X[m])	remove(u->rson, m, r, s);
        else	remove(u->lson, l, m, s), remove(u->rson, m, r, s);
    }
} mbrh;

int main(){
    int n, m, cmd, x, pos;
    double px, py;
    scanf("%d %d", &n, &m);
    mesh.read(n);
    mlist.init(), mbrh.init();
    for (int i = 0; i < m; i++) {
        scanf("%d", &cmd);
        if (cmd == 1) {
            scanf("%d %d", &x, &pos);
            mbrh.insert(mlist.insert(x, pos));
        } else if (cmd == 2) {
            scanf("%d", &x);
            mbrh.remove(mlist.remove(x));
        } else {
            scanf("%lf %lf", &px, &py);
            puts(mbrh.inside(px, py) ? "1" : "0");
        }
    }
    return 0;
}

PBDS + ROPE

#include<bits/stdc++.h>
using namespace std;

#include<ext/rope>
using namespace __gnu_cxx;

#include<ext/pb_ds/assoc_container.hpp>
#include<ext/pb_ds/tree_policy.hpp>
using namespace __gnu_pbds;

struct Mesh{
static const int MAXN = 1e5 + 5;
    int pt[MAXN][2];
    vector<int> X;
    void read(int n){
        for (int i = 0; i < n; i++)
            scanf("%d %d", &pt[i][0], &pt[i][1]);
        X.clear(); X.reserve(n);
        for (int i = 0; i < n; i++)
            X.push_back(pt[i][0]);
        sort(X.begin(), X.end());
        X.erase(unique(X.begin(), X.end()), X.end());
    }
} mesh;

class Rope{
    rope<int> r;
public:
    void init(){
        r.clear();
    }
    tuple<int, int, int> insert(int data, int pos) {	// make [pos] = data
        assert (pos <= r.size());
        r.insert(pos, data);
        int p = r.mutable_reference_at((pos-1+r.size())%r.size());
        int q = r.mutable_reference_at((pos+1+r.size())%r.size());
        return make_tuple(p, data, q);
    }
    tuple<int, int, int> remove(int pos) {
        assert(pos < r.size());
        int del = r.at(pos);
        int p, q;
        if (r.size() == 1) {
            p = q = del;
        } else {
            assert(r.size() > 0);
            p = r.at((pos-1+r.size())%r.size());
            q = r.at((pos+1+r.size())%r.size());
        }
        r.erase(pos, 1);
        return make_tuple(p, del, q);
    }
} mlist;

struct Pt{
    double x, y;
    Pt() {}
    Pt(int xy[]):Pt(xy[0], xy[1]) {}
    Pt(double x, double y):x(x), y(y) {}
    bool operator<(const Pt &o) const {
        if (x != o.x)	return x < o.x;
        return y < o.y;
    }
};

struct PtP{
    static double x;
    Pt p, q;
    PtP(Pt a, Pt b) {
        p = a, q = b;
        if (q < p)
            swap(p, q);
    }
    double interpolate(const Pt& p1, const Pt& p2, double& x)const {
        if (p1.x == p2.x) return min(p1.y, p2.y);
        return p1.y + (p2.y - p1.y) / (p2.x - p1.x) * (x - p1.x);
    }
    bool operator<(const PtP &o) const {
        return interpolate(p, q, x) < interpolate(o.p, o.q, x);
    }
};
double PtP::x = 1;

struct SegSeg{
    struct Node{
        Node *lson, *rson;
        tree<pair<PtP, int>, null_type, less<pair<PtP, int>>, rb_tree_tag, tree_order_statistics_node_update> segs;
        Node() {
            lson = rson = NULL;
        }
    };
    Node *root;
    int xn;
    Node*newNode(){
        return new Node();
    }
    void freeNode(Node *u){
        free(u);
    }
    void init(){
        root = NULL;
        xn = mesh.X.size();
    }
    void insert(tuple<int,int,int> e){
        int p, q, r; tie(p, q, r) = e;
        remove(0, xn-1, {PtP(Pt(mesh.pt[p]), Pt(mesh.pt[r])), r});
        insert(0, xn-1, {PtP(Pt(mesh.pt[p]), Pt(mesh.pt[q])), q});
        insert(0, xn-1, {PtP(Pt(mesh.pt[q]), Pt(mesh.pt[r])), r});
    }
    void remove(tuple<int,int,int> e){
        int p, q, r; tie(p, q, r) = e;
        remove(0, xn-1, {PtP(Pt(mesh.pt[p]), Pt(mesh.pt[q])), q});
        remove(0, xn-1, {PtP(Pt(mesh.pt[q]), Pt(mesh.pt[r])), r});
        insert(0, xn-1, {PtP(Pt(mesh.pt[p]), Pt(mesh.pt[r])), r});
    }	
    int inside(double x, double y){
        PtP::x = x;
        return count(root, 0, xn-1, x, y)&1;
    }
    int count(Node* u,int l, int r, double x, double y){
        if (u == NULL)
            return 0;
        int ret = 0;
        if ((mesh.X[l] > x) != (mesh.X[r] > x))
            ret += u->segs.order_of_key({PtP(Pt(x, y), Pt(x, y)), -1});
        if (l == r)
            return ret;
        int m = (l+r)/2;
        if (x <= mesh.X[m])
            ret += count(u->lson, l, m, x, y);
        if (x >= mesh.X[m])
            ret += count(u->rson, m, r, x, y);
        return ret;
    }
    void insert(int l, int r, pair<PtP, int> s){
        if (s.first.p.x != s.first.q.x)
            insert(root, l, r, s);
    }
    void remove(int l, int r, pair<PtP, int> s){
        if (s.first.p.x != s.first.q.x)
            remove(root, l, r, s);
    }
    void insert(Node* &u,int l, int r, pair<PtP, int> s){
        if (u == NULL)
            u = newNode();
        if (s.first.p.x <= mesh.X[l] && mesh.X[r] <= s.first.q.x) {
            PtP::x = (mesh.X[l] + mesh.X[r])/2.0;
            u->segs.insert(s);
            return;
        }
        if (l == r)
            return;
        int m = (l+r)/2;
        if (s.first.q.x <= mesh.X[m])			insert(u->lson, l, m, s);
        else if (s.first.p.x >= mesh.X[m])	insert(u->rson, m, r, s);
        else	insert(u->lson, l, m, s), insert(u->rson, m, r, s);
    }
    void remove(Node* u,int l, int r, pair<PtP, int> s){
        if (u == NULL)
            return;
        if (s.first.p.x <= mesh.X[l] && mesh.X[r] <= s.first.q.x) {
            PtP::x = (mesh.X[l] + mesh.X[r])/2.0;
            u->segs.erase(s);
            return;
        }
        if (l == r)
            return;
        int m = (l+r)/2;
        if (s.first.q.x <= mesh.X[m])			remove(u->lson, l, m, s);
        else if (s.first.p.x >= mesh.X[m])	remove(u->rson, m, r, s);
        else	remove(u->lson, l, m, s), remove(u->rson, m, r, s);
    }
} mbrh;

int main(){
    int n, m, cmd, x, pos;
    double px, py;
    scanf("%d %d", &n, &m);
    mesh.read(n);
    mlist.init(), mbrh.init();
    for (int i = 0; i < m; i++) {
        scanf("%d", &cmd);
        if (cmd == 1) {
            scanf("%d %d", &x, &pos);
            mbrh.insert(mlist.insert(x, pos));
        } else if (cmd == 2) {
            scanf("%d", &x);
            mbrh.remove(mlist.remove(x));
        } else {
            scanf("%lf %lf", &px, &py);
            puts(mbrh.inside(px, py) ? "1" : "0");
        }
    }
    return 0;
}