LinkCut木を実装してみた

LinkCut木について勉強する機会があったので@iwiwi先生のスライドと,Tarjanの本(邦訳は絶版っぽい.そして原著高いので誰か買ってください)を参考に実装してみた.

LinkCut木は,木をパスに分解して木にして,それに含まれるパスをSplay木で管理するというよくわからない構造のデータ構造.
いろいろな操作を O ( \log n )ぐらいで行える...らしい...

実装した機能は以下の通り

link(c,v)
根をcとする木をvに接続する(ただしcは必ず根)
cut(c)
木をcから親に向かう辺を削除して2つの木に分ける
findroot(v)
vを持つ木の根を見つける
addcost(v,x)
vから根に向かうパスの中にある頂点の重みにxを加える.
findcost(v)
vから根に向かうパスの中で最小の重みと,その重みをもつ最も根に近い頂点の組を返す.
evert(v)
vを持つ木の根をvに変更する

iwi先生の言う通りに実装したらSplay木の一番左のノードが示す頂点が根に最も近い頂点を表すことになった.Tarjanの本では逆だった気がする.

link,cut,evertはiwi大先生がスライドにソースを載せていらっしゃるのでそっちを参照したらかける.
問題は値の管理の方で,Tarjanの本にはさらっと書いてある割に考えることが多い.

cost(x)は自分が持つコスト,mincost(x)は自分とその子が持つコストの最小値を表す.Splay木のノードに以下で定義する値を持たせる.

<br />
\Delta cost(x) = cost(x) - mincost(x) \\<br />
\Delta min(x) = \left\{<br />
mincost(x) \text{    ($x$ is root of splay tree)} \\<br />
mincost(x) - mincost(parent(x)) \text{    (otherwise)}<br />
\right<br />

この値をノードを回転させたり,exposeしたり,linkしたり,cutしたりするときにがんばって更新してあげないといけない.
操作を行う前のノードの値を \Delta cost (x), \Delta min (x),操作を行った後のノードの値を \Delta cost^\prime (x), \Delta min^\prime (x)で表す.

初めに回転を行った時の更新を考える.
回転は右と左があるが,片方について考えれば他方は向きを反転して適用してやればいいので,右回転についてのみ考える.
現在参照しているSplayTreeのノードをp,その親をq,pの左子供をa,右子供をb,qの右子供をcとする.

f:id:blue_jam:20130627183653p:plain

この時,次の等式が成り立つ.

\Delta min^\prime(p)=\Delta min(q)
\Delta cost^\prime(p)=\Delta cost(p)+\Delta min(p)
\Delta min^\prime(q)=\min(\Delta cost(q), \Delta min(c), \Delta min(p) + \Delta min(a))
\Delta cost(q) + \Delta min(q)=\Delta min^\prime (p) + \Delta min^\prime(q) + \Delta cost^\prime(q)
\Delta min^\prime (a)=\Delta min(a)+\Delta min(p)
\Delta min (b) + \Delta min(p) = \Delta min^\prime (b) + \Delta min^\prime (q)
\Delta min^\prime (c) = \Delta min^\prime (c) + \Delta min^\prime (q)

次にexposeをした後の更新を考える.
exposeのコード(構造に関するところだけ)を示す.

node *expose(node *x){
    node *rp = NULL;
    for(node *p = x; p != NULL; p = p->pp){
        p->splay();
        p->cp[R] = rp;
        rp = p;
    }
    x->splay();
    return x;
}

f:id:blue_jam:20130627183638p:plain

この時,次の等式が成り立つ.
 \Delta min^\prime (b) = \Delta min(b) + \Delta min(p)
 \Delta min^\prime (p) = \min (\Delta min(p) + \Delta cost(p), \Delta min(p)+\Delta min(a), \Delta min(rp))
 \Delta min^\prime (rp) = \Delta min(rp) + \Delta min^\prime(p)
 \Delta min^\prime (a) + \Delta min^\prime(p) = \Delta min(a) + \Delta min(p)
 \Delta cost^\prime (p) + \Delta min^\prime (p) = \Delta cost(p) + \Delta min(p)

次にlink(c, p)をした後の更新を考える.
引数の前提よりcは根であるので,expose(c)を行った後,cを表すSplayTreeのノードは子を持たない.
また,expose(p)を行っているので,pを表すSplayTreeのノードは右子供を持たない.

f:id:blue_jam:20130627183648p:plain

この時次の等式が成り立つ.
 \Delta min^\prime (p) = \min(\Delta min(p), \Delta min(c))
 \Delta min^\prime (c) = \Delta min(c) - \Delta min^\prime(p)
 \Delta min^\prime (a) + \Delta min^\prime(p) = \Delta min(a) + \Delta min(p)
 \Delta cost^\prime (p) + \Delta min^\prime(p) = \Delta cost(p) + \Delta min(p)

最後にcut(c)をした後の更新を考える.
expose(c)を行った後では,cを表すSplayTreeのノードは右子供を持たない.

f:id:blue_jam:20130627183622p:plain

この時,次の等式が成り立つ.

 \Delta min^\prime(p) = \Delta min(c) + \Delta min(p)
 \Delta min^\prime(c) = \Delta min(c) + \Delta cost(c)
 \Delta cost^\prime (c) = 0

以上をプログラム中にうまく書いてやれば正しく値の更新ができる.
findrootやfindcost,addcostの実装方法はTarjanの本を参照するか推して知るべし.

ソースの全文は次のようになる.

struct LinkCutTree{
    static const int L = 0, R = 1;
    struct node{
        node *pp, *cp[2];//*lp, *rp;
        bool rev;
        int dcost, dmin;
        int id;
        node(){
            id = -1;
            pp = cp[L] = cp[R] = NULL;
            rev = false;
            dcost = dmin = 0;
        }
        // Splay木の根であるか
        bool isRoot(){
            return !pp || pp->cp[L] != this && pp->cp[R] != this;
        }
        // ノードの状態(反転など)を子に伝搬させる
        void push(){
            if(rev){
                swap(cp[L], cp[R]);
                if(cp[L]) cp[L]->rev ^= true;
                if(cp[R]) cp[R]->rev ^= true;
                rev = false;
            }
        }
        // dがRなら右回転,dがLなら左回転を行う
        void rot(int d){
            node *q = pp, *r = q->pp;
            int e = 1 - d;
            if(q->cp[e] = cp[d]) cp[d]->pp = q;
            cp[d] = q; q->pp = this;

            int qd = q->dmin;
            q->dmin = q->dcost;
            if(q->cp[e]) q->dmin = min(q->dmin, dmin + q->cp[e]->dmin);
            if(q->cp[d]) q->dmin = min(q->dmin, q->cp[d]->dmin);
            q->dcost -= q->dmin;
            if(q->cp[e]) q->cp[e]->dmin += dmin - q->dmin;
            if(q->cp[d]) q->cp[d]->dmin -= q->dmin;
            if(cp[e]) cp[e]->dmin += dmin;
            dcost += dmin;
            dmin = qd;
            
            if(pp=r){
                if(r->cp[L] == q) r->cp[L] = this;
                if(r->cp[R] == q) r->cp[R] = this;
            }
        }
        // Splay操作
        void splay(){
            push();
            while(!isRoot()){
                node *q = pp;
                if(q->isRoot()){
                    q->push(); push();
                    if(q->cp[L] == this) rot(R);
                    else rot(L);
                }
                else{
                    node *r = q->pp;
                    r->push(); q->push(); push();
                    for(int i = 0; i < 2; ++i) if(r->cp[i] == q){
                        if(q->cp[i] == this){ q->rot(1-i); rot(1-i); }
                        else{ rot(i); rot(1-i); }
                        break;
                    }
                }
            }
        }
    };
    vector<node> nodes;
    LinkCutTree(int n) : nodes(n) {
        for(int i = 0; i < n; ++i) nodes[i].id = i;
    }
    node *expose(node *x){
        node *rp = NULL;
        for(node *p = x; p != NULL; p = p->pp){
            int ndmin;
            p->splay();

            if(p->cp[R]) p->cp[R]->dmin += p->dmin;
            ndmin = p->dmin + p->dcost;
            if(rp) ndmin = min(ndmin, rp->dmin);
            if(p->cp[L]) ndmin = min(ndmin, p->dmin + p->cp[L]->dmin);
            if(rp) rp->dmin = rp->dmin - ndmin;
            if(p->cp[L]) p->cp[L]->dmin += p->dmin - ndmin;
            p->dcost += p->dmin - ndmin;
            p->dmin = ndmin;

            p->cp[R] = rp;
            rp = p;
        }
        x->splay();
        return x;
    }
    int findroot(int n){
        node *v = expose(&nodes[n]);
        while(v->cp[L]) v->push(), v = v->cp[L];
        v->splay();
        return v->id;
    }
    pair<int, int> findcost(int n){
        node *v = expose(&nodes[n]);
        for(;;){
            v->push();
            if(v->cp[L] && v->cp[L]->dmin == 0) v = v->cp[L];
            else if((!v->cp[L] || v->cp[L]->dmin > 0) && v->dcost > 0) v = v->cp[R];
            else break;
        }
        v->splay();
        return make_pair(v->dmin, v->id);
    }
    void addcost(int n, int x){
        node *v = expose(&nodes[n]);
        v->dmin += x;
    }
    void cut(int v){
        node *c = &nodes[v];
        expose(c);
        node *p = c -> cp[L];
        c->cp[L] = NULL;
        p->pp = NULL;

        p->dmin += c->dmin;
        c->dmin += c->dcost;
        c->dcost = 0;
    }
    void link(int v, int w){
        node *c = &nodes[v], *p = &nodes[w];
        expose(c);
        expose(p);
        c->pp = p;
        p->cp[R] = c;

        p->dcost += p->dmin - min(p->dmin, c->dmin);
        if(p->cp[L]) p->cp[L]->dmin += p->dmin - min(p->dmin, c->dmin);
        p->dmin = min(p->dmin, c->dmin);
        c->dmin = c->dmin - p->dmin;
    }
    void evert(int v){
        node *r = expose(&nodes[v]);
        r->rev ^= true;
    }
};

vectorで管理しているノードのポインタを扱っているので,vectorのサイズが変わった時にほぼ確実に変な挙動を示すはず.なので,大きさは絶対に変えないこと.どうしても変えたい場合はlistにするか,ポインタの代わりにvectorのindexを使って実装してください.
オンラインジャッジにLinkCut木を使う問題があるかどうか知らないのでVerifyはしてない.*1
代わりに,link,cut,findroot,addcost,findcost,evertを愚直に実装したものを用意して,ランダムで選んだ操作に対して正しく動くか一晩かけて試したところ,問題なく動いたのでたぶん大丈夫.
でも,これを使って何か被害があったとしても僕は責任を負わない.

何か間違いがあったら教えてください.

*1:@gachizei_tcさんにJOI2013年春合宿4日目のSpaceShipsを紹介していただき,AtCoderにてlink,cut,findrootをVerifyしました.