Toptree 導入編
みんな日本語記事を待っていたはず....!
toptreeがどんな感じで動いているのかを書いてみます
実装はここにあります
https://github.com/niuez/toptree-rust
0. toptree is なに
toptreeはlink-cut treeの上位互換です. 木を切ったりつなげたり, パスのクエリを処理したり, 木上の二分探索ができたりします
今回はそのベースとなる構造の話です
1. Compress Rake
木をまとめる
ここで言う木は, toptreeが表す木のことです. 曖昧にならないようにこのことをreal treeと呼ぶことにします. (木を木で表現するの文章が曖昧になりがち)
toptreeでは, つながっている2つの辺をまとめる操作を繰り返したものを表現した木です. 1つの辺にまとめ上げることでパスを表現します.
まとめる操作は２つあり, それぞれCompress, Rakeといいます.
具体的に
こんな感じでまとめていきます
このまとめていく操作を木で表現するのがtoptreeです.
Compress Tree
例えば, a - c - bという一直線のreal treeを扱う時, 辺acとcbをcompressをします.
これをtoptreeで表現すると
四角は辺を表すノードです. Edge Nodeといいます. toptreeでは, Edge Nodeを葉にします.
丸はcompressした後の辺を表すノードです. Compress Nodeといいます.
Compress Nodeが節, Edge Nodeが葉のこの木をCompress Treeといいます.
重要なのは, Compress Treeの根がパスの端点を結ぶ, Compressされた辺を表しているということです. またそのcompressされた辺の端点は必ず次数が1です.
0 - 1 - 2 - 3 - 4 - 5 というreal treeをtoptreeで表すと, 一例としては以下のようになります.
辺の向き付け
ここで葉のEdge Nodeの順番がパスの辺の順番になっている点に注意してください.
toptreeでは辺の向きに注意して操作しないとダメです.
僕の実装では, 0 - 1 - 2 - 3 - 4 - 5のtoptreeを0 -> 1 -> 2 -> 3 -> 4 -> 5と向き付けると解釈しています. 以後のreal treeの図では向き付けしたものを用います.
compress, rakeを, 同じ向きのものをまとめる操作と解釈することにしましょう. すると, compress, rakeを向き付けたreal treeについて改めて考えると以下のようになります.
上の具体的にで示したreal treeはこんな感じで向き付けすると同じようにまとめる操作ができるはずです.
Rake Tree
では一直線ではないreal tree, 例えばこれはどうやってtoptreeにするのでしょうか.
ここでrakeを使います.
?????????????????????
辺14を追加したreal treeをtoptreeにすると
ひし形のノードはcompressと同じように察せるはずです. 31と41をrakeしたものを表現しており, ひし形のノードをRake Nodeといいます.
また, Rake Nodeが節, Compress Treeの根が葉の木をRake Treeといいます.
図では, Compress Nodeに今までの左右の子と, 赤の線でつながった子があります.
赤の線でつながった子は, Compress Nodeの右の子とrakeされるノードです.
このように3分木にしてreal treeの情報を持ちます.
具体例のtoptree
こんな感じになります
四角はEdge Node, 丸はCompress Node, ひし形はRakeNodeです.
Compress Tree(青の点線で囲った部分), Rake Tree(赤の点線で囲った部分)はそれぞれここです
Splice
重要なのは, Compress Treeの根がパスの端点を結ぶ, Compressされた辺を表しているということです.
この木のtoptreeをもう一度見てみます.
0 - 1 - 2がこのtoptreeの主役のパスになっています.
でも, 3 - 1 - 2を主役にしたいときもあるはずです. それは, 31と01を入れ替えることで達成できます.
0 - 1 - 3にを主役にしたいときもあるはずです. それは, 31と12の向きを反転させてから入れ替えることで達成できます.
この, Rake Treeの葉のノードと, Compress nodeの子を入れ替えてCompress Nodeの表すパスを変える操作をspliceといいます.
Splay
Splay Treeを知っていますか? wikipediaを見て
Splay Treeでは, splayという木を回転させてノードを根に持ってくるという操作をします. まあwikipediaみて
toptreeで扱っている木は, 葉木です. 葉はsplayできないことに注意しましょう.
Handle
spliceをするとパスを変形できることはわかりましたが, 具体的にどのノードをspliceすると良いのでしょうか?
それを示すのがhandleという概念です. handleは各頂点に対して割り振られるもので, toptree上のCompress/Edge Nodeが割り振られます.
具体的には, 下のルールで構成します.
- Compress Node
abの左右の子がac,cbのとき, 頂点cのhandleはCompress Nodeab - Compress/Edge Node
abの親が - いない(toptreeの根): 頂点
a,bのhandleはCompress Nodeab - それ以外(Rake Treeの葉になっている): 頂点
aのhandleはCompress Nodeab
具体例を見たほうが早い気がします.
頂点0, 5はルール2., 8, a, b, cはルール3., それ以外は1.です.
今はとりあえずこういうものとしておくのがいいと思います.(あとで大活躍します.)
また頂点vに対して, N_vをvのhandleのNodeとします. (上の図で言えばN_2はtoptreeの根です)
Expose
exposeという操作を導入したいと思います. (これが超本質)
任意の頂点vのhandleをtoptreeの根にするのがexposeです.
先にどうやってexposeするか書いてしまいます.
expose(v)
1. N_vをCompress Tree上でsplayする
2. N_vの親が
- いない: N_vはtoptreeの根になったのでexpose終了
- Compress Node: そのCompress Nodeをnとおく
- Rake Node: そのRake Nodeをrとおく, rをRake Tree上でsplayし, rの親をnとおく(nはCompress Nodeになります)
3. nをCompress Tree上でsplay
4. nの左のノードとN_vを入れ替える
5. N_vがEdge Nodeのとき, N_vをnにする
6. 1に戻る.
1. splay(N_v)
N_vの属しているCompress Tree上でN_vを根にします.
2.
これがちょっとむずかしいです.
親がいない場合は目的達成なので終了です.
Compress Nodeの場合, こんな状態です.
Rake Nodeの場合, 例えばこんな状態です.
splay(r)をすると, こうなります.
3. splay(n)
4で行う操作を簡単にするために行います. なんで簡単になるかはsoft_exposeで解説したいと思います.
4. splice(N_v)
入れ替えます.
5.
これは何かというと, N_vがEdge Nodeのとき, spliceするとvのhandleの位置が変わります. これに対応するためです.
6.
これで, N_vをtoptreeの根にすることができます.
Soft Expose
soft_exposeは任意の頂点v, u間のパスのCompress Node vuを作る操作です！(やっとここまできた)
こんな形にtoptreeを変形します.
(8/5 なんか頭悪い画像になっていたので修正しました)
手順を先に言ってしまいます
soft_expose(v, u)
1. expose(v)
2. N_vとN_wが
- 同じ: toptreeの根はvuかuvなので, uvであれば反転する. soft_exposeおわり
- 違う: 続く
3. N_vをguardする(????)
4. expose(w)
5. N_vのguardを外す
6. N_vから見てN_uが右側なら, 反転させる.
7. おわり
toptreeの根をguardするとは, splay操作があってもtoptreeの根を変えさせないようにすることです.
これは, N_vをtoptreeの根にした後, N_vの左側にN_wを持ってくる必要があり, N_vが根であり続ける必要があるからです.
guardされているときのspliceの操作が少し違います.
nの親がguardされていて, 親から見て左側にある場合, spliceはnの左の子と交換しないといけません.
しかし, 親から見て右側にある場合, spliceはnの右の子と交換しないといけません.
これはtoptreeの, 葉がパスの辺の順番になっているルールに違反するからです.(toptree壊れる)
Path Query
soft_exposeができるようになると, パスに関するクエリを処理することができます.
パスの長さとか, パス中の辺の長さの最大値とかです.
各ノードに情報をもたせて, セグ木みたいに左の子の情報と右の子の情報を演算するみたいな感じです. これをすると, soft_expose(v, w)をしてvwを見た時にパスv-wについての演算結果が求まるはずです. やったね.
ひとまず終了...
link, cut, select, 各種クエリとかは後日にします... 疲れた...
niuez.hatenablog.com link cutかいた
僕が書くのサボった厳密な話をnoshi91さんが書いています こちらも読んでください
top-tree実装体験木
読者層が限定されすぎていませんか？
Link Cut Treeを書いたことがない人はこちら！
部分木クエリについてはこちら！
最遠点クエリについてはこちら！
え〜っと, 多分toptreeが書けました. 多分! verifyしたいんですが, Rustのバージョンで困っています...
実装
論文を読むとかける！(素振り) 間違ってたらごめんなさい
LinkCutTreeがならしO(logN)らしいし, これもO(logN)だと信じている(解析できねえ...)
読んだもの
- Maintaining Information in Fully-Dynamic Trees with Top Trees
- Design and Analysis of Data Structures for Dynamic Trees
- noshi91さんのツイート
toptree?
LinkCutTreeでは, Heavy-edgeのパスを頂点の順番でSplayTreeで管理していました. また, 上のLinkCutTreeでの部分木クエリでは, Light-edgeでつながったものを, multisetで管理しています.
これを, Heavy-edgeのパスを辺の順番で, 辺を葉とするSplayTreeで管理し, Light-edgeでつながったものをSplayTreeでまとめたデータ構造が, TopTreeの雑な説明です. このSplayTreeにクエリを処理させれば, 部分木クエリなどができます. link, cutもできます.
ここでもっと雑な説明をしてしまうとTopTreeが分からなくなってしまいそうなので, 他の記事にゆっくりまとめたいと思います... まとめる時間をください...(学業がなくなれば)
実際にできたクエリ
TopTreeにはClusterを載せてクエリを処理させます. ClusterのTraitは以下のようになっています.
pub trait Cluster: Clone { fn identity() -> Self; fn compress(left: Self, right: Self) -> Self; fn rake(left: Self, right: Self) -> Self; fn reverse(&mut self); }
v-uパスの長さ
impl Cluster for usize { fn identity() -> Self { 0 } fn compress(left: Self, right: Self) -> Self { left + right } fn rake(a: Self, _: Self) -> Self { a } fn reverse(&mut self) {} }
というふうにClusterを定義すると
pub fn path_length_test() { println!("path_length"); let v: Vec<_> = (0..13).map(|i| Vertex::new(i)).collect(); let edges = [ (0usize, 1usize, 1usize), (1, 2, 10), (1, 3, 3), (1, 4, 4), (0, 5, 3), (5, 9, 4), (9, 10, 7), (10, 11, 9), (10, 12, 1), (0, 6, 3), (6, 7, 3), (7, 8, 7), ]; let mut es = Vec::new(); for (a, b, w) in edges.iter() { es.push(link(v[*a], v[*b], *w)); } assert!(path_query(v[1], v[0]) == 1); assert!(path_query(v[0], v[4]) == 5); assert!(path_query(v[1], v[9]) == 8); assert!(path_query(v[3], v[11]) == 27); assert!(path_query(v[6], v[12]) == 18); assert!(path_query(v[12], v[6]) == 18); assert!(path_query(v[2], v[4]) == 14); assert!(path_query(v[5], v[6]) == 6); }
木の直径クエリ
#[derive(Clone, Debug)] struct Diameter { diam: usize, max_dist_left: usize, max_dist_right: usize, length: usize } impl Diameter { fn new(l: usize) -> Self { Diameter { diam: l, max_dist_left: l, max_dist_right: l, length: l, } } } impl Cluster for Diameter { fn identity() -> Self { Diameter { diam: 0, max_dist_left: 0, max_dist_right: 0, length: 0, } } fn compress(a: Self, b: Self) -> Self { Diameter { diam: *[ a.diam, b.diam, a.max_dist_right + b.max_dist_left].into_iter().max().unwrap(), max_dist_left: std::cmp::max(a.max_dist_left, a.length + b.max_dist_left), max_dist_right: std::cmp::max(b.max_dist_right, b.length + a.max_dist_right), length: a.length + b.length } } fn rake(a: Self, b: Self) -> Self { Diameter { diam: *[ a.diam, b.diam, a.max_dist_right + b.max_dist_right ].into_iter().max().unwrap(), max_dist_left: std::cmp::max(a.max_dist_left, a.length + b.max_dist_right), max_dist_right: std::cmp::max(a.max_dist_right, b.max_dist_right), length: a.length, } } fn reverse(&mut self) { std::mem::swap(&mut self.max_dist_left, &mut self.max_dist_right); } }
というふうにClusterを定義すると
pub fn diameter_cut_test() { println!("diameter cut"); let v: Vec<_> = (0..13).map(|i| Vertex::new(i)).collect(); let edges = [ (0usize, 1usize, 1usize), (1, 2, 10), (1, 3, 3), (1, 4, 4), (0, 5, 3), (5, 9, 4), (9, 10, 7), (10, 11, 9), (10, 12, 1), (0, 6, 3), (6, 7, 3), (7, 8, 7), ]; let mut es = Vec::new(); for (a, b, w) in edges.iter() { es.push(link(v[*a], v[*b], Diameter::new(*w))); } cut(v[0], v[5]); println!("0 diameter = {}", expose(v[0]).fold().diam); // -> 24 println!("5 diameter = {}", expose(v[5]).fold().diam); // -> 20 }
AOJの直径のやつ適当にいくつか通しました(提出できない)
最遠点クエリ
#[derive(Clone, Debug)] struct Farthest { ans: usize, max_dist_left: usize, max_dist_right: usize, length: usize } impl Farthest { fn new(l: usize) -> Self { Farthest { ans: l, max_dist_left: l, max_dist_right: l, length: l, } } } impl Cluster for Farthest { fn identity() -> Self { Farthest { ans: 0, max_dist_left: 0, max_dist_right: 0, length: 0, } } fn compress(a: Self, b: Self) -> Self { Farthest { ans: std::cmp::max(a.max_dist_right, b.max_dist_left), max_dist_left: std::cmp::max(a.max_dist_left, a.length + b.max_dist_left), max_dist_right: std::cmp::max(b.max_dist_right, b.length + a.max_dist_right), length: a.length + b.length } } fn rake(a: Self, b: Self) -> Self { Farthest { ans: 0, max_dist_left: std::cmp::max(a.max_dist_left, a.length + b.max_dist_right), max_dist_right: std::cmp::max(a.max_dist_right, b.max_dist_right), length: a.length, } } fn reverse(&mut self) { std::mem::swap(&mut self.max_dist_left, &mut self.max_dist_right); } }
というふうにClusterを定義すると
が解けます. サンプルは通りました(AtCoderさんverifyさせてくださいおねがいします)
pub fn farthest_test() { println!("farthest"); let mut buf = String::new(); std::io::stdin().read_to_string(&mut buf).unwrap(); let mut iter = buf.split_whitespace(); let q: usize = iter.next().unwrap().parse().unwrap(); let mut v: Vec<_> = (0..1).map(|_| Vertex::new(())).collect(); let edges :Vec<(usize, usize, usize)>= (0..q).map(|_| { ( iter.next().unwrap().parse().unwrap(), iter.next().unwrap().parse().unwrap(), iter.next().unwrap().parse().unwrap(), ) }).collect(); let mut es = Vec::new(); for (t, a, c) in edges.iter() { if *t == 1 { let new_v = Vertex::new(()); v.push(new_v); link(v[*a], new_v, Farthest::new(*c)); es.push((*a, v.len() - 1)); } else if *t == 2 { let p = es[*a - 1].0; let q = es[*a - 1].1; cut(v[p], v[q]); link(v[p], v[q], Farthest::new(*c)); } else if *t == 3 { println!("farthest from {} = {}", *a, expose(v[*a]).fold().ans); } } }
Nearest Marked Vertex Queryはまだやってない.
何が載るのかはさっぱりわかりません, だれか解明して.
ダイクストラとポテンシャルのはなし
はじめまして, niuezといいます. 競プロを少ししています.
最近勉強したことのメモ書きをしておきます.
ダイクストラ法
ダイクストラ法(Dijkstra)は負の長さの無いグラフで始点からの最短距離を求めるアルゴリズムです.
具体的には
距離が未確定の頂点の中で一番小さいものを選び, 距離を確定させる.
選んだ頂点から距離が未確定の頂点に伸びる辺で, 未確定な距離をより短いものに更新する.
を繰り返します. これを実装すると $O(N)$ですが, よく知られるダイクストラの計算量は $O((E+ V) \log E)$ です(heapとかを使う).
#include <set> #include <queue> #include <vector> struct edge { int u,v; int dist; }; std::vector<int> dijkstra(const std::vector<std::vector<edge>>& g, int s) { std::vector<int> dist(g.size(), 1e9); using node = std::pair<int,int>; std::priority_queue<node,std::vector<node>, std::greater<node>> Q; dist[s] = 0; Q.push(node(dist[s], s)); while(!Q.empty()) { int v = Q.top().second; int d = Q.top().first; Q.pop(); if(dist[v] < d) continue; for(const auto& e: g[v]) { if(dist[e.u] + e.dist < dist[e.v]) { dist[e.v] = dist[e.u] + e.dist; Q.push(node(dist[e.v], e.v)); } } } return std::move(dist); }
ベルマンフォード法
ベルマンフォード法(Bellman-Ford)は任意の長さのグラフで始点からの最短距離を求めるアルゴリズムです. 負の長さの閉路があるときはもちろん求められませんが, この記事では考えないことにします.
$O(VE)$ で直感的にもわかりやすいアルゴリズムですね.
#include <vector> struct edge { int u,v; int dist; }; std::vector<int> dijkstra(const std::vector<std::vector<edge>>& g, int s) { std::vector<int> dist(g.size(), 1e9); dist[s] = 0; for(int c = 0;c < g.size();c++) { for(int v = 0;v < g.size();g++) { if(dist[v] == (int)1e9) continue; for(const auto& e: g[v]) { if(dist[e.u] + e.dist < dist[e.v]) { dist[e.v] = dist[e.u] + e.dist; } } } } return std::move(dist); }
負の重みがあるときはベルマンフォード法しか無い?
ダイクストラ法のほうが定数倍が早かったりするので, できるだけベルマンフォード法よりはダイクストラ法を使いたいですよね?
一回だけ最短経路を求めるときはベルマンフォード法を使うしかありません.
複数回最短経路を求めるときはどうでしょうか?
実はこの場合, ベルマンフォード法を最初に一回だけしておくことで, 後の複数回はダイクストラ法を使うことが出来ます.
ダイクストラ法を使うとすると, 辺の長さをうまいことして正の長さにする必要があります.
最短経路とは??
始点を頂点 $s$ とした最短経路を数式に落とし込むと, こういう定義になります.
$d_s = 0$ とする.
すべての辺 $(i,j)$ において $d_i + dist(i,j) \ge d_j$ が成り立つときの, $d$ のそれぞれの取れる最大値.
これを頭に入れておくと次がわかります.
ポテンシャル
ここで天才をします. 先人は天才です.
ある$p_i$という値を用意して, 距離を $dist'(i,j) = dist(i,j) + p_i - p_j$ としたグラフを考えます.
長さ $dist'$ のグラフで, 頂点 $s$ を始点とした最短距離を計算して, ${d'}_i$ を求めたとしましょう.
${d'}_i + dist'(i, j) \ge {d'}_j$
${d'}_i + dist(i, j) + p_i - p_j \ge {d'}_j$
${d'}_i + p_i + dist(i, j) \ge {d'}_j + p_j$
よく見ると
$d_i = {d'}_i + p_i$
とすれば, $d_i$は最短距離の定義を満たしているように見えますね.しかし
$d_s = {d'}_s + p_s = p_s`$
なので
$d_s = 0$ を満たしていません.
なので, すべての頂点 $i$ について $ans_i = {d'}_i + p_i - p_s$ を計算すれば, $ans$ は最短経路を示しています.
このとき, $p_i$ のことをポテンシャルと呼びます.
では, $dist'$を正の長さにしたい気持ちになります.
$dist(i, j) + p_i - p_j \ge 0$
$p_i + dist(i, j) \ge p_j$
これは何かな. 最短距離の定義そのままですね(天才).
これを使うと色々なものが効率的に求めることが出来ます.
負の重みがあるグラフでの全点間最短距離問題
全点間最短距離問題とは, すべての頂点の間での最短距離を計算する問題のことです.
よく知られているのはワーシャルフロイド法の $O(V3)$ ですが, これを $V$ 回のダイクストラに置き換えることが出来て, $O(V(E + V) \log V)$ になります.
疑似コード
proc all_pair_shortest_path(G, dist)
let potential = bellman_ford(G, dist, 0) //引数は グラフ, 距離, 始点 です
for e = (i, j) in G
dist2(i, j) = dist(i, j) + potential[i] - potential[j]
for s in [0, |V| - 1]
result[s] = dijkstra(G, dist2, s)
for j in [0, |V| - 1]
result[s][j] += potential[j] - potential[s] // result[s][j]... s -> jの距離
return result
最小費用流
最小費用流はたぶん皆さんなら, 最短路反復法で実装していると思いますが...(RHS-algorithmなんて強多項式計算量知らない)
このとき負の辺があるときはダイクストラが使えないと思われがちですが, 同じように最初にベルマンフォードでポテンシャルを求めておけば, 高速で計算が可能です.
しかし, 逆辺が負の重みを持つので, ポテンシャルは, その時求めた最短距離を加算して行くことで, 更新をし続けなければなりません.
libalgoが参考になります.
スケーリングアルゴリズム
スケーリングを用いたダイクストラのアルゴリズムは重みが非負整数のときに使える高速化手法です.
簡潔に言うと, 「辺の重みを半分にしたものでダイクストラをして, その結果の二倍をポテンシャルに使ってダイクストラをする.」 を再帰的に行うことで, ダイクストラを高速化するテクを使うというアルゴリズムです.
下に例を示します.
このようなグラフがあったとします.
このグラフの重みを半分にした(小数点以下切り捨て) グラフでダイクストラをします.
最短距離は赤色に示した通りです.
この値を二倍した値を, 半分にする前のグラフのポテンシャルに使います.
辺の重みをポテンシャルによって置き換えると以下のようになります.
この置き換えた重みでダイクストラをします.
最短距離は青色で示しました.
それぞれの頂点で赤色と青色の値を足すと, ポテンシャルの性質により半分にする前のグラフの最短距離が求まります.
このグラフの重みは二進数にしたとき高々2桁なので1回半分にするだけで済みましたが, 一般に $ \log W$ 回再帰的に「重みを半分にして二倍してポテンシャルに使う」という動作をすれば求まります.
高速化
正直こんなことしなくてもこのままのアルゴリズムであれば, 大元のグラフをダイクストラすればいいだけの話です.
ですが, このポテンシャルで変更を加えた後のグラフに性質があります.
重みを半分にしたグラフでの, 頂点$s$から頂点$g$の最短経路($P$とします)に含まれる辺の数を $L$とします.
このとき, 半分にする前の重みをポテンシャルで変更を加えたグラフでの $s$から$g$の最短経路の重みは $L$以下です.
なぜなら, 二進数を考えると$P$上の辺の重みはポテンシャルで変更を加えると $0$ か $1$にしかならないからです.
つまりこのスケーリングアルゴリズムで行うダイクストラは, $V$個のQueueを用意してダイクストラをするものを使えば $O(m + n \log W)$ で計算できます.
〆
実はダイクストラの定数倍が速すぎてスケーリングはそんなに速くなりません
