CodeChef CHEFFILT

はじめに

この問題を本番中解いていた時に，$O(4^W)$から$O(2^W)$に落とせることに気づいて面白かったのでメモ．ぶっちゃけ気づいたのは，提出前に手元でテストしててDPの表を出力してみたからで，完全に偶然だったりする．

問題

参照

https://www.codechef.com/DEC15/problems/CHEFFILT

要約

$S~(0 \leq S < 2^W-1)$ と $n~(1 \leq n \leq 100,000)$ 個の数列$F_i(1 \leq i \leq n, 0 \leq F_i < 2^W)$ が与えられる．$F_i$ の部分集合のうち全要素の xor を取った結果が $S$ となるものの個数を$1,000,000,007$で割った余りを求めよ．（$F$ には重複した値が含まれるが，異なるものとして扱う）

想定解法

参照

https://discuss.codechef.com/questions/77624/cheffilt-editorial

要約

$F_i$ に含まれる相違なる値は $2^W$ 個以下である．そこで，$F$ 中で値が $x$ である要素の個数を $C(x)~(0 \leq x < 2^W)$とおく．

同じ値同士の xor は0(すなわち $x \mathrel{\mbox{xor}} x = 0$) であることから，ある値 $x$ が $k$ 個存在する場合，

偶数個を集合に含めた場合は，含めない場合と同じ結果
奇数個を集合に含めた場合は，$x$ 以外の全ビット列と $x$ の xor を取ったものが結果

となる．ここで，($k$ 個のものの中から偶数個選ぶ選び方）＝($k$ 個のものの中から奇数個選ぶ選び方）＝$2^{k-1}$ であることに注意すると，以下のDPを考えることができる．

$D(x, y) := $ $F$ で $x$ 以下の値のみからなる部分集合のうち，その全要素の xor が $y$ となるものの個数．

$$ D(x, y) = \begin{cases} D(x - 1, y) & C(x) = 0 \\ (D(x - 1, y) + D(x - 1, y \mathrel{\mbox{xor}} x)) \cdot 2^{C(x)-1} & C(x) > 0 \\ \end{cases} $$

この解法の計算量は $O(4^W)$ となる．

O(2^w) 解法

実は以下の性質が成り立つ．

任意の $y$ について，$D(x, y) \neq 0$ ならばその値は $x$ によって一意に定まる $$ D(x, y) = \begin{cases} 0 & D(x, y) = 0 \\ E(x) & D(x, y) \neq 0 \\ \end{cases} $$
$E(x)$ は $D(x - 1, x)$ の値から一意に計算可能 $$ E(x) = \begin{cases} E(x - 1) \cdot 2^{C(x)-1} & D(x - 1, x) = 0 \\ E(x - 1) \cdot 2^{C(x)} & D(x - 1, x) \neq 0 \end{cases} $$
$E(x)$の構成手順から，解は $D(2^W-1, S)=0$ ならば0, そうでなければ $2^{n-m}$．ただし，$m$は$D(x, x - 1) \neq 0$ なる $x$ の個数．

補題

$$ \begin{cases} D(x - 1, x) = 0 \wedge D(x - 1, y) \neq 0 \Rightarrow D(x - 1, y \mathrel{\mbox{xor}} x) = 0 \\ D(x - 1, x) \neq 0 \wedge D(x - 1, y) \neq 0 \Rightarrow D(x - 1, y \mathrel{\mbox{xor}} x) \neq 0 \end{cases} $$

補題の証明

$D(x - 1, x) = 0 \wedge D(x - 1, y) \neq 0 \wedge D(x - 1, y \mathrel{\mbox{xor}} x) \neq 0$ となる $y$ が存在すると仮定する．このとき，全要素の xor が $y$ になる集合の1つを $U$ ，$y \mathrel{\mbox{xor}} x$ になる集合を $V$ とおくと，$V$ から共通集合 $U \cap V$ を取り除き差集合 $U \setminus V$ を加えた集合の全要素のxorは $x$ となるため，矛盾する．$D(x - 1, x) \neq 0 \wedge D(x - 1, y) \neq 0 \Rightarrow D(x - 1, y \mathrel{\mbox{xor}} x) \neq 0$も同様に示せる．

証明

$x=0$ のとき， $D(0, 0) = 2^{C(0)}$ かつ $D(0, y) = 0~(y \neq 0)$ なので，成り立つ($E(0)=2^{C(0)}$)

$x-1$のとき，成り立つと仮定すると，$x$のとき，

$D(x - 1, x)=0$ならば，補題から，$D(x - 1, y) \neq 0$ なる $y$ について，$D(x - 1, y \mathrel{\mbox{xor}} x) = 0$であり，$y \neq z \Rightarrow y \mathrel{\mbox{xor}} x \neq z \mathrel{\mbox{xor}} x$であることから， $$ D(x, y \mathrel{\mbox{xor}} x) = D(x, y) = D(x - 1, y) \cdot 2^{C(x)-1} = E(x-1) \cdot 2^{C(x)-1} $$ であり，$y$ および $y \mathrel{\mbox{xor}} x$ 以外の値 $z$ については $D(x, z) = 0$ となる．
$D(x - 1, x)\neq0$ならば，補題から，$D(x - 1, y) \neq 0$ なる $y$ について，$D(x - 1, y \mathrel{\mbox{xor}} x) = E(x - 1)$であり， $$ D(x, y) = (D(x - 1, y) + D(x - 1, y \mathrel{\mbox{xor}} x)) \cdot 2^{C(x)-1} = E(x-1) \cdot 2^{C(x)} $$ であり，$y$ 以外の値 $z$ については $D(x, z) = 0$ となる．

以上から，$x$のときにも性質1,2が成り立つことがわかる．

性質3は，実際に$E$の漸化式を展開してみるとわかる．

実装

重要なのは$D(x, y) > 0$かどうかだけであり，個数は独立して計算できるので，$D(x, y) > 0$かどうかを $2^W$ の配列$A$で管理する．これで，$D(x - 1, x) > 0$かどうかを$O(1)$で判定できる．さらに，$D(x, y) > 0$ なる $y$ の集合をリスト$B$でも管理しておく．$D(x, x - 1) = 0$となった場合，リスト$B$の各要素$y$について$y \mathrel{\mbox{xor}} x$ を計算し，$A$および$B$を更新する．$D(x, y) > 0$ なる $y$ の集合は更新のたびに2倍になるため，更新は高々 $\log 2^W=W$ 回しか発生せず，$i$回目の更新の段階におけるリストの要素数は$2^i$のため，全体の計算量は$O(2^W)$

ソースコード

from collections import Counter
from sys import stdin
_data = iter(stdin.read().split('\n'))
input = lambda: next(_data)
MOD = 1000000007
L = 10
for _ in range(int(input())):
    s = int(input().replace('b', '0').replace('w', '1'), 2)
    n = int(input())
    ct = Counter(int(input().replace('-', '0').replace('+', '1'), 2)
                for _ in range(n))
    g = [False] * 2**L
    g[s] = True
    m = 0
    for f, x in ct.items():
        if not g[f ^ s]:
            for i in range(2**L):
                if g[i]: g[i ^ f] = True
            m += 1
    print(pow(2, n - m, MOD) if g[0] else 0)

追記

Editorialページのコメントに上がっていたのだけど，Gauss Jordan 使えば$O(\min(n, 2^W) \cdot W^2)$になるらしい．行空間の基底を求める問題に帰着できるっぽい？

JAG 夏合宿 2014 Day4 E - AI 解説＋裏話

はじめに

この記事は Competitive Programming （その2） Advent Calendar 2015 5日目の記事です．

1年半にわたって，自分が担当していた問題の解説を放置していたのでこの機会に書きます．

あと，単純に解法書くだけだとさみしい気がするので，作問体験記ということで，作問時の裏話も書きます．

問題

コンテスト: JAG 夏合宿 2014 Day4
問題文: AI | Aizu Online Judge
作問，解説: @mitaki28
テスター: @ustimaw

解法

この問題のプログラムにおいて，条件判定の結果はロボットの位置と向きにのみ依存します．
つまり，ロボットの位置と向きが全く同じ状態でプログラム中のある文が2回実行された場合，ロボットは以後，同じ動作を無限に繰り返すことがわかります．
- サンプル3でロボットがどう動くのか試してみるとわかると思います．
したがって，(ロボットのx座標)×(ロボットのy座標)×(ロボットの向き)×(プログラム中で実行中の文)について，1度発生した状態をすべて記憶していき，同じ状態に2回到達した時点でゴールに到達不可能なものとすることで解くことができます．
計算量は $O(hw|s|)$ となり，状態数は高々 $50 \times 50 \times 4 \times 1000=10000000$ で解けます．
- 状態数の上限がわかっているので，10000000ステップ回してゴールに到達しなかったら諦めるとかでも解けます．
- ただし，単純に動作文の数だけ数えていると後述する空のプログラムのパターンでTLEするため注意が必要です．

実装上の注意

if 文や while 文については条件によってプログラムを対応する括弧まで読み飛ばす必要があります．毎回，対応する括弧を見つけても良いですが，前処理で対応する括弧を計算しておくと，多少実装が楽になるかもしれません．
同じ状態に到達したかどうかの判定は，条件文の中で行う必要があります．動作文で判定した場合，while 文中のプログラムが空の場合に無限ループが発生してしまいます．以下のケースでTLEしている人が多かったです．
```
5 3
###
#g#
#.#
#s#
###
{T}
```
実装量的には1500〜3000byte程度になると思います．

ジャッジ解

@ustimaw, 1412 byte
@mitaki28, 2029 byte

#include <cassert>
#include <iostream>
#include <algorithm>

using namespace std;

const int MAX_W = 55;
const int MAX_H = 55;
const int MAX_L = 1100;
const int D = 4;
const int dx[] = {0, 1, 0, -1};
const int dy[] = {-1, 0, 1, 0};
const string ACT = "^v<>";
const string DIR = "NESW";

int H, W, K;
char a[MAX_H][MAX_W];
string s;
int y, x, d, p, gy, gx, ct;
bool used[MAX_H][MAX_W][D][MAX_L];
int jump[MAX_L];

void prog();
void stmt();
void if_stmt();
void while_stmt();
void action();
bool cond();

bool cond() {
  if (used[y][x][d][p]) throw false;
  used[y][x][d][p] = true;
  bool r = false;
  bool neg = false;
  if (s[p] == '~') {
    neg = true;
    p++;
  }
  if (s[p] == 'T') {
    r = true;
  } else if (s[p] == 'C') {
    r = (a[y + dy[d]][x + dx[d]] == '#');
  } else {
    r = (DIR.find(s[p]) == d);
  }
  p++;
  if (neg) r = !r;
  return r;
}

void action() {
  ct++;
  int nx, ny, nd;
  if (s[p] == '^') {
    ny = y + dy[d], nx = x + dx[d], nd = d;
  } else if (s[p] == 'v') {
    ny = y - dy[d], nx = x - dx[d], nd = d;
  } else if (s[p] == '<') {
    ny = y, nx = x, nd = (d + 3) % D;
  } else if (s[p] == '>') {
    ny = y, nx = x, nd = (d + 1) % D;
  }
  if (a[ny][nx] != '#') {
    y = ny;
    x = nx;
    d = nd;
    if (y == gy && x == gx) throw true;
  }
  p++;
}

void while_stmt() {
  int q = p;
  while (cond()) {
    prog();
    p = q;
  }
  p = jump[q];
}

void if_stmt() {
  int q = p;
  if (cond()) {
    prog();
  } else {
    p = jump[q];
  }
}

void prog() {
  for (;;) {
    if (ACT.find(s[p]) != string::npos) {
      action();
    } else if (s[p] == '[') {
      assert(s[p++] == '[');
      if_stmt();
      assert(s[p++] == ']');
    } else if (s[p] == '{') {
      assert(s[p++] == '{');
      while_stmt();
      assert(s[p++] == '}');
    } else {
      break;
    }
  }
}

void preprocess(char end) {
  int q = p;
  while (s[p] != end) {
    if (s[p] == '[') {
      assert(s[p++] == '[');
      preprocess(']');
      assert(s[p++] == ']');
    } else if (s[p] == '{') {
      assert(s[p++] == '{');
      preprocess('}');
      assert(s[p++] == '}');
    } else {
      p++;
    }
  }
  jump[q] = p;
}


int main() {
  cin >> H >> W;
  for (int i = 0; i < H; i++) {
    for (int j = 0; j < W; j++) {
      cin >> a[i][j];
      if (a[i][j] == 's') {
        y = i;
        x = j;
        a[i][j] = '.';
      }
      if (a[i][j] == 'g') {
        gy = i;
        gx = j;
        a[i][j] = '.';
      }
    }
  }
  cin >> s;
  s.push_back('$');
  
  p = 0;
  preprocess('$');
  
  p = 0;
  ct = 0;
  bool ok = false;
  try {
    prog();
  } catch(bool e) {
    ok = e;
  }
  cout << (ok ? ct : -1) << endl;
  return 0;
}

裏話

せっかくなので，作問時の裏話をつらつら書きます．

問題が作られた経緯

ICPCということで，強実装系の構文解析問題を入れようという話になった
過去の構文解析系の問題は式の解析が多かったので，if文やwhile文を使ってみた
プログラムで実際にロボットを盤面上で動かしてみると面白いんじゃないかと思いつく
- 発想自体は Misterious Gems | Aizu Online Judge とかと同じ感じがする．
単純な構文解析と見せかけてメモ化が必要というのがひねりがあって面白そうという話になり採用

テストケースについて

この問題は，むしろテストケースを作るのが大変だった
- ランダムケースを作るのが難しい
  - 適当に作ると，ほぼ100%ゴールまで辿りつけない
最終的にランダムケースは以下の手順で作成
1. ランダムな迷路を生成
2. ランダムなプログラムを生成
3. 1, 2を使って2つのケースを生成
  1. ランダムな位置にゴールをおいたもの(到達失敗ケース)
  2. ランダムなプログラムで移動する範囲の中でスタート地点から最も離れた位置にゴールを置いたもの(到達成功ケース)
とはいっても，今回の場合，ランダムケースはあまり当てにならないので，手動ケースも作成した
- 無限ループに関するもの: 4つ
- 空プログラムを含むもの: 1つ
- 50x50のほとんど空の盤面上を動き回るもの: 6つ
- if文やwhile文を限界までネストしたもの: 2つ

これでも，テストケースとしてはちょっと弱いなぁと思ってたら @ustimaw 君がマラソンマッチパワーを発揮してとんでもなく強力なケースを作ってくれていた

50 50
##################################################
#s..#..#.#.....#..............#...#.#.###...#....#
#...........#...#.............#...#.#.#...#..#..##
#....##.#......##..#.##...#.........#..##..#...###
#.#..#....#.....##....#.#..##....#...###..#.....##
#..###....#.##......#.#...#..###.###.............#
#...#..#.#.....#..#..#..#.#...#..#.#.......#.....#
#............#....#.##..#.....##...#...###....##.#
#..#.#....##.#.###..#...#.##..##......#..#.###.#.#
#.#..#.##............#............#.#.#...#....#.#
#..#........###...#......#.......#........##....##
#.....#..##.....#..#..#.#..####...#........#...#.#
#..#...#.##..#..#...#...#.......##...#.#....##..##
#..##.........#..#.........#...####......##.....##
#..#####...#..#..#.##..##.##........#...#........#
#..#........##.##..#.#.#..#..#.#.###.#.#..#..#.###
#..#.#..#..#.#...#.....#..#..##.....#.....#..#..##
#..#...##.#..#..#..##........#.#..######.....###.#
#....#.#.##.....#....#...##..#.##...#....###.#...#
##.......#.#...#..#....#.....#.#...#........#...##
#..#.##..#.........#..#.#...#......#....##..#....#
#.........#..#..#.#....#..####..##...#...#.#######
#.#.#..#..#..#..#....##....##.##..#...#.##.##..#.#
#...##.#..#.#......##...##...#....#..#..#.#...#..#
#...#..#..#..#####....##..##..#....#.#.....#.##..#
#.#....#..#........#..#.##...#...##..####...#.#..#
#.........#.#........#.........#..####..####.###.#
#..#..#..#....###..#....#..#...#.#..#...#...#....#
#.#...#..#..#.#.#..######.###.#...#.....#..##....#
##....#.#..#..#...##.#......#...#..####..#....#..#
#.###.......#......#.#.#....#...#...#...##.##.####
#......#.#.....#.#.....###.#..##...#..##..#...####
#....#.....#..#.#....##..##....#...#.#..##..##...#
#..#.#.##...#.....#......#....#....#..#..#...#...#
#..........##..#..#........##.#......#...###..#.##
#.###.#.......#.#.#...###...#...###...##.....#..##
#......#.###....##..#....###....#..#..##..#...#.##
#..#..#..#......#.#.#.....##....#...#..##........#
#...#.......###.........#.##.###...#..#.##....#..#
#.......#.#..#..##.#.##.....#............#.#.....#
#.#..#..#..........##........#..##.##....#.#...#.#
##..#..#.....#..#...##..#.....#..........##...#..#
#....#.#..#..#.#.#....#...##......#.####.......#.#
#....#.#..##....#.##....#.#..#.#.....#.....#.#...#
#...#.#..#.#.........##.#.#....##........#.....#.#
#.#..##......#...##..............#...###...##.##.#
#.......#.###....#.#.....#...##..##......#.#..####
#..#.#....#....#..#..##..#..#...#...##...#.....###
#...#..#...#....................................g#
##################################################
{T<><<><^><><<<<<<><<<<<<><<<<<<<<<><<<<<<><<<<<><<<<<>>><<<>^<<<<<<<<<<<<<><<<{Cv<<><<<<>><<<^v<<<><<><<<><<<^><><<<<<<><<<<><<<>>>>><<><<<<^<>>><<^<><<><<>><<<><<<<<>><<^><><><<<<<<<><<<<><<<<<<<<<<<>^>><>><<<><><<<><<<><<<<<vv>><<<<<<<<<<>>><<<<><<<<<>><<><<>>^<<<<<><<<><<<>><v><^<<<<^>>^<<<<<<<>^<<<><<<><<><<<v<<><<<<<<<><v<<<<><<<<><<<<<<><><<<<<<><<><<<>><<<<<<<<<>><<<<<<<><<<<<v<<><<<<<<<v<<<<><<<<<^<<<<><<<<<><<<<<<<<<<<<<<<<<>><<v<<>>><<<<<<<<^<<<<<<<<<<<<<><<<<<<^<<<<v<<<><<<<<<<<<<<<v<>><<<<<^<<<<<<<<<<><<<<><<<<v<><<><<><<<<><><<<<<<<<<<<<<><<<<>>>><<>><><<>><<<<>^<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<><<<<<><<<^<<<<<<<><<<<<<<<<<<<<<<<<<<<<><<<<<><<<<<<><<<<<<<<<<v<<<<<<v<<<^<<<<<<<^<<<<><<<<<<<<><<><<<<<<>><<<<<<<<<<<<><<><<<<><<<<><<<<<<<v<<<<<><<><<<<><<<<<<<><<^<<<<<<<><<<<<<<<<<><<<<<<><<<<<<<<<<<<<<<<<<<<>><<<<<<<<<><<<}<<<<<<<<<><<<>><<<<<<<<<<><>><<<<<<<<<<<>>><<<<<<<<<<<<<<>><><><<<<<<<<<<<<<<<<<<<<<<><<<<<<<<<<<<<<^>><<<><<><<<<<<>><<<<<^<<<><>><<><<>><<<<<<<<<><<><<>}

このケース，ゴールまでに1210310(動作文)かかる
- 理論上の最大値の1/10
これよりも長い動作文が必要なケースを作れたという方はご一報ください．

コンテスト中のジャッジルームにて

先ほど書いたように，「動作文のない無限ループが存在しうる」という罠があるにも関わらず，コンテスト開始時に問題文が『「プログラム」は複数の「文」からなり，...』となっており，あわてて修正した．
AtCoder の仕様で問題文を修正すると，"<" や ">" のエスケープが解除されてしまうというトラップがあり，Sample Input 中で使われていた"<" や ">" がぶっ壊れてしまい，更なる対応に追われた．
- 基本的にAtCoderで問題文を訂正するときは，AtCoder上から書き換えるのではなく，問題文を再生成してコピペするほうが安全
構文解析＋探索でそれほど内容的には難しくないと考えていたので，想像以上にWAやTLEが多くミスジャッジしてしまっているのではないかとヒヤヒヤしていた．
上記のトラブルや，WAやTLEの多さからジャッジルームでずっとヤバいヤバい言ってて周囲を苛つかせていたかもしれない．反省している．

Code Festival 2015 参加記

1日目

起床〜会場到着

朝5時起きで会場に向けて出発．

このとき，事前にアンケートで答えていたバスの予約を2日目のものと勘違いしていて，思いっきりすっぽかしてしまう．

おかげで，10時ごろに東京駅に到着したあたりで，携帯にスタッフの方から電話がかかってきていたり，Twitterでメンションが飛んできてたりすることに気づいて本番の開始時間を勘違いしてたんじゃないかと非常に焦っていた．

ご迷惑をおかけしてすみませんでした・・・

本戦

A, B はまあ解ける．

C は貪欲を書くも + と - を書き間違えてシャリが無限に増殖していたため 1WA

D はこの前のARCみたいに累積和使う賢い解法があるんだろうなぁと思いつつ，ぱっと思いつかなかったので LazySegmentTree でとっとと殴り倒す．

E は出力は高々6パターンしかないことがわかったので，場合分けした．

F は読んでみても方針が立たなかったのでスルーしてGを解く．

G は，以前，Markdownエディタを作る過程で，形の違う木の手抜きローリングハッシュが一致する木（違う形状で行きがけの順序が一致する木）にはどんなのがあるのだろうということを少しだけ考えていたので，問題はすんなり頭に入ってくる．その後，$O(n^3)$のDPをわりとスムーズに思いついてAC．途中 int 同士の掛け算がオーバーフローしてるのに気づかなくて15分ぐらい無駄にしたのがもったいなかった・・・

結局，残り1時間半を費やしても H が解けずにそのまま終了．最終的に思いついていた方針は，解説とは違っていたので，だめかなと思っていたけれど，本番後に書いてみたところ通ったので，かなり惜しいところまで解けていた．本番中に通したかった．

最終的な順位は27th．6完の中では2番目にペナルティが少なかったらしい．正直，凍結後に結構みんな F を投げていたし，31th以降に放り出されると思っていたので本当に運がよかった．

トークライブ

チーム Unagi のトークライブと @laycrs さんのトークライブを聞いた．

ICFPCは毎年参加したいと思いつつ，問題文の長さで心が折れてしまっているので，来年こそは参加したいと思った．

@laycrs さんのトークライブを聞いて，昔作問した時のことを思い出して懐かしんでいた．そういえば，あのときの解説まだ書いてないな・・・

エキシビション

途中， @Mi_sawa さんがエディタ上で解説し始めたり，@yosupot さんが画面使って散々遊びつつサラッとAを通していたりしておもしろかった．

そして，終了直前にりんごさんがまさかの全完．やばすぎる．

解散〜就寝

メインディッシュを食べそこねてしまって，お腹が空いていたのでホテルのロビーにいたメンバーでマクドナルドに行く．そして，就寝．

2日目

朝プロ

寝坊して結局10時に会場に着く．とりあえずHardのAはすぐ通したが，残り15分でBを解くのは無理ゲーな気がしたので，そのまま諦める．

トークセッション

ライトニングトークとりんごさんのトークセッションを聞く．ライトニングトークでは @uwi さんを初めてリアルで目撃した．りんごさんのトークセッションは去年にも増してカオスだった．

リレー

うちのチームは2人組を作って2問ずつ担当するという作戦を取った．僕は @imulan さんと組んでD, Iを担当．Dはすぐ解けたが，Iが解けない．$O(n^5)$のDPしか思いつかず，悩んでいるところを @qnighy さんに助けてもらう（ちょうど同じタイミングで @nikollson さんも @qnighy さんと同じ解法を話していた気がする）．結局，僕は聞いた解法を実装しただけだったけど，最終的な順位は7th．チームメンバーたちに感謝．

解散〜帰宅

すこし時間があったので，関西勢で適当にご飯を食べに行く．目黒周辺をうろうろしていて見つけたイタリアンの店でピザとサラダを軽く食べる．

帰りの新幹線の中で暇だったので本戦 H, I と朝プロ Hard の B を通す．

その後，無事家に帰る．

永続赤黒木を実装した時のメモ

はじめに

Copy&Paste の解説スライドを読んで，永続赤黒木を実装してみました．

merge-split ベースの赤黒木の実装は insert-erase ベースのものと比較するとかなりシンプルで，1重回転だけで実装できますし，行数的にもRBST+50〜60行で実装できます．

ただし，上のスライドだけだと，なんでこれで動くのかを結構考察する必要があったり，split に関しては擬似コードどおりに書いても正しく動かなかったりするので，備忘録も兼ねてメモしておきます．

実装前に理解しておくべきこと

P.11 の平衝条件はしっかり頭に叩き込んでおきましょう．
1. 各ノードには赤か黒の色がついている
2. 根は黒
3. 葉は黒
4. 赤いノードは隣接しない
5. 根から葉までのパス上の黒いノード数は一定
さらに，上記の条件から以下のことが導けます
1. 条件3,5 から，葉以外のノードは必ず2つの子を持つ
  - 持たないとすれば，条件5に反する
  - この性質を覚えておくと，実装上 null チェックの数が少なくて済みます
2. 変更によって根が赤になった場合，単純に根を黒に戻せば良い．
  - 条件4, 5は保たれるので問題ない
このアルゴリズムで値を持つのは葉だけ
- 一般に出回っているRBSTとは違うので注意しましょう．
各ノードはランクを持つ（P. 13）
- 葉から自身までの黒いノードの数．ただし自分自身は含まない．
- 条件5 から，a.l.rank + is_black(a.l.color) と a.r.rank + is_black(a.r.color) は一致する．
- したがって，計算時は左右のどちらの子を使って計算してもよい．

ノードの生成

葉と中間ノードの性質が全く異なるので，別々の生成関数を用意しておくと便利です．

function leaf(value)
    a := new Node.
    a.color := 黒.
    a.size := 1.
    a.rank := 0.
    a.value := value.
    a.l := null.
    a.r := null.
    return a.

function node(l, r, color)
    a := new Node.
    a.color := color.
    a.size := l.size + r.size.
    a.rank := l.rank + (1 if l.color = 黒 else 0).
    a.l := l.
    a.r := r.
    return a.

Merge

Merge 操作は P.17〜22 に書かれています．基本的には，P.18, 19 の擬似コードを実装すれば動きます．ただし，（上と同様）の部分は，左右反転して実装する必要があることに注意してください．全部書くと以下のようになります．

function mergeSub(a, b)
    if (a.rank < b.rank):
        c := mergeSub(a, b.l).
        if (b.color = 黒 and c.color = 赤 and c.l.color = 赤):
            if (b.r.color = 黒):
                return node(c.l, node(c.r, b.r, 赤), 黒).
            else:
                return node(node(c.l, c.r, 黒), node(b.r.l, b.r.r, 黒), 赤).
        else:
            return node(c, b.r, b.color).
    else if (a.rank > b.rank):
        c := mergeSub(a.r, b).
        if (a.color = 黒 and c.color = 赤 and c.r.color = 赤):
            if (a.l.color = 黒):
                return node(node(a.l, c.l, 赤), c.r, 黒).
            else:
                return node(node(a.l.l, a.l.r, 黒), node(c.l, c.r, 黒), 赤).
        else:
            return node(a.l, c, a.color).
    else:
        return node(a, b, 赤).

この節の以降では，実装時に，疑問に思ったことをメモしておきます．興味がなければ読み飛ばしてもらって構いません．

場合分けについて

P. 18 の実装では，a.rank < b.rank（a.rank > b.rank）が確定したあとは，P.20, 21の2パターンとそれ以外の場合しか場合分けをしていません．実は，これで，条件を満たすために必要な操作がすべて列挙されています．

null チェックを行っていない点

この実装では，b の子および c の子の色を参照していますが，null チェックは必要ありません．

まず，b の子についてですが，a.rank < b.rank が成り立っておりa.rank >= 0 であることから b.rank >= 1であることがわかります．したがって，bは葉ではなく，条件6から，葉以外のノードは必ず2つの子を持つので，null チェックは不要です．

また，c については，少なくとも1つ赤のノードを挿入することから c 自身が null になることはなく，c を作成する時点で，a も b も null になり得ない（少なくとも葉まで到達した段階で a.rank == b.rank となる）ことから，c の子についても null チェックは不要であることがわかります．

`c` が赤の場合は一旦そのままにしている点

b が赤のノードの場合に，そのまま c を左の子とした場合に，赤のノードが連続し，条件4が満たされないままになる可能性があるように思えるかもしれません．

しかし，実際には以下の2通りしかなく，問題ありません．

b が根のとき，条件7からbを黒にすれば平衝条件を満たす．
b が根ではないとき，b と b の親の間には条件4が成り立っていることから，bの親は必ず黒であり，この場合，b の親の段階で平衝条件を満たすように変更が加わるため，問題ない．

`c` が赤の場合， `c` の右の子は黒で確定する点

まず，重要なポイントとして，a と b の rank が突然入れ替わることはないため，最初に a.rank < b.rank が成り立てば，a.rank == b.rank が成り立つまでは，常に，a.rank < b.rank が再帰的に成り立ちます．つまり，a.rank < b.rank の場合は，再帰的に b の左側の子に対して，a をマージする操作が繰り返されることになります．

さて，c が赤である場合を列挙してみると，

c が a.rank == b.rank のケースで新たに挿入されたノードの場合， rank の定義から，a.rank == b.rank となるとき， b は常に黒のノードなので c の右の子は黒．
c が P.21 のケースで赤に変更されたノードの場合，右の子は黒．(実際には左の子も黒なのでこのケースが問題になることはない)
a < b.l.rank で，a が b.l.l にマージされ，なおかつb.l.color == RED のとき，c.r == b.l.r であり，条件4から b.l.r は必ず黒．

Split

split 操作は，P. 24 の擬似コードをそのまま実装しても正しく動きません． 正確には，以下の擬似コードを実装する必要があります．

function split(a, k)
    if (k = 0):
        return <null, asRoot(a)>.
    if (k = a.size):
        return <asRoot(a), null>.
    if (k < a.l.size):
        <L, R> := split(a.l, k).
        return <L, merge(R, asRoot(a.r))>.
    else if (k > a.l.size):
        <L, R> := split(a.r, k - a.l.size).
        return <merge(asRoot(a.l), L), R>.
    else:
        return <asRoot(a.l), asRoot(a.r)>. 
        
function asRoot(a)
    if (a is null):
        return null.
    if (a.color = 黒):
        return a.
    else:
        return node(a.l, a.r, 黒).

正しく動かない原因は，条件2を満たさないためです．a.l, a.r は a から切り離した段階で独立した木になるため，戻り値として戻す，または，mergeする前の段階で，条件7を使って，根を黒にしておく必要があります．

計算量について

merge: $O(|a.rank - b.rank|)$
split: $O(\log n)$

RBSTと比較すると，葉にしか値を持たなかったり，追加の情報を保持しないといけなかったりと，定数倍が重そうなイメージですが，merge の計算量が$O(|a.rank - b.rank|)$なので，実際には insert や erase の操作のように似たような rank の木をくっつける場合には RBST より高速に動作する・・・はず．また，merge が$O(|a.rank - b.rank|)$なおかげで，線形時間で初期構築できます．(P.26参照)

あと，split において，merge を何回も呼び出しているのに $O(\log n)$ になるのは，非直感的かもしれません．split の操作は，再帰の各段階において，

左の部分木を L に merge する
右の部分木を R に merge する

のいずれかの操作を繰り返し行っているものと考えることができます．このとき，部分木は rank の小さい順に merge されていくため，マージの計算量が$O(|a.rank - b.rank|)$であることから，$O(\log n)$になります．（$1,2,5, \ldots \log n$ のようなソート済みの列があった時に隣接する要素の差の和は $O(\log n)$ にしかなりません．）

永続 RBST を撃墜するケース

背景

RBSTはコピー可能は嘘という疑惑が興味深かったので僕もいろいろ試してみたところ，要素数 $n$ に対してそこそこの確率で木の高さが $n/4$ 程度になるクエリのパターンを見つけました．

テストケース

RBSTはコピー可能は嘘という疑惑のテストケースに対して，一旦 split してそのまま merge するという一見なんの意味もない操作を付け加えただけです．なんでこれで劇的にバランスが悪くなるのか僕にもわかりません．

int main() {
  RBST<Value> tr; 
  tr = tr.insert(0, 1); //trは1ノード

  srand(time(nullptr)); //このブロックは乱数を適当に消費するためにあります
  int rn = rand() % 10000;
  for (int i = 0; i < rn; i++) {
    int sz = tr.count();
    tr = tr.merge(tr); 
    tr = tr.split(sz).first; //同じものをくっつけて分離するだけの意味のない動作
  }
  cout << tr.count() << endl;
  assert(tr.count() == 1);

  for (int i = 0; i < 20000; i++) {
    int sz = tr.count();
    cout << sz << endl;
    RBST<Value> tr1, tr2;
    tie(tr, ignore) = tr.split(tr.count() / 2 + 1);
    tie(tr1, tr2) = tr.split(tr.count() / 2 + 1); //いったん分離して
    tr = tr1.merge(tr2);  // そのままくっつける
    tr = tr.merge(tr);
    printf("%d\n", tr.max_depth());
  }

  //ここでtrはバランスが崩れまくっている(ことがある)はず
  printf("%d\n", tr.max_depth());
    
  return 0;  
}

可視化

試しに，木のバランスが露骨に悪くなった時の形状を出力してみました．2本の長いパスのいずれかに頂点が集中してしまっている状態のようです．

rbst

CodeChef CHEFFILT

はじめに

問題

参照

要約

想定解法

参照

要約

O(2^w) 解法

補題

補題の証明

証明

実装

ソースコード

追記

JAG 夏合宿 2014 Day4 E - AI 解説＋裏話

はじめに

問題

解法

実装上の注意

ジャッジ解

裏話

問題が作られた経緯

テストケースについて

コンテスト中のジャッジルームにて

Code Festival 2015 参加記

1日目

起床〜会場到着

本戦

トークライブ

エキシビション

解散〜就寝

2日目

朝プロ

トークセッション

リレー

解散〜帰宅

永続赤黒木を実装した時のメモ

はじめに

実装前に理解しておくべきこと

ノードの生成

Merge

場合分けについて

null チェックを行っていない点

c が赤の場合は一旦そのままにしている点

c が赤の場合， c の右の子は黒で確定する点

Split

計算量について

永続 RBST を撃墜するケース

背景

テストケース

可視化

`c` が赤の場合は一旦そのままにしている点

`c` が赤の場合， `c` の右の子は黒で確定する点