🕸️

Rustで動的計画法の実装：Grid

2023/04/13に公開

Rust

AtCoder

動的計画法

tech

はじめに

動的計画法を実装してみて、Rustの勉強をやってみる。問題としてはEducational DP Contestという動的計画法の練習を目的としたコンテストのものを使用。

https://atcoder.jp/contests/dp/

AからZまで問題が設定されているが、今回はHのGrid1とYのGrid2、壁のあるグリッドでの最短経路の数を求める問題。同じ問題であるがGrid1はグリッドのサイズに制限があり、壁の数に制限がない。Grid2は壁の数に制限がある代わりにグリッドのサイズの制限がGrid1より緩い。このため、解き方が異なる。

利用したライブラリ

標準入力からデータを読み取る部分はproconioを利用。入力のノード番号が1オリジン(1-indexed)なので、0オリジン(0-indexed)に変更するためにmarker::Usize1も使う。

Grid1: 完成したコード

動的計画法っぽく、順番を意識しながら(0,0)から数え上げれば解ける。オーバーフローを避けるためにMODで割った余りを出力するという指定があるので、計算の度に引いている。

grid_counting.rs

use proconio::input;
use proconio::marker::Chars;

const MOD:i64 = 1000000007;

fn main(){
    input!(
        h: usize, w: usize, n:usize,
        grid: [Chars; n]
    );

    let mut dp: Vec<Vec<i64>> = vec![vec![0; w]; h];
    dp[0][0] = 1;
    for i in 0..h {
        for j in 0..w {
            if i < h - 1 && grid[i+1][j] == '.' { dp[i+1][j] += dp[i][j]; if dp[i+1][j] >= MOD { dp[i+1][j] -= MOD; } }
            if j < w - 1 && grid[i][j+1] == '.' { dp[i][j+1] += dp[i][j]; if dp[i][j+1] >= MOD { dp[i][j+1] -= MOD; } }
        }
    }
    println!("{}", dp[h-1][w-1]);
}
//

Rustの&&はLazy処理

Rustの&&はLazy Operatiorである。つまり、&&の前の項がfalseであれば後ろの項は処理されない。このため、以下のように配列のindexのチェックと配列の参照を&&で繋げるとエラーは発生しない。

grid_counting.rs

            if i < h - 1 && grid[i+1][j] == '.' { dp[i+1][j] += dp[i][j]; if dp[i+1][j] >= MOD { dp[i+1][j] -= MOD; } }

つまり、以下のような処理と同じになることが保証されている。

grid_counting.rs

            if i < h - 1 {
                if grid[i+1][j] == '.' {

Perlでは&&, ||と別にLazyな処理が保証された演算子であるandとorが導入されていた。以下のような書き方はPerlの定番である。

orのサンプル

fopen(INFILE, "foo.txt") or die "Cannot open file: $!";

昔は実装依存とかそんな話があったような気がするが、Wikipediaによると、いまどきはだいたいLazy Operatorsであるようだ。

Grid2: 解いてみる

最初Gird2をみた時、入力フォーマットの違いだけかと思った。つまり、以下のようにフォーマットを変換してあげれば解けるのではないかと。実際は入力例の4、 $10^5\times 10^5$ のグリッドの問題において処理が終わらない。

grid_counting_2.rs

fn main(){
    input!(
        h: usize, w: usize, n:usize,
        wall: [(usize, usize); n]
    );

    let mut grid = vec![vec!['.'; w]; h];
    for (x, y) in wall {
        grid[x - 1][y - 1] = '#';
    }

包除原理で解く

数え上げる方法での計算オーダーは $O(WH)$ 、つまりグリッドの面積に比例してしまうので、計算が終わらない。そこで、組み合わせ計算と包除原理で解くのが正しいらしい。この方法だと、壁の数 $N$ に対して $O(N^2)$ の計算量で済む。

剰余演算も必要なので実装も面倒ではあるが、AtCoder Liburary ac_libraryに便利なライブラリが用意されているので、これを使って実装した。

Closure

Closureは無名関数だが、scope内の変数を参照できるので、インラインマクロのようにも使える。
階乗をfactorialとfactorial_invで予め計算して、それを使ってconbinationを計算する関数combination(a,b)を以下のように定義している。

grid_pie.rs

let combination = |a, b| factorial[a] * factorial_inv[b] * factorial_inv[a - b];

two-phase borrow

以下の部分。+=を使って書きたくなるがエラーがでてしまう。dp[i][j]がModInt<1000000007>型であるのが問題のようで、これまでi64やusizeのような組込の型では問題にならなかった。two-phase borrowが適用されるかどうかという話らしい。以下のように展開してしまえば問題ない。

grid_pie.rs

dp[j][m] = dp[j][m] + dp[i][1 - m] * combination(dx + dy, dx);

詳しくは以下の文章が分かりやすい。

計算量の測定

問題を作成するスクリプトを作って計算量の測定もやってみた。入力値のチェックをしていないので、グリッドサイズに対して大量の壁を要求すると、処理が終わらないとは思うが、テスト用なので問題はない。

gen_grid.py

#!/user/bin/env python
import random
import sys

(h, w, n) = (int(sys.argv[1]), int(sys.argv[2]), int(sys.argv[3]))

s = set()
while len(s) < n:
    pos = (random.randint(1,h), random.randint(1,w))
    if pos != (1,1) and pos != (h, w) and pos not in s:
        s.add(pos)

print('{} {} {}'.format(h, w, n))
for (y, x) in s:
    print('{} {}'.format(y,x))

数え上げ（counting）と包除原理（PIE; Principle of Inclusion-Exclusion) での計算時間の比較。PIEの方は階乗計算を予め計算しており、その時間も含む。

上は壁の数を10個に固定して $10 \times 10$ から $10^4\times 10^4$ まで変化させた時の計算時間。Countingは線形に延びているが $10^4\times 10^4$ の時に1秒になっている。このため、単純に線形で推定すると $10^5\times 10^5$ の時は100秒になる。PIEはグリッドサイズに依存しないので直線になる。

下のグラフはグリットサイズを $10^3 \times 10^3$ に固定して、壁の数 $N$ を $10$ から $10^4$ まで変化させた時の計算時間。PIEは $O(N^2)$ の曲線になっているので $N$ が1000を超えるとCountingより計算量が大きくなってしまう。二乗曲線なので、 $N$ が大きくなるとさらに問題になるかと思う。

PIEでは今回、階乗だけでなく階乗の逆元も予め求めており、Conbinationを求める時は掛け算だけで処理できるようにしている。当然、元の式である以下のようにしても計算可能である。

combination

let combination = |a, b| factorial[a] / (factorial[b] * factorial[a - b]);

逆元を予め求める場合（PIE）と予め求めない場合（PIE not inv)の比較を以下に示す。予め計算するための時間にも差があるはずであるが、それよりも計算量の増加が大きい。メモリ量や計算規模の問題もあるが、逆元は予め計算しておいた方が速い。

Rustで動的計画法の実装：Grid

はじめに

利用したライブラリ

Grid1: 完成したコード

Rustの&&はLazy処理

Grid2: 解いてみる

包除原理で解く

Closure

two-phase borrow

計算量の測定

関連記事

Discussion