PentominoをBitboardを使って解く

子どもが百均で買ってきたパズルをやってみたら、全然うまく出来なくて悔しかったのでプログラムで解を探すことにした。

Pentominoとは

調べるまで知らなかったが、このような 5つの正方形を繋げて作られる12種類のピースを使って敷き詰める、というパズルを Pentomino(ペントミノ) というらしい。

en.wikipedia.org

最も一般的なものはこれらを1種類ずつ使って長方形に敷き詰めるもので、5マス x 12種類なので面積60のものが考えられる。回転や反転によって重複するものを同一視すると、それらの解の数は以下の通り、とのこと。

  • 6x10 には 2,339通り
  • 5x12 には 1,010通り
  • 4x15 には 368通り
  • 3x20 には 2通り

そして冒頭の写真のものは2x2の正方形ピースが使われていて、合計面積64なので 8x8 の正方形が作れる。一般的にはこの2x2の正方形は中央に固定され、周囲のドーナツ状の部分に12種類を敷き詰める問題、となる。このときの解は 65 通りになる、とのこと。

探索アルゴリズムと実装

ということで、この問題を解くプログラムをRustで書いてみた。

計算量概算

まず、どれくらいの計算量がかかるのかを概算してみる。

ピースは12種類だが、反転と 90°, 180°, 270°の回転を考えると、それぞれに対し8通りの置き方がある。

12種類のピースの中から未使用のものを順番に選び
  選んだものに対して8通りの反転/回転のどれかを選び
    可能であれば(既に置かれている他のピースとぶつからなければ)置く

12種類の選ぶ順番だけで  12! (= 479,001,600) 通りあり、置き方を考慮するとその8倍、そして置ける場所を探すのにも計算が必要で愚直にやると盤面の全探索になったりして、とにかくかなりの時間がかかる、というのは容易に想像できる。

効率的な探索の方針

無駄な計算を省くために、まずは「端っこから順番に埋めていく」という方針を考える。 例えば最上段の左端から右端へ、次はその下の段の左端から右端へ、と。次に埋めるべきマスが決まっているので、現在選択しているピースではどうやってもそのマスを埋めることができない、となればその先の探索は必要なくなり枝刈りされる。

backtracking

こういう探索でよく使われるのがおそらくbacktrackingだろう。ざっくり以下のような感じ。

fn backtrack(board: &mut <盤面の状態>) {
    if すべてのマスが埋まっている(すなわち12種類のピースをすべて使いきっている) {
        答えを出力
        return;
    }
    for ピース in12種類のピース.filter(|ピース| 未使用) {
        for 反転/回転 in8通り {
            if 現在の盤面で最も上段・左端の空きマスを埋めるように置ける {
                boardに選択したピースを置いた状態にする
                backtrack(board);
                boardを1つ前の状態に戻す
            }
        }
    }
}

可能なものを順番に置いていき、置けなくなったら戻って次のものを置いていく、という処理を、盤面の状態を更新/復元しながら再帰的に探索していく。

あとはこれをどのようなデータ構造で実現するか、というところ。 長方形の盤面なので Vec<Vec<Option<u8>>> のように2次元配列で表していっても良いが、埋めるべき左上の空きマスを探すときや ピースを置けるか否かの判定などに毎回多くのループ処理が必要になってしまう。

Bitboardによる検索と判定

ここで登場するのがBitboard。 将棋ライブラリを自作 したときに学んだものが役に立つ。

u64 を使って64マスの盤面を表現する。0が空きマス、1が置けないマス、を表す。 8x8の中央に2x2の正方形が固定されていると考えると、初期盤面は以下のようになる。

00000000
00000000
00000000
00011000
00011000
00000000
00000000
00000000

この盤面のそれぞれのマスがu64の何bit目に対応するか、を割り当てていく。上段の左端から順に割り振っていくと、以下のようなレイアウトになる。

00 01 02 03 04 05 06 07
08 09 10 11 12 13 14 15
16 17 18 19 20 21 22 23
24 25 26 27 28 29 30 31
32 33 34 35 36 37 38 39
40 41 42 43 44 45 46 47
48 49 50 51 52 53 54 55
56 57 58 59 60 61 62 63

なので前述の初期盤面は 103_481_868_288 (0x0000_0018_1800_0000) という値で表現される。

こうして表される盤面Bitboardに対し、例えばT型のピースを一番左上に置くとき、

11100000
01000000
01000000
00000000

...

のようなBitboardを作る。盤面とピースのそれぞれのBitboardとして表現されるu64同士をAND演算(&)すると、重なっている部分が1になるので、(board & piece) == 0 で盤面とピースが重なるマスが無いかどうか(すなわち空きマスだけにピースを置くことができるか否か)、を判定できる。 また、その後に board | piece とOR演算(|)を適用することで、盤面にピースを置いた状態を得ることができる。

また、盤面の左上から順番に埋めていくという方針を考えると、前述のレイアウトにしておくことで、 「最も上段・左端の空きマス」の位置を u64::trailing_ones() を使って一発で取得できる。 これは x86_64ならbsf命令、aarch64ならclz命令を使って計算されるので非常に高速に動作する。

あとはすべてのピース、その反転/回転、それらを盤面内の可能な範囲で移動したものをすべてBitboardで表現できれば、探索が進められそうだ。

候補の事前計算

まずピースについて、もう少し深く考えてみる。 前述したように12種類のピースにそれぞれ反転と4種類の回転があるので、全部で 12 * 2 * 4 = 96 通りありそうだが、実際には線対称・点対称な形もあるのでそこまで多くはならない。例えばT型のピースは線対称な図形なので反転を考慮する必要がなく4通りだけ。I型の細長いピースは縦か横の2通りだけだし、X型のピースにいたっては回転しても同じ形になるので1通りしかない。これらを考慮すると、実際には合計で63種類の形しか存在しないことが分かる。

で、これらをそれぞれ各方向に移動したものを考えていくが、この移動によって配置可能な位置の数はピースの横幅・縦幅に依存する。 ちなみに前述のレイアウトにしておくことで、最も左上に配置したときのBitboardを基準として、1つ右に移動したものは piece << 1 で、1つ下に移動したものは piece << 8 で表現できる。 ただしこれは8x8の左上から並べるレイアウトのときに限る話であり、例えば6x10の横長盤面で同様に左上から並べるレイアウトでは、1つ下に移動する操作は piece << 10 になる。

全63種類の形を盤面上の全位置に置く、のを列挙するのは大変ではあるが、探索時にはこれらを毎回すべて確認する必要はない。

上述の探索方針で次に置けるピースを探すとき、埋めるべきマスの位置は分かっていて、盤面でそのマスよりも左上(下位ビット)のものは既にすべて埋まっている、という前提になっている。 ということは、次に置く試行をすべきピースは必ず「埋めるべきマスを含む」ことが分かっているし、さらに「そのマスより下位のビットを含まない」、つまりそのピースを表すBitboardでの最下位ビットが埋めるべきマスに一致する配置のものしか有り得ない、ということになる。

.O##.   ..O#.   ..O..
..#..   .##..   .###.
..#..   ..#..   ..#..

例えば上図のように.以外で表現される形のピースの場合、Oの位置が埋めるべきマスになっている必要があり、そうでないとピース内の他の箇所が既に埋まっているはずの盤面にぶつかってしまう。

この条件を考慮して「埋めるべきマスの位置」(00-63)のそれぞれに対して「そのマスを埋めることができるピースの候補」を引けるようにしておく。対象1箇所あたり高々63種類で、すべて事前に計算しておくことができる。

実装と実行結果

というわけで、ここまでを踏まえて以下のような実装。

use std::array;

const NUM_PIECES: usize = 12;

pub type Bitboard = u64;

pub struct SimpleSolver {
    table: [[Vec<Bitboard>; NUM_PIECES]; 64],
}

impl SimpleSolver {
    pub fn new(rows: usize, cols: usize) -> Self {
        assert!(rows * cols <= 64);
        let shapes = calculate_shapes();
        let mut table = array::from_fn(|_| array::from_fn(|_| Vec::new()));
        for (n, shape) in shapes.iter().enumerate() {
            for s in shape {
                if s.iter().any(|&(x, y)| x >= cols || y >= rows) {
                    continue;
                }
                let v = s.iter().map(|p| (1 << (p.0 + p.1 * cols))).sum::<u64>();
                let (w, h) = s
                    .iter()
                    .fold((0, 0), |(xmax, ymax), &(x, y)| (xmax.max(x), ymax.max(y)));
                for y in 0..rows - h {
                    for x in 0..cols - w {
                        let offset = x + y * cols;
                        table[s[0].0 + offset][n].push((v << offset).into());
                    }
                }
            }
        }
        Self { table }
    }
    fn backtrack(
        &self,
        current: Bitboard,
        solutions: &mut Vec<[Bitboard; NUM_PIECES]>,
        s: &mut [Bitboard; NUM_PIECES],
    ) {
        if !s.iter().any(|u| *u == 0) {
            return solutions.push(*s);
        }
        let target = current.trailing_ones() as usize;
        for i in 0..NUM_PIECES {
            if s[i] == 0 {
                for &b in self.table[target][i].iter() {
                    if (current & b).is_empty() {
                        s[i] = b;
                        self.backtrack(current | b, solutions, s);
                        s[i] = Bitboard::default();
                    }
                }
            }
        }
    }
    pub fn solve(&self, initial: Bitboard) -> Vec<[Bitboard; NUM_PIECES]> {
        let mut solutions = Vec::new();
        self.backtrack(
            initial,
            &mut solutions,
            &mut [Bitboard::default(); NUM_PIECES],
        );
        solutions
    }
}

通常backtrackingというとVecをstackとして使ってpush()pop()による更新/復元を使って探索していくことが多いと思うが、この問題の場合は必ず全部使い切るまで探索を続けるものなので固定長の配列を使っている。

手元の環境(MacBook Pro 14-inch, 2021 Apple M1 Pro 32 GB)でrelease buildを使うことで以下のような結果を得た。

  • 8x8-2x2: Found 520 solutions in 325.862ms
  • 6x10: Found 9356 solutions in 11.493s
  • 5x12: Found 4040 solutions in 24.804s
  • 4x15: Found 1472 solutions in 50.362s
  • 3x20: Found 8 solutions in 37.092s

解の数は、反転や回転による同一視をしていないので長方形では4倍、正方形では8倍の値になっている。 それぞれ割れば前述した通りの数になる。

Bitboardによる判定と候補の事前計算によってある程度は効率的に探索ができているはずで、8x8-2x2の盤面なら1秒以内に全解答を列挙できているが、他の長方形はかなりの時間がかかっている…

反転・回転での重複の除外

全解答を探索することはできたが、反転や回転で同一になる解を除外することはできていない。 これらを除外するためには、探索できた解とそれを反転・回転してできるものをすべてチェックし、既に発見済みの解と同一かどうかを判定する必要がある。

そのためには盤面に合わせてBitboardの反転・回転を実装する必要があるが、一旦これは後回しにして、とりあえず盤面のグリッド表現を愚直に作ってしまうことにする。

#[derive(Debug, Copy, Clone, PartialEq, Eq, Hash, FromPrimitive)]
pub enum Piece {
    O = 0,
    P,
    Q,
    R,
    S,
    T,
    U,
    V,
    W,
    X,
    Y,
    Z,
}

impl Solver {
    ...

    fn represent_solution(&self, solution: &[Bitboard; NUM_PIECES]) -> Vec<Vec<Option<Piece>>> {
        let mut ret = vec![vec![None; self.cols]; self.rows];
        for (i, b) in solution.iter().enumerate() {
            let p = Piece::from_usize(i);
            for (y, row) in ret.iter_mut().enumerate() {
                for (x, col) in row.iter_mut().enumerate() {
                    if b & (1 << (x + y * self.cols)) != 0 {
                        *col = p;
                    }
                }
            }
        }
        ret
    }
}

これによって2次元配列として表現される盤面を得られるので、単純に反転・回転を実装することができる。 8x8-2x2の場合はXY軸を反転させたものも考慮する必要があるので、それも用意しておく。

#[derive(Debug, Clone, PartialEq, Eq, Hash)]
struct Board(Vec<Vec<Option<Piece>>>);

impl Board {
    fn is_square(&self) -> bool {
        let (rows, cols) = (self.0.len(), self.0[0].len());
        rows == cols
    }
    fn flip_x(&self) -> Self {
        let mut ret = self.clone();
        for row in ret.0.iter_mut() {
            row.reverse();
        }
        ret
    }
    fn flip_y(&self) -> Self {
        let mut ret = self.clone();
        ret.0.reverse();
        ret
    }
    fn flip_xy(&self) -> Self {
        let mut ret = self.clone();
        for row in ret.0.iter_mut() {
            row.reverse();
        }
        ret.0.reverse();
        ret
    }
    fn transpose(&self) -> Self {
        let mut ret = self.clone();
        for (y, row) in ret.0.iter_mut().enumerate() {
            for (x, col) in row.iter_mut().enumerate() {
                *col = self.0[x][y];
            }
        }
        ret
    }
}

あとはこれを使ってHashSetなどに同一解を登録しておき、本当に未発見の解だけを記録するようにすれば良い。

use std::collections::HashSet;

struct UniqueSolutionStore {
    solutions: Vec<[Bitboard; NUM_PIECES]>,
    set: HashSet<Board>,
}

impl UniqueSolutionStore {
    fn add_solution(&mut self, board: Board) {
        if !self.set.contains(&board) {
            self.solutions.push(*pieces);
        }
        self.set.insert(board.flip_x());
        self.set.insert(board.flip_y());
        self.set.insert(board.flip_xy());
        if board.is_square() {
            self.set.insert(board.transpose());
        }
    }
}

impl Solver {
    fn backtrack<S: SolutionStore>(
        &self,
        current: Bitboard,
        pieces: &mut [Bitboard; NUM_PIECES],
        store: &mut S,
    ) {
        if !pieces.iter().any(|u| *u == 0) {
            return store.add_solution(pieces);
        }
        ...
    }
    pub fn solve(&self, initial: Bitboard) -> Vec<[Bitboard; NUM_PIECES]> {
        let store = UniqueSolutionStore::new(...);
        self.backtrack(
            initial,
            (1 << NUM_PIECES) - 1,
            &mut [Bitboard::default(); NUM_PIECES],
            &mut store,
        );
        store.get_solutions()
    }
}

これによって、

  • 8x8-2x2: Found 65 solutions in 362.164ms
  • 6x10: Found 2339 solutions in 11.707s
  • 5x12: Found 1010 solutions in 25.003s
  • 4x15: Found 368 solutions in 51.320s
  • 3x20: Found 2 solutions in 39.203s

と、全列挙していたときと探索時間はほぼ変わらず、正しく重複を除外した解を得ることができたようだ。

高速化

さて、とりあえずは正しく解を列挙できるようになったので、次は探索の高速化を考えてみる。 これには様々な手法やテクニックが古くから研究されていたようだ。

短辺から埋める

まず重要なのが、「最上段から順番に左端から右端へと埋めていく」と考えたときに、盤面が長方形の場合は横長の盤面より縦長の盤面の方が探索が早く終わる、ということ。

これは例えばU型のピースを最初に変な向きで置いてしまったときのことを考えると分かりやすい。

##.............
.#.............
##.............
...............

4x15の盤面で左上から埋めていくと、最上段の残りの13マスをまず埋める探索を繰り返していくが、頑張って埋めたとしても埋める対象が2段目に来た時点で左端の1マスの空間はどうしても埋められないことが分かりそこで探索が止まる。

これが15x4の縦長盤面だと

##..
.#..
##..
....

...

となり、もっと早く2段目の不可能なマスに到達できるので、そこまでの無駄な探索が少なくなる。

なので、この探索方針の場合は横幅は短ければ短いほど効率が良くなる、ということになる。 ので、6x10(6行10列)の横長盤面の場合は10x6(10行6列)の縦長盤面として、5x1212x5に、と縦と横を入れ替えて探索することで効率が良くなる。 解の表示が横長の方が良かったら表示時にだけまた入れ替えて戻してやれば良い。

当然8x8-2x2の正方形の盤面に対しては効果が無いが、それ以外の長方形たちは

  • 6x10: 約12s → 約1.02s
  • 5x12: 約25s → 約344ms
  • 4x15: 約51s → 約82ms
  • 3x20: 約39s → 約3.7ms

と何百倍、何千倍も高速化される。ものすごい効果だ。

X を最初に配置する

次に大きな効果を得られるのが、X型のピースを最初に配置して探索範囲を絞る、というもの。

4. Xピースによる解の重複排除 の記事が詳しいが、X型のピースは反転しても回転しても同じ形になるので、盤面左上の1/4の領域にだけ最初に配置することで探索範囲を絞ることができる。それ以外の領域にX型ピースを置くことで得られる解は必ずその左上領域に置いて得られた解を反転/回転することで得られるものになるからだ。

さらにX型ピースは最も左上に配置すると左上隅に隙間ができてしまうので有り得ない、ということも分かるので、実際に初期配置として有り得る箇所は  \lfloor\frac{(W-1)}{2}\rfloor\times\lfloor\frac{(H-1)}{2}\rfloor-1 通りになる。最も多くても5x12のときの 9通り程度だ。それらを配置した状態を初期状態として、X型以外の11種類のピースで今まで同様のbacktrackingで探索していけば良い。

こうして得られる解は必ずX型ピースが左上の領域に存在することになり、前述した重複除外の手順と逆に、これらが重複を含まない解となる。そしてそれらを反転/回転したものを含めるようにすることで全探索した場合の解を得ることができるようになる。実際にはX型ピースが奇数マス長の辺の中央にいた場合は反転での重複が有り得るし、8x8-2x2の正方形の場合はXY軸を反転させたものへの重複チェックなどが必要になる。

最初の1ピースの配置箇所候補が1/4になるので、単純に探索規模が約1/4に減ることになる、つまり約4倍の高速化が見込まれる。

結果として

  • 8x8-2x2: 約22ms
  • 6x10: 約103ms
  • 5x12: 約43ms
  • 4x15: 約9ms
  • 3x20: 約380μs

と、実際に約4〜10倍程度の高速化が確認された。

github.com

Bitboardの反転/回転

重複除外などの計算で盤面を反転/回転したものを得る必要があるが、簡単にするためにまずは盤面の2次元配列表現を愚直に作って処理してしまっていた。これを高速化するためにBitboardの反転/回転を実装する。 これによってHashSetなどで重複チェックする際にも、キーとして [Bitboard; NUM_PIECES] がそのまま使用でき、2次元配列を作る必要がなくなる。

とはいえ盤面によってビット表現のレイアウトは違うし、それぞれに対する反転や回転を計算するのは簡単ではない。 勿論64回のループを回して1ビットずつ反転後/回転後のビット位置を計算しながらコピーしていくこともできるが、それでは2次元配列作るのとあまり変わらず、旨味がない。

Delta Swap による手法

参考にしたのは以下の記事。 qiita.com

Delta Swap というXOR演算を利用した数回の操作で、maskで指定した範囲のビットを他の箇所のビットと入れ替えることができる。 これを利用することで、例えば隣り合う列同士、のように複数のビットの交換を同時に行うことができる。次は隣り合う2列同士、その次は隣り合う4列同士、とすれば8列の反転が3回のDelta Swapで実現できる。行・列のそれぞれの長さに対し  O(\ln{n}) 程度で盤面上での反転が機能する。 また、8x8-2x2におけるXY軸の反転も、上述の記事に書かれている通り斜めに並ぶmaskを指定することで同様に3回のDelta Swapで実現できる。

盤面の縦横サイズによってビットのレイアウトが変わるので、必要なmaskの値もシフトすべきdeltaの値も異なるものを用意する必要があるが、盤面が決定した時点で事前に計算できるので、それらをVecで用意しておくことで、実際にはDelta Swapの計算は高速に行うことができる。

Const Genericsなんかを使って盤面サイズごとに定数値を埋め込んでおくともっと効率的かな、とも思ったが、試してみた感じではそれほど差異は無さそうだった。実際にBitboardの反転/回転の計算が必要になるのは解が見つかったときだけなので、探索の計算量と比較すると無視できる程度の差にしかならないのかもしれない。

巨大テーブルによるループ効率化

とにかく探索のループが何度も実行されるので、この処理を高速化することが重要になる。

まず使用済みのピースをBitboard同様にビットで管理する。i番目のピースが使用済みか否かのチェックは used & (1 << i) == 0 で判定でき、探索の終了はこの used で使うビットがすべて1になっているか否かで判定できる。そうするとbacktrackingから戻ってきたときにpieces[i]をわざわざdefault値に戻したりする必要はなくなる。

impl Solver {
    fn backtrack<S: SolutionStore>(
        &self,
        current: Bitboard,
        used: usize,
        pieces: &mut [Bitboard; NUM_PIECES],
        store: &mut S,
    ) {
        if used == (1 << NUM_PIECES) - 1 {
            return store.add_solution(pieces);
        }
        let target = current.trailing_ones() as usize;
        for i in 0..NUM_PIECES {
            if used & (1 << i) == 0 {
                for &b in self.table[target][i].iter() {
                    if (current & b).is_empty() {
                        pieces[i] = b;
                        self.backtrack(current | b, used | (1 << i), pieces, store);
                    }
                }
            }
        }
    }
}

ここまでは特に効率化できていないが、ここでself.tableから候補を探索する処理を考えてみる。

このtableには target(埋めるべきマスの位置)それぞれに対してi番目の種類のピースで配置可能なBitboardのリストが格納されている。ピースがすべて未使用なら12個のVecをすべて調べるし、残り1種類なら1個だけのVecに対して判定をしていくことになる。いずれにしても2重ループでの処理となる。

ただ、table自体は変更されるものではないので実際にtargetusedの値の組み合わせに対して探索されるBitboardのリストは不変だ。 少し富豪的な発想だが「使用済みピースの状態数は高々4096(= 1<<12)程度なのだから、これも事前計算でそれぞれ1つのVecとしてすべて用意してしまえば良いのでは?」と考える。

    let mut table: [[Vec<Bitboard>; NUM_PIECES]; 64] = array::from_fn(|_| array::from_fn(|_| Vec::new()));
    for (i, shape) in shapes.iter().enumerate() {
        ...
        for target in 0..64 {
            let b: Bitboard = ...;
            table[target][i].push(b);
        }
    }

のようにVec<Bitboard>を格納する[64, 12]の2次元配列として準備していたものを、

    let mut table: [Vec<Vec<(usize, Bitboard)>>; 64] = array::from_fn(|_| vec![Vec::new(); 1 << NUM_PIECES]);
    for (i, shape) in shapes.iter().enumerate() {
        ...
        for target in 0..64 {
            let b: Bitboard = ...;
            for j in 0..(1 << NUM_PIECES) {
                if (j & (1 << i)) == 0 {
                    table[target][j].push((i, b));
                }
            }
        }
    }

Vec<(usize, Bitboard)>を格納する[64, 4096]の2次元配列に拡張する。こうすることで、usedに対応する全候補リストを一発で引けるようになる。

    fn backtrack(
        &self,
        current: Bitboard,
        used: usize,
        pieces: &mut [Bitboard; NUM_PIECES],
        store: &mut dyn SolutionStore,
    ) {
        ...
        let target = current.trailing_ones() as usize;
        for &(i, b) in &self.table[target][used] {
            if current & b == 0 {
                ...
            }
        }
    }

結局判定する内容は変わらないが、このようにループ処理を変えるだけでも意外に明確に速度が改善された。

ただし、初期化のコストは大幅に増大する…。tableの準備段階で1<<12回のループを回すことになるし、格納する値も冗長になりメモリ消費量が増える。手元の環境では6x10の盤面に対し巨大テーブルを使うようにすることで探索時間が約70msから約50msへと短縮されたが、そのぶん初期化にかかる時間が約0.1msから約70msへと激増していて、トータルの実行時間でいうとむしろ遅くなっている。このあたりは探索を速くしたいのかトータルでの実行時間を短くしたいのか、によって使い分けることになるだろう。

github.com github.com

もう一方のループ効率化

巨大テーブルを使わないとどうしても速くならないのか、というとそうでもなく、usedの判定をしながらループを回すのも多少は効率化できる。

毎回12回のループを回して使用済みの場合はskip、としていると無駄な計算が多くなる。

    for i in 0..NUM_PIECES {
        if used & (1 << i) != 0 {
            continue;
        }
        ...
    }

ので、ここでもビット演算によって未使用ビットに該当するものだけを列挙するように変更する。

    let mut u = !used & ((1 << NUM_PIECES) - 1);
    while u != 0 {
        let i = u.trailing_zeros() as usize;
        ...
        u &= u - 1;
    }

これも将棋Bitboardで学んだテクニック。u &= u - 1で最下位の1を消していくことができるので、これをすべて0になるまで繰り返すだけ。少数のビットしか未使用のものがない場合に無駄なループ処理を省くことができる。

これによっても速度は改善し、巨大なテーブルを用意しなくても90msかかっていたものが70ms程度になるくらいの効果があった。

孤立点検出による枝刈り

探索処理自体の高速化はこれくらいにして、最後に枝刈りについて考える。

探索途中のある状態において、この先はどう置いても解が得られないと分かっている場合は、その時点で探索を打ち切ることでその先での無駄な探索を省くことができる。この枝刈りをどれだけできるかでも大きな速度の差が出る。前述した縦長と横長の交換も枝刈りを多く発生させるための工夫だ。

まず一番分かりやすいのはU型ピースを壁際に配置するときで、

#.#.
###.
....
....

##..
.#..
##..
....

のように1箇所の穴が壁際に生じてしまう向きで配置すると、その穴のマスはどうやっても埋めることができないのでその先の探索はすべて無駄になる。

同様に四隅に空白ができるような配置も考えられる。左上から埋める限りはその左上に空間はできないように埋めていくが、それ以外の3箇所の角では起こりうる。

########
########

...

###.....
..##....

こうした配置は最初から考慮できるので、候補の事前計算時に除外しておくことはできる。

あとは探索途中の状態における孤立空間を検出して枝刈りできるかどうか。

######
###o..
.#....

...

のように#が埋まっているとき、次に埋めるべきマスがoになるが、ここで次に配置するピースによっては、その右側や左下に空間が生じ得る。

######
####..
.#.##.
....##

######
####..
.#.#..
..###.

のような感じ。少なくとも5の倍数の面積でない空間が生じていると、もうそれらを埋める手段は存在しないので、可能であれば即座に探索を打ち切ってしまいたい。

しかし、こういった空間の検出やその面積の計算まで毎回行うのは探索の高速化には逆に計算負荷が増大して逆効果になったりする。検出のコストと枝刈りの効果のバランスを考える必要がある。

Bitboardを使用して簡単に検出することができるのは、孤立点の空間だ。あるマスが空いていて、その上下左右がすべて壁もしくは他のピースによって埋まっている、という状態。これは

..#..
.###.
..#..

のようなX型のmaskを用意し、そのmaskと現在の盤面BitboardのAND演算の結果が

..#..
.#.#.
..#..

の値になるか否か、という判定で簡単に実現できる。盤面上のマスすべてに対してこのmaskと値を事前に計算しておけば、いつでもこういった「穴」のマスを検出することができる。

とはいえ実際にこのような「穴」が生じ得るのは、targetとなる埋めるべきマスにピースを置いた後の、その周辺だけ。なので毎回の探索ごとに穴の出現を判定するのは最小限にしたい。

targetを埋めるようにピースを置いたときに穴が生じやすいのはtargetのすぐ右隣のマス、もしくはすぐ下の段の左隣のマス、だ。先述の図でいうとoに対してXの位置。この2点だけに絞ってピースを置いた後に穴が生じるかどうかを判定して枝刈りを行うことにした。 本当はtargetの場所によって穴の発生が起こりやすい位置が違うかもしれないので一律に同じ相対位置にしないで調整した方が良いかもしれないが、そこまではやっていない。

######
###oX.
.#X...

...
    fn backtrack(
        &self,
        current: Bitboard,
        used: usize,
        pieces: &mut [Bitboard; NUM_PIECES],
        store: &mut dyn SolutionStore,
    ) {
        ...

        let target = current.trailing_ones() as usize;
        for &(i, b) in &self.table[target][used] {
            if current & b == 0 {
                let next = current | b;
                // `target`に対してチェックすべき`holes`は初期化時に計算しておく
                if self.holes[target].iter().any(|&(u, v)| next & u == v) {
                    continue;
                }
                pieces[i] = b;
                self.backtrack(next, used | (1 << i), pieces, store);
            }
        }
    }

最終的にこれらの工夫を施すことで、巨大テーブルを使用し初期化コストを無視することで探索自体は

  • 8x8-2x2: 約22ms → 約9ms
  • 6x10: 約103ms → 約44ms
  • 5x12: 約43ms → 約18ms
  • 4x15: 約9ms → 約3ms
  • 3x20: 約380μs → 約147μs

くらいまで短縮することができた。 6x10ではどうしてもまだ40ms以上かかってしまうが、まぁ10秒以上かかっていた頃からすれば200倍以上は速くなっているし、じゅうぶん速くはなったと思う。

github.com

あとはマルチスレッド化したりすればもっと速くなるのかなぁ…。

その他の実装・手法

ここまでまったく触れなかったが、Pentomino Solverの実装として"Knuth's Algorithm X"というものを使う手法があるらしい。

en.wikipedia.org

Pentominoを含むPolyominoによる敷き詰め問題をExtract cover problemとして考え、その解法としてDancing Linksというデータ構造を利用して効率的に計算していく、というものらしい。 以下の記事が詳しい。

matsu7874.hatenablog.com

実際に試していない(RustでLinked List的なものを扱うのが面倒そう、というのもあり…)が、これを使ってもそこまで劇的に速くなるわけでもないっぽい。 参考記事でも触れているが全パターンの探索には効率的かもしれないが枝刈りのような処理がないはず?なので、今回のようにX型を最初に置いたり孤立点チェックをしたりといった高速化の工夫の方が探索時間の短縮には効果的なように思える。

WebAssembly化してWebAppで使う

RustでSolverを書けたので、WASM化してWebAppから使えるようにしてみた。 最近はViteを使ってReact+TypeScriptの環境をさくっと作れると聞いたので、それで。 探索は100msかからないくらいだしメインスレッドで実行しても良いかな、とも思ったが念のためWebWorkerで実行するようにしている。

sugyan.com

手元の環境ではCLIで実行するときと遜色ないパフォーマンスがブラウザでも実現できているが、どうもスマホ(iPhone SE)で試してみたところ 巨大テーブルを使う方は最初の初期化に10秒以上の時間がかかる。メモリ確保がすごい重いんだろうか… 一度読み込んでしまえばその後はそこそこの速度になるようだけども。ともかくこれだとWebWorker経由じゃないと使いものにならない。

100msかからずに全解答を列挙できるようなブラウザ上で動く実装は他では見かけなかったので現時点では最速のWebAppになっているかもしれない。 もっと速いものを知っている方は教えてください。

参考

とても参考にしたページ。

まとめ

最大限にチューニングするとどこまで速くできるものなのかは気になるところだけど、とりあえず自分が満足できるくらいのパフォーマンスでの探索は実現できた。

AT Protocol(Bluesky)のためのRustライブラリを作った

前記事OCamlやってくぞ、と書いたけど結局Rustです。

Bluesky

Twitter代替の候補として噂される(?)分散型SNS Bluesky。 現状ではまだprivate betaということで招待コードが無いと使えないのだけど、先月運良くコードをいただくことができて使い始めてみている。

一通りのSNS的な動きはしているものの、まだまだ鋭意開発中ということで足りていない機能があったり 日々新機能が追加されたりと目紛しく動いている世界のようだ。

AT Protocolとエコシステムの現状 (〜2023/04)

AT Protocolについてはこちら。

atproto.com

日本語では以下の記事が詳しいでしょう。

ざっくりとした理解では、新しい分散型SNSのために「AT Protocol」が開発されていて、いま使っているBlueskyというSNSもこのAT Protocolに従って実装されているSmall-worldの一つである、という感じなのかな。

「AT Protocol」は様々なテクノロジーの総称という感じで、その中の汎用通信レイヤーとして「XRPC」というものがある。これを使ってサーバーと通信する、ということらしい。 BlueskyはWeb UIも提供されているが、DevToolsを覗いてみる限りでは確かにすべてのデータ取得や更新などこのXPRCのリクエストによって行われているようだ。

GitHubでは以下のリポジトリが公開されている。

主要なプロトコル実装としてTypeScriptで書かれた atproto があり、そのGo版ということで indigobluesky-social organizationによって管理されている。

その他の非公式の言語実装やClient, Applicationなどが以下のrepositoryにまとめられている。

言語としてTypeScriptがメインであり それを使ってすべての現存機能が動かせるので、Webフロントエンドで完結する、ということでWeb Clientが多く作られているようだ。

Rust版の実装

上記で紹介されているライブラリでRustのものは、(2023/04時点で)以下の adenosine のみ。

このrepository ownerはBlueskyの中の人ではあるが、公式としてGitHubのorganizationに入れられはしなかった、ようだ。

Lexiconとコード生成

前述したXRPCによるメッセージングは、すべて Lexicon というスキーマシステムによって型が定義されている。 JSON Schemaのように記述されるJSONによって、すべてのリクエストやレスポンス、その中のデータ型やフィールド名などが定義される。

ので、これがAPI通信の仕様のすべて、ということになる。

前述した atprotoindigoAPIライブラリは、このLexicon schemaを示すJSONファイルからコード生成によって作成されているようだ。

この仕組みで作っていない限り、どうしてもclient libraryは最新のAPIに追従できない、ということになる。AT Protocolはまだまだすごい頻度で変更が加えられていて、lexicon schemaも毎日、とはいかないまでも毎週くらいのペースでは何らかの変更がある。

上述の adenosine は、このコード生成によって作られたものではないため、実際既に最近追加されたAPIは使えないようだった。

その他にも reqwest の blocking client に依存している、とか 結局内部のxprc clientには serde_json::Value を渡していて何でも自由にできてしまう、など 微妙だなと思う点もあったので、自分でまったく別のRustライブラリを作ってみよう、ということになった。

(一応他にもそれっぽいlibraryやrepositoryを検索して先行事例を調査してみたりもしたが、ちゃんとそういったコード生成によって作られているものはなさそうだった。)

自作Rust版実装: ATrium

github.com

AT Protocolのための小さなライブラリが集まり青空を望む中庭、ということでATrium、という名前に。最初は atprs とか雑な名前で書き始めていたけど ChatGPT に「いいかんじのプロジェクト名ないですかね」と相談したところ、ATriumという名前を提案してくれたので採用した。最近 YAPC::Kyoto の会場でも聞いた単語だったし丁度いいな、と。

Lexicon schema

まずは必要になるのはLexicon schemaのJSONを読み込むためのライブラリ。 これ自身のJSON Schemaなどが存在するわけではなく、実際のvalidationは zodを使って書かれている ので、これを必死に翻訳する形で。

このJSON自体もかなり柔軟性があり、色んな入力が有り得るので大変…。いちおう type の tag を使って判別可能なようにはなっているので、それを使って enum に振り分けることでどうにか serde で正しく読み込めるようになった。 Option なfieldが多く、一度 deserialize してから元のJSONに同等に戻せるようにすると記述が大変だったので、serde_withskip_serializing_none を使った。

unionで新しい定義を作られるところが、単純に enum に変換すると一段深くなってしまうしtagがずれても困るし…ということで仕方なくいちいちコピーしてしまっているが、もっといい方法あるだろうか。…macroを使えば良いのか?

コード生成

ともかくLexicon schemaを正しく読み込むことができれば、あとはそこからRustコードに落とし込んでいける。 細かいところはまだ実装できていないが、基本的な object などはだいたい書けている、と思う。

Lexiconではdefsの名前は camelCase だが Rustコードでは 型名は PascalCase で field名やmodule名は snake_case にしたい、などの変換が必要だったので、 heck というライブラリを使用した。

正しくnamespaceをマッピングさせていれば ref はそのまま定義済みの型として参照できるが、 union はまぁ厄介で、それ用の enum を追加で定義して Lexicon解析のとき同様に $type tagで判別して振り分けるようにする必要があった。

とりあえず正しいコードとして生成されているかどうかは cargo build で確認できて 通ればまぁ動くだろう、という判断ができて便利だった。

API設計

型の定義だけならまだラクだが、Lexiconには query もしくは procedure のXRPCリクエストの定義が含まれる。 これは関数の実装を提供することになるが、そのためにはHTTP requestを実行するためのclientが必要になる、がAPIライブラリとしてはそういったものは含めたくない、と考えた。

実はRustでHTTP clientをマトモに使ったことが無かったが、 hyper だったり isahc だったり surf だったり、いろいろあるらしい。ライブラリでClientも含めることになると、どれを使うか判断する(もしくはすべてサポートするなどの)必要がでてきてしまう。

ので、それらは回避して HttpClient というTraitを定義して、それを実装するものを外部から注入する形にした。 各XRPCリクエストは「HttpClient をSupertraitとする、デフォルト実装を持つTrait」となり、外部で実装したClientがそれらの実装を宣言すれば各リクエストが使えるようになる、ということになる。

#[async_trait::async_trait]
pub trait HttpClient {
    async fn send(&self, req: http::Request<Vec<u8>>) -> Result<http::Response<Vec<u8>>, Box<dyn Error>>;
}

最初こういうのもライブラリとして存在しているのでは、と探して http-client というのがあったので使おうとしたが、同梱されている HyperClient を使ってみたら依存のtokioのバージョンが古すぎて動かず、またこのTraitで使われている http-types も良いけど http を使ったほうがより汎用的で良いですよ、とChatGPTが教えてくれたので自前で簡単に定義することにした。

ともかく、この形にすることで、ライブラリのユーザーは一手間はかかるけど自分の好きなHTTP clientを使って各種XRPCリクエストを送ることができるようになっている。はず。

とりあえずこうして作った APIライブラリをまず atrium-api という名前でpublishした。

docs.rs https://crates.io/crates/atrium-api

Lex

ちなみに Lexicon schema解析のための型定義だけのもので一つのライブラリとして切り出して作っていて、まだ publishはしていないけどやる予定。 上記の atrium-api の設計が気に入らない、という人はこれを使って自分でコード生成を作って別のAPIライブラリを作ってくれれば良いと思います。

CLI

作成したAPIの動作実験も兼ねて、 mattn/bsky のようなCLIも作ってみている。 これはこれで便利、なのかな? バイナリ配布だけしてみても良いかもしれない。

今後

まだまだ未完成で未実装のところも多いので埋めていく予定。 あとは本家のLexicon更新を検知して自動でコード生成してpushして publishまで…が自動化できても良さそうだけど そこまでは無理か? ちょっと調べてやってみたいところ。

あとは…折角Rustで動かせるようになるならTauriとか使ってGUIアプリまで作れると良いのかな… さすがにそんな余裕は無さそうだけども。

ともかくちゃんと使ってもらえるかもしれない(現時点で22 starsとかだけど)OSSっぽいものを初めて作った気がするので、まずは色々整備していきたいところ。

40歳から始める関数型言語、OCaml

というわけで、今年に入ってからのここ数ヶ月、OCamlを勉強し始めている。

動機

これまで「関数型言語をちゃんと触ったことがない」ということが若干コンプレックスになっていた。40歳になった今、唐突に始めてみる気になったので、やることにした。

Why OCaml

そもそも関数型にどんな言語があるのか、それぞれがどんな特徴であり、どういったところで使われているのか、という情報すらほとんど知らなかったので、「関数型言語やってみたいんですけど、どれがいいんすかね?」と某所で雑に相談したところ、OCamlを挙げてもらったので、素直にそれに従ってみることにした。

そんな適当で良いのか?とも思うけど、とにかくどれであっても触ってみないことには何も分からないし。実際、OCamlを選んで間違いは無かったと今は思う。

学習方法

Real World OCaml

とにかくOCamlの学習にはこれが良いらしい、という話をきいた。

dev.realworldocaml.org

ので、これを最初の数章は頑張って読んでみて、あとは実際にコードを書いてみつつ必要になったら参照して勉強していく、という感じで。 実際とても詳しく丁寧に書かれていて、とても良い本だと思う。

Github Copilot と ChatGPT

これらはもう現代のプログラマにとって欠かせないものになっていると思う。

最近のメジャーな言語と比較してOCamlに関する学習データがどの程度の量と質で、どの程度優れたコード出力能力であるのかは分からないが、少なくとも初心者から始める人間にとっては十分に助けになる。

特にChatGPTは、環境構築からライブラリの使い方などは勿論、「Rustでこう書いてたのはOCamlではどう書くの」とか「こういうのできないの」「こう書いてみたけど他にも書き方あるのかな」みたいな、独り言で呟くようなことでもChatで気軽に書いて質問できて体験が良かった。 めちゃくちゃ親しくて時間を奪うこと気にせずに何でも質問できる友人、が隣に居ればそういう人に頼れるかもしれないが、そういう友達が居ない人間にとってはこうした心理的なハードル低く幾らでも質問できる存在はとても有り難いと思う。

とはいえやはり正確性は保証されなくてChatGPTも平気で実在しないライブラリ関数をでっちあげて使用例を提示してきたりするし、油断せずにちゃんと公式ドキュメントなど探して確認しながら書く必要はある。 このへんは仕方ないと腹を括って、「俺が学習データを提供してやるぜ」くらいの気概でできるだけ良い(と自分が思う)コードを書いてとにかく公開していく。

オンラインジャッジ (競プロ)

簡単なコード片を書いて正しい入出力を得られたか確認できる、また別の人の書いた同じ目的のコードを気軽に参照できる、という点で競技プログラミングの過去問などは良い学習材料になると思う。とりあえず他人と競うことは考えず、オンラインジャッジシステムとして利用させてもらう。

個人的には LeetCode でやりたかったが、LeetCodeにはOCamlの選択肢が無い…。ので、 AtCoder を使うことにした。 AtCoder Problems の「Boot camp for Beginners」 というのがあることを @naoya_ito さんのTweet で知ったので、真似してこれらを解いてみることにした。

とはいえAtCoderOCamlのランタイムや使えるライブラリが少し古くて(2023/04時点)、手元で Core 最新を使って動くコード書けて意気揚々とSubmitしてみたら CE (Compilation Error) 出てガックリすることも。 ここは今後 言語のアップデート で更新されてより使いやすくなることを期待している。

簡単な問題を解くことで基礎的なデータ構造やアルゴリズムの扱いは練習できるが、難易度が上がると当然 数学的な思考力やより高度なアルゴリズムの知識が主眼になってくる。 そこは言語習得とはまた別の話になってくるのでそこまで頑張る必要はない(本当に競プロを頑張るならC++とかRustとか他の言語でやったほうが良いと思う)、ということでほどほどにして止めるつもり。

Advent of Code

より実践的に「プログラミングによる問題解決」をできるようにする練習として、 Advent of Code は非常に良い題材だと思う。 AtCoderのような競技プログラミング的な要素もあるが、それだけではなく

  • 入力データが必ずしも扱いやすいものばかりではない
    • → 様々な形式に対応したparse処理を書く練習になる
  • part1/part2 で同じ入力から異なる解を出す、など出力も様々な型・形式になり得る
    • → interfaceを考える練習になる
  • アルゴリズムだけで解決しないこともある
    • → 力技での膨大な探索や泥臭い実装などが必要になることも

といった点でまた別の実装力を鍛えることができると考えられる。

2022を完走した記事 にも書いたが、2022年の問題は特にバランス良く様々な題材のものがあって言語習得の練習に十分に適している、と感じた。

というわけで、OCamlAdvent of Code 2022 を解くチャレンジをしている。

github.com

まだ18日目までだが、どうにか最後まで頑張りたい…!

その次?

やはり「問題を解く」だけではなく、最終的には「何らかのアプリケーションを自作する」ことができるくらいにはなりたいので、それを考えていく。

自分の場合は「過去に他の言語で作ったことがあるもの」「自分が使いたいもの」を作るのが最もモチベーションが高くなるので、そういうもので考えると、、、今はコンピュータ将棋関連だろうか。 幸い(?) opam を検索してみた限りではビットボード実装や将棋関連のライブラリなどは存在していないようなので、自分で作ってみるのは面白そうだ。関数型言語での将棋プログラムなんて全然需要は無さそうだけどw とりあえずOCamlでの実装がC++やRustと比較してどの程度の速度になるのかは知りたいので perft できるくらいまではやってみたい。

あとは今の流行でいうと Bluesky で使われている AT Protocol の実装をしてみる、とかだろうか。これもわざわざOCamlでやろうとする人はあんまり居なそう…

所感

とりあえず数週間触ってみての感想。

|> は多用することが分かったので自作キーボード(Claw44)のキーマップにマクロとして追加した。

関数型という概念

結局、現状どれくらい理解できているだろうか…。まだまだ理解が浅いとは思う。

基本的に Base ライブラリを使用しているが、 List に対する様々な処理が再帰で表現できるというのを見て驚かされた。というか List というデータ構造がこんなにも使いやすくて色々できるんだなぁ、と少し親しみを持てるようになった。 とりあえず、 List などを使って |> でパイプを繋げてデータを変換していく、という操作は慣れてきた。

最初はどうしても「手続き型の言語だと普通こう書くので、それを関数型で表現すると…」と変換していく感じにはなりがちだったが、だんだん「この関数とこの関数を組み合わせてこういう関数を作って、そしてこの入力にこういう出力を返す関数を作る、」という感じで考えられるようになってきた気がする。

高階関数の部分適用とか使ってみるとめっちゃ便利だな〜と感心したり。

プログラミングHaskell という本を何年か前に買ったものの序盤ですぐに挫折してしまっていたのだけど、OCamlを通じて多少なりとも理解が深まってきたおかげか、今は再チャレンジして楽しく読めている。

OCamlの書き味

思っていた以上にツールチェインなどが整っていて、環境構築からエディタ設定などほとんど問題なくできて良かった。 VSCodeでは OCaml Platform 入れるだけだし、 ocamlformat を使って自動整形もできる。このあたりの「揃っていて欲しい」ものは大抵ちゃんと揃っている。

コンパイラのエラーメッセージは正直あまり親切な感じはなくて、型が合わないのはどこをどうすれば… みたいなのがなかなか詰まりやすかったりした。とはいえ「ちゃんとコンパイルが通ればだいたいちゃんと動く」という感覚がRustに近いものがある、と感じる。

あとは dune がRustでいう cargo と同じような感覚で使えるので、これも便利。

Rust, Python の経験

Rustを書いたことがあると、色んなところで「似ているな」と感じるところがあった。例えば option とかそれに対するパターンマッチとか。 Rustのiteratorをメソッドチェーンで繋いで処理するのと同じような感覚で |> で繋ぐ処理を書けるな、とか。 Rustがまったくの未経験だったらもっとOCamlに戸惑っていたかもしれない。

Pythonはそんなに似てないかな…? Haskellみたいに内包表記があるともっと違う感想になったかな。 ただOCamlの影響で自分が書くPythonコードがちょっと変わってきたような気がする。関数の使い方をより意識するようになったというか、なんというか。上手く言語化できないけど。

関数型言語をやるとプログラミングが上手くなる」という言説を聞いたことがあるが、こういうところで少しは効果がでるのかな…?気のせいかもしれない。

AIとの親和性

ところで Github Copilot にコードの補完をしてもらうときに、関数のシグネチャの情報はわりと重要だと思っていて。 「こういう名前の関数でこの型のこの名前の引数を受け取って、この型の値を返す」という情報が書いてあると当然補完も正確になりやすいだろう、と。

一方でOCaml型推論が強力なので、型はあるもののtype annotationを書く必要が殆どなく、またちょっとした処理は無名関数を書くことが多いし、そうでなくとも let rec loop acc x = ... とか let f x y = ... とか、簡潔な書き方になりがち。 そうすると、補完しようとするAIとしては情報が足りなすぎて難しいよなぁ、ということを思った。

明瞭な関数名や型注釈をいちいち書くのは面倒で人間としては省略できる方が有り難いが、それはそれでAIに助けてもらいながら書く場合には意図や情報を渡せるという点で良いのかもしれない。 今後新しく作られる言語の設計もそういう観点が考慮されるようになってくるんだろうか、とか。

まとめ

コードを書く仕事はAIに奪われていくかもしれないが、プログラミングは楽しいので好きなように書きたいだけ書いたら良い。

YAPC::Kyoto 2023 に参加した #yapcjapan

yapcjapan.org

おそらく YAPC::Tokyo 2019 以来?4年ぶりにオフラインのイベントに参加しました。地元開催ということで日帰りで行けてありがたい…! 最初だけちょろっと家族で参加。

子ゃーんの相手をしてくださった方々、ありがとうございました。

久しぶりすぎて誰と会ってもキョドってしまうくらいだったけど、とにかく色んな人たちとリアルで再会できて嬉しかった!!

そんなに長居はしてなかったのでトークは多くは訊けなかったけど、

id:ar_tama さんの『あの日ハッカーに憧れた自分が、「ハッカーの呪縛」から解き放たれるまで』がとても良かったです。 このトークを聴いた若手のエンジニアがこれを心に刻んで大きく成長し、また10年後のYAPCで良い発表してくれるといいな、と思いました。

(…と思ったら本当にベストトーク賞だったようで。おめでとうございます!)

とにかく楽しい1日でした。 運営の皆様、今回もありがとうございました!!

Advent of Code 2022 を完走した

毎年12月に開催されている Advent of Code に、2019年から参加している。

過去記事:

2022年のAdvent of Codeにも挑戦していて、年が明けてしまったが先日ようやく25日すべての問題に解答して 50 個のスターを集めることができた。

2022年の問題もどれもとても面白かった。

データ構造を扱う問題、ある程度のアルゴリズム知識が求められる問題、ちょっとした数学的な考え方が必要な問題、ひたすら実装が大変な問題、視覚化が面白い問題、などが揃っているのは毎年だが、特に2022年はバランス良くバラエティに富むものが出題されていて取り組み甲斐があったように思う。

工夫のしどころも多くあって、自分では上手く解いたつもりでも Reddit で他の人の解答を見ると目から鱗のものがたくさん発見できたり。

基本的にはヒント無しで自力で正解まで辿り着いたが、day19のpart2だけはどうにもならなくて他の人の解答を見たりしてようやく、という感じになってしまった。悔しい。

全部Rustで解いていて、あとでPythonでも解こうとしている。

github.com

とても楽しかったし勉強になったので、もうちょっと同じ問題に取り組んだ人達でわちゃわちゃと「ここが難しかった」「これはこんな手法があって」「こんな書き方もあって」など議論できると嬉しい…。 それはRedditで英語で頑張れば良いのだろうけど、ともかくもう少し日本語圏の開発者の間でも流行って仲間が増えてくれるといいな、と思っている。

という想いもあり、2020のときにもやったけど 日本語の解説本を書いてみている。

github.com

まだ半分もいってないけど… 少しでも多くのプログラマの人たちに届くといいな

2023パズル をRustで解いてみる

tkihiraさんの問題が面白そうだったので挑戦してみた。

既に解説記事が出ているので解答はこちらをどうぞ。

nmi.jp

結局自分は自力では解けなくて 他の人の解法や上記の解説記事を読んでようやくできた、のだけど… 自分なりに理解して改めてRustで実装してみた。

RPN(逆ポーランド記法)の backtracking

まずは型定義。

#[derive(Clone, Copy, PartialEq, Eq)]
enum Op {
    Add,
    Sub,
    Mul,
    Div,
}

#[derive(Clone)]
enum ExpressionElement {
    Operand(i32),
    Operator(Op),
}

この ExpressionElement を stack に詰んでいって、末尾が Operator だったらさらに Operand を 2つ pop() して、演算結果を push() していくことで計算が実現できる。

常に Operand の方が多い状態でないと Operator は挿入できないので、各々の出現回数をカウントしながら条件を満たしているときだけ操作するようにする。

様々な実装が試せるように Trait を定義。

trait Rpn {
    fn traverse(
        &mut self,
        expr: &mut Vec<ExpressionElement>,
        (i, j): (usize, usize),
        results: &mut Vec<Vec<ExpressionElement>>,
    ) {
        if i == j + 1 && self.get(i).is_none() {
            if self.evaluate(expr) {
                results.push(expr.clone());
            }
            return;
        }
        if let Some(&n) = self.get(i) {
            self.backtrack(expr, (i + 1, j), results, ExpressionElement::Operand(n));
        }
        if i > j + 1 {
            for op in &[Op::Add, Op::Sub, Op::Mul, Op::Div] {
                self.backtrack(expr, (i, j + 1), results, ExpressionElement::Operator(*op));
            }
        }
    }
    fn evaluate(&self, expr: &[ExpressionElement]) -> bool;
    fn get(&self, i: usize) -> Option<&i32>;
    fn backtrack(
        &mut self,
        expr: &mut Vec<ExpressionElement>,
        (i, j): (usize, usize),
        results: &mut Vec<Vec<ExpressionElement>>,
        e: ExpressionElement,
    ) {
        expr.push(e);
        self.traverse(expr, (i, j), results);
        expr.pop();
    }
}

self.get(i)i番目に使う数字を呼び出す。すべて使い切っていて Operatorも適切な回数使われていたら self.evaluate(expr) で計算結果を確認し、条件に合致していた場合のみ その exprresults に格納する。

exprpush / pop して再帰的に self.traverse() を呼ぶ backtracking の部分は別の実装で上書きできるようにあえて独立のmethodで定義。

何も難しいことを考えずにこれを満たす実装を作ると、以下のようになる。

struct Fraction(i32, i32);

struct DefaultSearcher {
    nums: Vec<i32>,
    target: i32,
}

impl DefaultSearcher {
    pub fn new(nums: Vec<i32>, target: i32) -> Self {
        Self { nums, target }
    }
}

impl Rpn for DefaultSearcher {
    fn evaluate(&self, expr: &[ExpressionElement]) -> bool {
        let mut stack = Vec::new();
        for e in expr {
            match e {
                ExpressionElement::Operand(n) => stack.push(Fraction(*n, 1)),
                ExpressionElement::Operator(op) => {
                    if let (Some(n0), Some(n1)) = (stack.pop(), stack.pop()) {
                        if let Some(n) = op.apply(&n1, &n0) {
                            stack.push(n);
                        } else {
                            return false;
                        }
                    }
                }
            }
        }
        stack
            .last()
            .map(|n| n.1 * self.target == n.0)
            .unwrap_or(false)
    }
    fn get(&self, i: usize) -> Option<&i32> {
        self.nums.get(i)
    }
}

impl Op {
    fn apply(&self, lhs: &Fraction, rhs: &Fraction) -> Option<Fraction> {
        match self {
            Self::Add => Some(Fraction(lhs.0 * rhs.1 + rhs.0 * lhs.1, lhs.1 * rhs.1)),
            Self::Sub => Some(Fraction(lhs.0 * rhs.1 - rhs.0 * lhs.1, lhs.1 * rhs.1)),
            Self::Mul => Some(Fraction(lhs.0 * rhs.0, lhs.1 * rhs.1)),
            Self::Div if rhs.0 != 0 => Some(Fraction(lhs.0 * rhs.1, lhs.1 * rhs.0)),
            _ => None,
        }
    }
}

解説記事では浮動小数点でも十分ということでやっていたけど、自分的には整数値だけで処理したいので Fraction(i32, i32) を定義して分数で計算するようにした。約分のことはとりあえず考えなくてよくて、最終的な値の計算結果の確認は 分子が分母の 2023 倍になっているか否か、で判定できる。 各Opの適用結果を Op::apply で計算。ゼロ除算の可能性があるので返り値は Option<Fraction> に。

self.evaluate() 内で stackを用意し、exprを順に見ながら Operand だったら Fraction に変換してpushし、 Operator だったら 2つを pop して演算結果をまた格納。理論上 正しい形式のRPN式であれば最終的にはstackには1要素だけ残るはずなので、その値を確認すれば良い。

これでとりあえず、(98の間の割り算を無視して)計算結果が 2023 になる式の組み合わせを列挙できる。

let nums = (1..=10).rev().collect();
let mut results = Vec::new();
DefaultSearcher::new(nums, 2023).traverse(&mut Vec::new(), (0, 0), &mut results);

手元の環境で実行すると 3370 件の結果が約2分弱で列挙された。

...

10-(9-(((8*(7*6))+(5-4))*(3*(2*1))))
10-(9-(((8*(7*6))+(5-4))*(3*(2/1))))
10-(9-(((8*(7*6))+(5-4))*(3+(2+1))))
10-(9-(8*(((7*(6+(5/4)))*(3+2))-1)))
10-(9-(8*((7*((6+(5/4))*(3+2)))-1)))
Completed 3370 results in 111.933449708s

探索の高速化

計算結果を保持

前述の実装だと、式が完成するたびにstack用意して計算結果を確認して… というのが非常に無駄になる。解説記事によると 1,274,544,128 通り生成されるようなので、ここはもっと計算量を削りたい。

途中までの計算結果は変わらないものなので、それを保持しながら探索した方が効率的になる。 内部に別のstackを持って、 Operator を受け取った時点で演算処理結果を格納しておくようにしておく。

前述の self.evaluate() でやっていた計算を self.backtrack() の中で逐次やっておく、という感じ。

pub struct FastSearcher {
    nums: Vec<i32>,
    target: i32,
    stack: Vec<Fraction>,
}

impl FastSearcher {
    pub fn new(nums: Vec<i32>, target: i32) -> Self {
        let stack = Vec::with_capacity(nums.len());
        Self {
            nums,
            target,
            stack,
        }
    }
}

impl Rpn for FastSearcher {
    fn evaluate(&self, _: &[ExpressionElement]) -> bool {
        self.stack[0].1 * self.target == self.stack[0].0
    }
    fn get(&self, i: usize) -> Option<&i32> {
        self.nums.get(i)
    }
    fn backtrack(
        &mut self,
        expr: &mut Vec<ExpressionElement>,
        (i, j): (usize, usize),
        results: &mut Vec<Vec<ExpressionElement>>,
        e: ExpressionElement,
    ) {
        match e {
            ExpressionElement::Operand(n) => {
                self.stack.push(Fraction(n, 1));
                expr.push(e);
                self.traverse(expr, (i, j), results);
                expr.pop();
                self.stack.pop();
            }
            ExpressionElement::Operator(op) => {
                if let (Some(n0), Some(n1)) = (self.stack.pop(), self.stack.pop()) {
                    if let Some(n) = op.apply(&n1, &n0) {
                        self.stack.push(n);
                        expr.push(e);
                        self.traverse(expr, (i, j), results);
                        expr.pop();
                        self.stack.pop();
                    }
                    self.stack.push(n1);
                    self.stack.push(n0);
                }
            }
        }
    }
}

これによって self.evaluate() に到達した時点で計算結果は出ているので、その値だけを確認すれば良い、ということになる。また、途中でゼロ除算が生じた場合にその先の探索をしなくなるので無駄を省く効果もありそうだ。

これで、前述の 3370 件の列挙が 約30秒に短縮された。約4倍の高速化。

探索結果のキャッシュ

また、この内部stackは同じ状態を経過することも多い。例えば

  • 10*(9-8)+...
  • 10/(9-8)+...

8 までの計算結果はどちらも変わらない。どちらか一方で全探索して1件も見つからなかったなら、もう片方でも1件も見つからないことが分かる。

ので、探索して条件に合致する式が見つかったか否かを返すようにし、それが false だった場合には同様の状態からの探索をスキップするように変更。

use std::collections::HashSet;

type Key = (Vec<Fraction>, (usize, usize));

pub struct FastSearcher {
    nums: Vec<i32>,
    target: i32,
    stack: Vec<Fraction>,
    seen: HashSet<Key>,
}

impl FastSearcher {
    pub fn new(nums: Vec<i32>, target: i32) -> Self {
        let stack = Vec::with_capacity(nums.len());
        Self {
            nums,
            target,
            stack,
            seen: HashSet::new(),
        }
    }
}

impl Rpn for FastSearcher {
    fn traverse(
        &mut self,
        expr: &mut Vec<ExpressionElement>,
        (i, j): (usize, usize),
        results: &mut Vec<Vec<ExpressionElement>>,
    ) -> bool {
        if i == j + 1 && self.get(i).is_none() {
            let found = self.evaluate(expr);
            if found {
                results.push(expr.clone());
            }
            return found;
        }
        let key = (self.stack.clone(), (i, j));
        if self.seen.contains(&key) {
            return false;
        }
        let mut ret = false;
        if let Some(&n) = self.get(i) {
            ret |= self.backtrack(expr, (i + 1, j), results, ExpressionElement::Operand(n));
        }
        if i > j + 1 {
            for op in &Op::ALL {
                ret |= self.backtrack(expr, (i, j + 1), results, ExpressionElement::Operator(*op));
            }
        }
        if !ret {
            self.seen.insert(key);
        }
        ret
    }
}

これで、約12秒くらいに短縮された。さらに2倍以上速くなった。

実際にはこの key に使っている self.stack は約分されていない計算結果なので、同じ状態を経由しているはずなのに 値が異なっているという扱いになってしまうケースがある。 gcd を使って既約分数を格納するようにしてみる。

impl FastSearcher {
    fn gcd(a: i32, b: i32) -> i32 {
        if b == 0 {
            a
        } else {
            Self::gcd(b, a % b)
        }
    }
}

impl Rpn for FastSearcher {
    ...

    fn backtrack(
        &mut self,
        expr: &mut Vec<ExpressionElement>,
        (i, j): (usize, usize),
        results: &mut Vec<Vec<ExpressionElement>>,
        e: ExpressionElement,
    ) -> bool {
        let mut ret = false;
        match e {
            ExpressionElement::Operand(n) => {...}
            ExpressionElement::Operator(op) => {
                if let (Some(n0), Some(n1)) = (self.stack.pop(), self.stack.pop()) {
                    if let Some(n) = op.apply(&n1, &n0) {
                        let gcd = Self::gcd(n.0, n.1);
                        self.stack.push(Fraction(n.0 / gcd, n.1 / gcd));

                        ...
                    }
                    self.stack.push(n1);
                    self.stack.push(n0);
                }
            }
        }
        ret
    }
}

これで、3370 件の列挙が 約7.5秒程度まで短縮された。 このキャッシュ手法は探索目標の値などによって効率が変わってくるので一概には言えないかもしれないが、少なくともこの例においては数倍の高速化が実現できた。

割り算を考慮した条件つき生成

ここまでは計算結果が 2023 になるものを全列挙することだけを考えていたが、件のパズル問題では「98の間には既に÷が入っています」という制約がある。

前述した 3370 件から条件に合致するものを抽出しても良いが、そもそも探索する時点で条件に合致しない式になるものを除外していくほうが効率が良い。

これは自分ではまったく思いつかなくて解説記事を読んで目から鱗だったのだけど、 8の数字が出た後、数字と演算子の数が一致していた場合」の演算子を割り算に固定 することによって実現できる。

解説記事を参考に、専用のSearcherを実装してみる。

struct Div8Searcher {
    nums: Vec<i32>,
    target: i32,
    stack: Vec<Fraction>,
    eight_depth: i32,
}

impl Div8Searcher {
    pub fn new(nums: Vec<i32>, target: i32) -> Self {
        let stack = Vec::with_capacity(nums.len());
        Self {
            nums,
            target,
            stack,
            eight_depth: -1,
        }
    }
}

impl Rpn for Div8Searcher {
    fn traverse(
        &mut self,
        expr: &mut Vec<ExpressionElement>,
        (i, j): (usize, usize),
        results: &mut Vec<Vec<ExpressionElement>>,
    ) -> bool {

        ...

        if i > j + 1 {
            for op in &Op::ALL {
                if self.eight_depth == 0 && op != &Op::Div {
                    continue;
                }
                ret |= self.backtrack(expr, (i, j + 1), results, ExpressionElement::Operator(*op));
            }
        }
        ret
    }
    fn backtrack(
        &mut self,
        expr: &mut Vec<ExpressionElement>,
        (i, j): (usize, usize),
        results: &mut Vec<Vec<ExpressionElement>>,
        e: ExpressionElement,
    ) -> bool {
        let mut ret = false;
        match e {
            ExpressionElement::Operand(n) => {
                let orig = self.eight_depth;
                if n == 8 {
                    self.eight_depth = 0;
                } else if self.eight_depth >= 0 {
                    self.eight_depth += 1;
                }
                self.stack.push(Fraction(n, 1));
                expr.push(e);
                ret |= self.traverse(expr, (i, j), results);
                expr.pop();
                self.stack.pop();
                self.eight_depth = orig;
            }
            ExpressionElement::Operator(op) => {
                if let (Some(n0), Some(n1)) = (self.stack.pop(), self.stack.pop()) {
                    if let Some(n) = op.apply(&n1, &n0) {
                        self.eight_depth -= 1;
                        self.stack.push(n);
                        expr.push(e);
                        ret |= self.traverse(expr, (i, j), results);
                        expr.pop();
                        self.stack.pop();
                        self.eight_depth += 1;
                    }
                    self.stack.push(n1);
                    self.stack.push(n0);
                }
            }
        }
        ret
    }
}

eight_depth の値を保持し、8 が出現したら 0 にセットして OperandOperator かによって値を増減。値が 0 のときには次に push する OpeatorOp::Div だけに制限される。

前述したような探索結果のキャッシュも使いたいところだが、新たな状態としてこの self.eight_depth もキーに含める必要があり、その割にはそれほど効果も無さそうなのでここでは省いた。

これによって 「98の間には必ず/が入る」式だけが探索されるようになる。

条件に合致する 530 件が、約4秒強で列挙されるようになった。全列挙してから抽出するよりも明らかに速いことが分かる。

...

(10*(9/((8-7)/((6*5)/(4/3)))))-(2/1)
(10*(9/((8-7)/(6*((5/4)*3)))))-(2*1)
(10*(9/((8-7)/(6*((5/4)*3)))))-(2/1)
(10*(9/((8-7)/(6*(5/(4/3))))))-(2*1)
(10*(9/((8-7)/(6*(5/(4/3))))))-(2/1)
Completed 530 results in 4.187013041s

中置記法、正規化

ここまでで探索して収集された resultsVec<Vec<ExpressionElement>> であり、これを中置記法に書き換える操作もまたstackを使って実現していくことになる。

impl Display for Op {
    fn fmt(&self, f: &mut Formatter<'_>) -> Result {
        f.write_char(match self {
            Self::Add => '+',
            Self::Sub => '-',
            Self::Mul => '*',
            Self::Div => '/',
        })
    }
}

fn solve(searcher: &mut impl Rpn) {
    let mut results = Vec::new();
    searcher.traverse(&mut Vec::new(), (0, 0), &mut results);
    for result in results {
        let mut stack = Vec::new();
        for e in &result {
            match e {
                ExpressionElement::Operand(n) => stack.push((n.to_string(), None)),
                ExpressionElement::Operator(op) => {
                    if let (Some((mut s0, o0)), Some((mut s1, o1))) = (stack.pop(), stack.pop()) {
                        if o1.is_some() {
                            s1 = format!("({s1})");
                        }
                        if o0.is_some() {
                            s0 = format!("({s0})");
                        }
                        stack.push((format!("{s1}{op}{s0}"), Some(*op)));
                    }
                }
            }
        }
        println!("{}", stack[0].0);
    }
}

連結していくための String と、どの Operator で得られたものか(もしくは Operand か)、を保持しながら処理していく。簡単にすべて括弧をつけてしまうのなら上述のように書ける。

不要な括弧を除去して出力する場合は、左operandの場合と右operandの場合で処理が変わる(引き算・割り算は右operandに注意する必要がある)。matches! で判定してやってみる。

        ...

        ExpressionElement::Operator(op) => {
            if let (Some((mut s0, o0)), Some((mut s1, o1))) = (stack.pop(), stack.pop()) {
                if matches!((op, o1), (Op::Mul | Op::Div, Some(Op::Add | Op::Sub))) {
                    s1 = format!("({s1})");
                }
                if matches!(
                    (op, o0),
                    (Op::Sub | Op::Mul, Some(Op::Add | Op::Sub)) | (Op::Div, Some(_))
                ) {
                    s0 = format!("({s0})");
                }
                stack.push((format!("{s1}{op}{s0}"), Some(*op)));
            }
        }

これで、得られた結果を HashSet などで重複排除することで、解説記事と同様の 81 件を約4秒強で列挙することができた。

((10+(9/8+7)*6*5)*4-3)*2-1
(10*9/((8-7)/(6*5))/(4/3)-2)*1
(10*9/((8-7)/(6*5))/(4/3)-2)/1
(10*9/((8-7)/(6*5))/4*3-2)*1
(10*9/((8-7)/(6*5))/4*3-2)/1
(10*9/((8-7)/(6*5)*4)*3-2)*1
(10*9/((8-7)/(6*5)*4)*3-2)/1
(10*9/((8-7)/(6*5)*4/3)-2)*1
(10*9/((8-7)/(6*5)*4/3)-2)/1
(10*9/((8-7)/(6*5/(4/3)))-2)*1
(10*9/((8-7)/(6*5/(4/3)))-2)/1
(10*9/((8-7)/(6*5/4))*3-2)*1
(10*9/((8-7)/(6*5/4))*3-2)/1
(10*9/((8-7)/(6*5/4)/3)-2)*1
(10*9/((8-7)/(6*5/4)/3)-2)/1
(10*9/((8-7)/(6*5/4*3))-2)*1
(10*9/((8-7)/(6*5/4*3))-2)/1
(10*9/((8-7)/6)*5/(4/3)-2)*1
(10*9/((8-7)/6)*5/(4/3)-2)/1
(10*9/((8-7)/6)*5/4*3-2)*1
(10*9/((8-7)/6)*5/4*3-2)/1
(10*9/((8-7)/6/(5/(4/3)))-2)*1
(10*9/((8-7)/6/(5/(4/3)))-2)/1
(10*9/((8-7)/6/(5/4))*3-2)*1
(10*9/((8-7)/6/(5/4))*3-2)/1
(10*9/((8-7)/6/(5/4)/3)-2)*1
(10*9/((8-7)/6/(5/4)/3)-2)/1
(10*9/((8-7)/6/(5/4*3))-2)*1
(10*9/((8-7)/6/(5/4*3))-2)/1
(10*9/((8-7)/6/5)/(4/3)-2)*1
(10*9/((8-7)/6/5)/(4/3)-2)/1
(10*9/((8-7)/6/5)/4*3-2)*1
(10*9/((8-7)/6/5)/4*3-2)/1
(10*9/((8-7)/6/5*4)*3-2)*1
(10*9/((8-7)/6/5*4)*3-2)/1
(10*9/((8-7)/6/5*4/3)-2)*1
(10*9/((8-7)/6/5*4/3)-2)/1
(10*9/(8-7)*6*5/(4/3)-2)*1
(10*9/(8-7)*6*5/(4/3)-2)/1
(10*9/(8-7)*6*5/4*3-2)*1
(10*9/(8-7)*6*5/4*3-2)/1
10*9/((8-7)/(6*5))/(4/3)-2*1
10*9/((8-7)/(6*5))/(4/3)-2/1
10*9/((8-7)/(6*5))/4*3-2*1
10*9/((8-7)/(6*5))/4*3-2/1
10*9/((8-7)/(6*5)*4)*3-2*1
10*9/((8-7)/(6*5)*4)*3-2/1
10*9/((8-7)/(6*5)*4/3)-2*1
10*9/((8-7)/(6*5)*4/3)-2/1
10*9/((8-7)/(6*5/(4/3)))-2*1
10*9/((8-7)/(6*5/(4/3)))-2/1
10*9/((8-7)/(6*5/4))*3-2*1
10*9/((8-7)/(6*5/4))*3-2/1
10*9/((8-7)/(6*5/4)/3)-2*1
10*9/((8-7)/(6*5/4)/3)-2/1
10*9/((8-7)/(6*5/4*3))-2*1
10*9/((8-7)/(6*5/4*3))-2/1
10*9/((8-7)/6)*5/(4/3)-2*1
10*9/((8-7)/6)*5/(4/3)-2/1
10*9/((8-7)/6)*5/4*3-2*1
10*9/((8-7)/6)*5/4*3-2/1
10*9/((8-7)/6/(5/(4/3)))-2*1
10*9/((8-7)/6/(5/(4/3)))-2/1
10*9/((8-7)/6/(5/4))*3-2*1
10*9/((8-7)/6/(5/4))*3-2/1
10*9/((8-7)/6/(5/4)/3)-2*1
10*9/((8-7)/6/(5/4)/3)-2/1
10*9/((8-7)/6/(5/4*3))-2*1
10*9/((8-7)/6/(5/4*3))-2/1
10*9/((8-7)/6/5)/(4/3)-2*1
10*9/((8-7)/6/5)/(4/3)-2/1
10*9/((8-7)/6/5)/4*3-2*1
10*9/((8-7)/6/5)/4*3-2/1
10*9/((8-7)/6/5*4)*3-2*1
10*9/((8-7)/6/5*4)*3-2/1
10*9/((8-7)/6/5*4/3)-2*1
10*9/((8-7)/6/5*4/3)-2/1
10*9/(8-7)*6*5/(4/3)-2*1
10*9/(8-7)*6*5/(4/3)-2/1
10*9/(8-7)*6*5/4*3-2*1
10*9/(8-7)*6*5/4*3-2/1
Completed 81 results in 4.229359166s

CLIアプリケーション化

ここまで幾つかの実装ができたので、折角なので引数でパラメータを変えながら全列挙できる CLI アプリケーションを作成してみた。

$ ./puzzle2023 --help
Solvers for https://twitter.com/tkihira/status/1609313732034965506

Usage: puzzle2023 [OPTIONS]

Options:
  -t, --target <TARGET>      [default: 2023]
  -m, --max-num <MAX_NUM>    [default: 10]
  -r, --rev <REV>            [default: true] [possible values: true, false]
  -s, --searcher <SEARCHER>  [default: fast] [possible values: default, fast, div8]
  -n, --normalize
  -v, --verbose
  -h, --help                 Print help information
  -V, --version              Print version information

全パターン列挙

$ ./puzzle2023 -v

...

10-(9-(((8*(7*6))+(5-4))*(3*(2*1))))
10-(9-(((8*(7*6))+(5-4))*(3*(2/1))))
10-(9-(((8*(7*6))+(5-4))*(3+(2+1))))
10-(9-(8*(((7*(6+(5/4)))*(3+2))-1)))
10-(9-(8*((7*((6+(5/4))*(3+2)))-1)))
Completed 3370 results in 7.338064208s

9と8の間の割り算に限定

$ ./puzzle2023 -v --searcher div8

...

(10*(9/((8-7)/(6*((5/4)*3)))))-(2/1)
(10*(9/((8-7)/(6*(5/(4/3))))))-(2*1)
(10*(9/((8-7)/(6*(5/(4/3))))))-(2/1)
Completed 530 results in 4.219031208s

正規化

$ ./puzzle2023 -v --searcher div8 --normalize

...

10*9/(8-7)*6*5/(4/3)-2/1
10*9/(8-7)*6*5/4*3-2*1
10*9/(8-7)*6*5/4*3-2/1
Completed 81 results in 4.239637041s

計算結果の確認

$ ./puzzle2023 --searcher div8 --normalize | python3 -c 'for s in open(0): print(f"{eval(s):.1f}")' | sort | uniq -c
  81 2023.0

使うのを 1 から 8 の昇順にして 2024 を作るようにしてみる

$ ./puzzle2023 -m8 -rfalse -t2024 -nv
(((1+2)*3*4+5)*6+7)*8
(1*2*(3+4*5*6)+7)*8
(1+(2+3*4)*(5+6+7))*8
(1+(2+3-(4-5))*6*7)*8
(1+(2+3-4+5)*6*7)*8
(1+2*(3+4)*(5+6+7))*8
(1+2*(3+4+5+6)*7)*8
(1+2*3+4)*(5*6-7)*8
(1+2/(3/((4+5)*6))*7)*8
(1+2/(3/((4+5)*6)/7))*8
(1+2/(3/((4+5)*6*7)))*8
(1+2/(3/(4+5))*6*7)*8
(1+2/(3/(4+5)/(6*7)))*8
(1+2/(3/(4+5)/6)*7)*8
(1+2/(3/(4+5)/6/7))*8
(1+2/3*(4+5)*6*7)*8
(1-(2-3*4))*(5*6-7)*8
(1-2*(3*4-5*6)*7)*8
(1-2*3*(4-5)*6*7)*8
(1-2*3/((4-5)/(6*7)))*8
(1-2*3/((4-5)/6)*7)*8
(1-2*3/((4-5)/6/7))*8
(1-2*3/(4-5)*6*7)*8
(1-2+3*4)*(5*6-7)*8
1*(2*(3+4*5*6)+7)*8
Completed 25 results in 131.605458ms

「8の前を割り算に限定」にしているのも、コマンドラインオプションで違う値で指定できるようにしたかったけど、複雑になりすぎる感じがしたのでそれは断念…

改善?

tkihiraさんの最適化ではstackを push / pop するのではなく固定長配列を用意して見るべきindexだけをズラしていく、といった工夫をしていた。真似して実装してみたが本実装においてはそこはあまり効果なさそうだった。

あとはこうして再帰的に探索した結果を一度 results という Vec に格納してしまっているのが実は非効率なのでは、と思ったが… Iterator として使いたかったが そうするためには探索の途中の状態を保持しておいて next() が呼ばれるたびに続きを探索していく、という処理にする必要があり、そのあたりの実装がうまくできなかった。何か良い方法あるのだろうか…

Repository

github.com

Rubyでバイナリデータに対するrindex検索の挙動でハマったので調べたことメモ

自分の手元の環境でこんなことが起きた。

$ ruby -v
ruby 3.1.2p20 (2022-04-12 revision 4491bb740a) [arm64-darwin21]
$ irb
irb(main):001:0> "\x01\x80\x00\x00".index("\x01")
=> 0
irb(main):002:0> "\x01\x80\x00\x00".rindex("\x01")
=> 1

\x010 番目にしかないのだから、 .index でも .rindex でも 0 が返ってくるはずではないの??

先に結論

バイナリデータを扱うときには必ずEncodingを ASCII-8BIT に指定しておくこと!

きっかけ

roo というgem (記事時点で 2.9.0)を使って、Excelファイルを開こうとした。

require 'roo'

p Roo::Excelx.new("hoge.xlsx")

これは問題ないが、 #initialize の引数は file_or_strem ということでファイルを open したものを渡しても良いはず、と

require 'roo'

File.open("hoge.xlsx") do |f|
  p Roo::Excelx.new(f)
end

ということをすると以下のような謎のエラーを吐いて落ちる。

/Users/sugyan/.rbenv/versions/3.1.2/lib/ruby/gems/3.1.0/gems/rubyzip-2.3.2/lib/zip/entry.rb:365:in `check_c_dir_entry_static_header_length': undefined method `bytesize' for nil:NilClass (NoMethodError)

      return if buf.bytesize == ::Zip::CDIR_ENTRY_STATIC_HEADER_LENGTH
                   ^^^^^^^^^

で困っていたので他の人に聞いてみたところ「自分のところではそんな問題は起きていない」と言われ。 試しにDockerでRuby環境作って実行してみると、確かに問題なく動く…何故? と調べはじめた。

String#rindex の謎挙動

自分の環境と問題なく動く環境とで比較しながら見てみると、roo が使っている rubyzip (記事時点で 2.3.2) の中でファイル内容を読み取った結果の値が違っていた。

module Zip
  class CentralDirectory
    include Enumerable

    END_OF_CDS             = 0x06054b50

    ...

    def get_e_o_c_d(buf) #:nodoc:
      sig_index = buf.rindex([END_OF_CDS].pack('V'))
      ...

buf の中から 0x06054b50 をpackしたもの つまり "PK\x05\x06" のデータを末尾から探すために .rindex() を呼んでいるのだが、手元の環境と 問題なく動く環境ではその値が 1 ズレている。。

何故?と思って色々試しているうちに冒頭の例のようなものに辿り着いた。

もう少し深く追う

少なくとも \x80 などを含まないASCIIだけのものであればおかしな挙動にはならない。

irb(main):001:0> "\x01\x7F\x00\x00".rindex("\x01")
=> 0

また、これが増えるとズレもどんどん大きくなる。

irb(main):001:0> "\x01\x80".rindex("\x01")
=> 0
irb(main):002:0> "\x01\x80\x80".rindex("\x01")
=> 1
irb(main):003:0> "\x01\x80\x80\x80".rindex("\x01")
=> 2
irb(main):004:0> "\x01\x80\x80\x80\x80".rindex("\x01")
=> 3

逆から走査する際に異常な文字が検出された際に読み飛ばし、その結果の走査数を文字列長から引いたものを返したことにより値がおかしくなる、というのが予想される。

ソースを読んでみる。 rstring 関連はこのへんのようだ。

#ifdef HAVE_MEMRCHR によって実装が分かれている。 memrchr というのがglibcに入っているもので、自分のmacOS環境では使えなくて Docker内のLinux環境などでは使える、これによって環境によって挙動が変わったりする、のかもしれない。

で、使えない環境の方での実装は。 対象が見つかるまで while loop の中で rb_enc_prev_char で前の文字を走査している、という感じのようだ。encodingの影響を受けそう…?

    while (s) {
        if (memcmp(s, t, slen) == 0) {
            return pos;
        }
        if (pos == 0) break;
        pos--;
        s = rb_enc_prev_char(sbeg, s, e, enc);
    }

Encodingと実行環境

Encodingが関係しているようだったので色々調べてみた。

Rubyではバイナリデータ(バイト列)も String で扱う。

バイナリの取扱い

Ruby の String は、文字の列を扱うためだけでなく、バイトの列を扱うためにも使われます。しかし、Ruby M17N には直接にバイナリを表すエンコーディングは存在しません。このため、バイナリを String で扱う際には、ASCII 互換オクテット列を意味する ASCII-8BIT を用います。これにより、ASCII 互換であるこの String は 7bit クリーンな文字列と比較・結合が可能となります。

https://docs.ruby-lang.org/ja/latest/doc/spec=2fm17n.html

ということで、バイナリの場合は本来は ASCII-8BIT のEncodingを持つ文字列として検索しなければならないのに、手元の環境では UTF-8 のままで検索しようとしていた、というところに「使い方の問題」があったようだ。

irb(main):001:0> "\x01\x80\x00\x00".encoding
=> #<Encoding:UTF-8>

Encodingが ASCII-8BIT になっていれば、rindex でも正しい値を得ることができそうだ。 全然知らなかったが String#bASCII-8BIT の複製を得ることもできるらしい。

https://docs.ruby-lang.org/ja/latest/method/String/i/b.html

irb(main):001:0> "\x01\x80\x00\x00".rindex("\x01")
=> 1
irb(main):002:0> "\x01\x80\x00\x00".force_encoding(Encoding::ASCII_8BIT).rindex("\x01")
=> 0
irb(main):003:0> "\x01\x80\x00\x00".b.rindex("\x01")
=> 0

なるほど〜。

ではこの Encoding は何で決まるのか?というのも上記リンクに書いてある。

リテラルエンコーディング

文字列リテラル正規表現リテラルそしてシンボルリテラルから生成されるオブジェクトのエンコーディングスクリプトエンコーディングになります。

またスクリプトエンコーディングが US-ASCII である場合、7bit クリーンではないバックスラッシュ記法で表記されたリテラルエンコーディングは ASCII-8BIT になります。

さらに Unicode エスケープ (\uXXXX) を含む場合、リテラルエンコーディングUTF-8 になります。

複雑。。。

とにかく通常はスクリプトエンコーディングがまず大事のようだ。

スクリプトエンコーディング

スクリプトエンコーディングとは Ruby スクリプトを書くのに使われているエンコーディングです。スクリプトエンコーディングは マジックコメントを用いて指定します。スクリプトエンコーディングには ASCII 互換エンコーディングを用いることができます。 ASCII 非互換のエンコーディングや、ダミーエンコーディングは用いることができません。

現在のスクリプトエンコーディング__ENCODING__ により取得することができます。

さらに、magic commentによってこのスクリプトエンコーディングを決定でき、それが無い場合の挙動も書かれている。

マジックコメントが指定されなかった場合、コマンド引数 -K, RUBYOPT およびファイルの shebang からスクリプトエンコーディングは以下のように決定されます。左が優先です。

magic comment(最優先) > -K > RUBYOPTの-K > shebang

上のどれもが指定されていない場合、通常のスクリプトなら UTF-8、-e や stdin から実行されたものなら locale がスクリプトエンコーディングになります。 -K オプションが複数指定されていた場合は、後のものが優先されます。

通常のスクリプト-e かどうか、でも変わったりするのね…。そして特に指定ない場合は最終的にはlocaleが使われる。 ので、 LC_CTYPE などでも挙動が変わり得るようだ。

$ ruby -e 's="\x01\x80\x00\x00"; p __ENCODING__, s.encoding, s.rindex("\x01")'
#<Encoding:UTF-8>
#<Encoding:UTF-8>
1
$ ruby -Kn -e 's="\x01\x80\x00\x00"; p __ENCODING__, s.encoding, s.rindex("\x01")'
#<Encoding:ASCII-8BIT>
#<Encoding:ASCII-8BIT>
0
$ RUBYOPT="-Kn" ruby -e 's="\x01\x80\x00\x00"; p __ENCODING__, s.encoding, s.rindex("\x01")'
#<Encoding:ASCII-8BIT>
#<Encoding:ASCII-8BIT>
0
$ LC_CTYPE=C ruby -e 's="\x01\x80\x00\x00"; p __ENCODING__, s.encoding, s.rindex("\x01")'
#<Encoding:US-ASCII>
#<Encoding:ASCII-8BIT>
0

という具合に、何もしないと UTF-8 文字列として扱ってしまうが、オプションや環境変数によって リテラルの Encoding を変えることもできる。

こうやって正しく ASCII-8BIT のEncodingを持つ文字列として検索すれば、 String#rindex の値がズレるということはなさそうだ。

つまり再現条件は

手元の環境で String#rindex の値がズレたのは2つの要因があって

  1. memrchr が使えない環境でビルドしたRubyで実行していて
  2. ASCII-8BIT でないEncodingの文字列に対して rindex をかけていた

ということになる。2つ揃っていないと起きないものと思われる。

Rooの問題

前述の、Rooでエラーが起きる件について考える。

そもそもEncodingが不明で渡ってくるものに対して安易に rindex を使うものではない、ということで rubyzip 側に落ち度がありそうではある。が 現在の master branch では リリースされている 2.3.2 とは大きく変わっていて、もう同じことは起きないのかもしれない。詳しくは追っていないので分からないが、今回のと関連するissueがあった。

既にCloseされているが、今回調べた Zip::CentralDirectory は直接使用するものではなく Zip::File を使ってください、ということのようだ。

つまりこの場合、 Roo 側の使い方が悪い、ということになりそう。 で、Rooの方も調べてみると ちゃんとそれに対応しようとしていると思われる修正があった。

これが取り込まれれば手元の環境でも問題なく動くようになるかな…?

もしくは現時点でも、使う側が渡す引数のEncodingを明示的に指定してやることで一応回避できそうではある。

require 'roo'

File.open("hoge.xlsx", "r:ASCII-8BIT") do |f|
  p Roo::Excelx.new(f)
end

Rubyのバグではないの?

しかし Encodingを正しく指定していなかったために起きているとしても、 String#indexString#rindex も検索した対象の開始indexを返すものであるはずなのだから、それがズレるのはおかしいのではないのか…?という気はする。

さらに memrchr が使えるか否かのビルド環境によってだけで結果が変わる(ちゃんと確かめてないけど おそらくそう…)、というのも気持ち悪い。

かなり昔から現在の実装になっているようだし 今からでは挙動を変えづらそうではある。 仕様といってしまえばそれまでだけど…。 このあたりはRuby開発陣の方々の見解も聞いてみたいところではあります。

3.2

ちなみに Ruby 3.2 からは String#byteindexString#byterindex が追加されるそうです。今回のようなバイナリデータに対する検索にピッタリ使えそうですね。