Rustで蟻本その2 〜動的計画法〜

の続きです。動的計画法の章から。

ナップサック問題

重さw, 価値vの品物のリストが与えられたとき、総重量をtwに抑えて総価値を最大化にしたときの総価値を求めよ。

再帰関数（メモ化を利用）

fn main() {
    println!(
        "{}",
        Napsack::solve(5, vec![(2, 3), (1, 2), (3, 4), (2, 2)])
    );
}

struct Napsack {
    wv: Vec<(i32, i32)>,
    memo: HashMap<(i32, usize), i32>,
}

impl Napsack {
    fn solve(tw: i32, wv: Vec<(i32, i32)>) -> i32 {
        Napsack {
            wv,
            memo: HashMap::new(),
        }
        .rec(tw, 0)
    }
    fn rec(&mut self, tw: i32, i: usize) -> i32 {
        if let Some(&ans) = self.memo.get(&(tw, i)) {
            ans
        } else if i == self.wv.len() {
            0
        } else {
            let (w, v) = self.wv[i];
            let ans = if tw < w {
                self.rec(tw, i + 1)
            } else {
                let dont_use_i = self.rec(tw, i + 1);
                let use_i = v + self.rec(tw - w, i + 1);
                std::cmp::max(dont_use_i, use_i)
            };
            self.memo.insert((tw, i), ans);
            ans
        }
    }
}

• メモ化をすることで計算量はtw×品物数になる（recの中にループがないことに注目）。
• struct Napsackに情報をまとめてみた。

メモテーブルを計算して解く方法

再帰関数のメモ化とやっていることは結果的に同じだが、コードは見やすくなる

const MAX_N: usize = 1000;
const MAX_W: usize = 1000;
fn dynpro(tw: usize, wv: &[(usize, i32)]) -> i32 {
    // dp[i][x] = wv[..i]を使ってx以下の重さに抑えたときの総価値
    let mut dp = [[0; MAX_W + 1]; MAX_N + 1];
    for (i, &item) in wv.iter().enumerate() {
        for x in 0..=tw {
            dp[i + 1][x] = if x < item.0 {
                dp[i][x]
            } else {
                std::cmp::max(dp[i][x], item.1 + dp[i][x - item.0])
            };
        }
    }
    dp[wv.len()][tw]
}

• 競プロのように変数の範囲が指定されていて、最悪の場合(変数がその最大値をとる場合)を想定すれば良い(変数が実際には小さい値の時の効率化はあまり考慮しなくて良い)ときは、2重配列は事前に固定サイズで与えてしまうのがよさそうだ。

最長共通部分列問題

2つの文字列の共通部分列として最長のものの長さを求めよ。

const MAX_N: usize = 1000;
const MAX_M: usize = 1000;
fn mss(a: &str, b: &str) -> i32 {
    // dp[i][j] = a[..i], b[..j]のmaximal subsequence
    let mut dp = [[0; MAX_M + 1]; MAX_N + 1];
    let a = a.as_bytes();
    let b = b.as_bytes();
    let (n, m) = (a.len(), b.len());
    for i in 0..n {
        for j in 0..m {
            // a[..=i], b[..=j]のmaximal subsequenceを求める
            dp[i + 1][j + 1] = if a[i] == b[j] {
                1 + dp[i][j]
            } else {
                std::cmp::max(dp[i][j + 1], dp[i + i][j])
            };
        }
    }
    dp[n][m]
}

• 計算時間はO(nm)

重複ありナップサック問題

const MAX_N: usize = 1000;
const MAX_W: usize = 1000;
fn napsack2(tw: usize, items: Vec<(usize, i32)>) -> i32 {
    // dp[i][x] = items[..i]を重複ありで使い、総重量≦xとした時の最大価値
    // dp[n][tw]を求めたい。
    let mut dp = [[0; MAX_W]; MAX_N];
    let n = items.len();
    for i in 0..n {
        for x in 0..=tw {
            dp[i + 1][x] = if x < items[i].0 {
                dp[i][x]
            } else {
                // items[i]を使わない、または1つ以上使う
                std::cmp::max(dp[i][x], dp[i + 1][x - items[i].0] + items[i].1)
            }
        }
    }
    dp[n][tw]
}

// 配列を使い回す方法
fn napsack2_2(tw: usize, items: Vec<(usize, i32)>) -> i32 {
    // i回ループがまわった時点で、items[..i]を使った場合の最大価値
    let mut dp = [0; MAX_W];
    for i in 0..items.len() {
        for x in items[i].0..=tw {
            dp[x] = std::cmp::max(dp[x], dp[x - items[i].0] + items[i].1);
        }
    }
    dp[tw]
}

重さが非常に大きい場合のナップサック問題

すべての重さに対して最大価値を求める方法は不可能になるため、総価値に対して、それを与える最小重量を求める。

const MAX_N: usize = 100;
const MAX_W: usize = 1000_000_000;
const MAX_V: usize = 100;
const INF: usize = MAX_N * MAX_W + 1;
fn napsack3(tw: usize, items: Vec<(usize, usize)>) -> usize {
    // dp[i][v] = items[..i]を使い、総価値をちょうどvにするようにした時の最小の総重量
    // 解がなければINF(実際に取り得ないほど大きい値)
    // dp[items.len()][v] ≦ twとなるような最大のvを求めたい
    let mut dp = [[INF; MAX_N * MAX_V + 1]; MAX_N];
    let tv: usize = items.iter().map(|(_, v)| v).sum(); // 全部の品物が使えた場合の総価値
    dp[0][0] = 0; // 0個の時、価値をゼロにすることはできるが、他は無理(INFのまま)
    for i in 0..items.len() {
        for v in 0..=tv {
            dp[i + 1][v] = if v < items[i].1 {
                dp[i][v]
            } else {
                std::cmp::min(dp[i][v], items[i].0 + dp[i][v - items[i].1])
            };
        }
    }
    dp[items.len()].iter().rposition(|&w| w <= tw).unwrap()
}

• 計算時間はO(n×tv)
• 最大のvを求めるにはrpositionメソッドが便利

個数制限付き部分和問題

fn subsum(k: usize, a: &[usize], m: &[i32]) -> bool {
    const MAX_K: usize = 100_000;
    // ループがi回回った時点で、dp[x] = a[..i]を使ってxを作ったときに、余っているa[i-1]の個数
    // xを作れないときは負の数。0以上なら作れる。
    // ループをa.len()回回して、dp[k] >= 0を返せば良い。
    let mut dp = [-1; MAX_K + 1];
    // a[..0] = []でも、0ならば作れる。
    dp[0] = 0;
    for i in 0..a.len() {
        for x in 0..=k {
            dp[x] = if dp[x] > -1 {
                m[i] // a[..i]で作れるのでa[i]はm[i]個余る
            } else if x < a[i] {
                -1 // a[..i]で作れず、a[i]を使っても作れない
            } else {
                dp[x - a[i]] - 1 // a[i]を1つ使う
            }
        }
    }
    dp[k] > -1
}

メモ

• 2次元配列にするとメモリを使い切ったらしくstack overflowになってしまうので、配列更新法にした。

最長増大部分列問題

Longest increasing subsequence

fn lis(a: &[i32]) -> i32 {
    // dp[i] = 最後がa[i]であるような最長増大部分列の長さ
    let mut dp = vec![0; a.len()];
    for i in 0..a.len() {
        dp[i] = a[..i]
            .iter()
            .zip(&dp)
            .filter_map(|(&x, &y)| (x < a[i]).then(|| y))
            .max()
            .unwrap_or(0)
            + 1;
    }
    *dp.iter().max().unwrap()
}

fn lis2(a: &[i32]) -> usize {
    const INF: i32 = 1000001;
    // dp[i] = k回回った時点で、a[..k]を使って長さi+1の増大部分列を作ったときの、もっとも小さい最終要素
    // 作れないときはINF
    let mut dp = vec![INF; a.len()];
    for &x in a {
        let i = dp.partition_point(|&y| y < x);
        dp[i] = x;
    }
    dp.partition_point(|&x| x < INF)
}

メモ

• Rust 1.50からboolをOptionに変換するbool::thenを使えるようになった。
• partition_pointメソッドでCのupper_bound, lower_boundと同じことができる。lis2の計算時間はO(nlogn)となる。

分割数

整数nをm以下の正整数の和に分割する(n, m ≧ 0)。m以下の正整数の和、というのは、m個の非負整数の和といっても同じである。

/// 整数nの分割
/// * `n` - 分割したい整数
/// * `m` - m個に分割する(0も使用可能)  
/// * `mm` - これで割った余を求める
fn partition(n: usize, m: usize, mm: usize) -> usize {
    const M_MAX: usize = 1000;
    // i回回った時点で、dp[j] = jをi個以下に分割する方法の総数
    let mut dp = [0; M_MAX + 1];
    dp[0] = 1; // 0は0個に分割できる(長さゼロの数列の和になる)
    for i in 1..=m {
        for j in i..=n {
            dp[j] = (dp[j] + dp[j - i]) % mm;
        }
    }
    dp[n]
}

重複組合せ

それぞれa[0], a[1], ..., a[n-1]個ずつあるn個の品物からm個選び出す組み合わせの総数

fn multi_choose(a: &[usize], m: usize, mm: usize) -> usize {
    const M_MAX: usize = 1000;
    const N_MAX: usize = 1000;
    // dp[i][j] = a[..i]を使ってmを作る方法の総数
    let mut dp = [[0; M_MAX + 1]; N_MAX + 1];
    // a[..0] = []を使って0を作る組み合わせは1通り
    dp[0][0] = 1;
    // dp[i+1][j] = dp[i][j-a[i]] + ... + dp[i][j] の実装
    // dp[i+1][j-1]から足したり引いたりして計算する
    for i in 0..a.len() {
        dp[i + 1][0] = dp[i][0];
        for j in 1..=m {
            dp[i + 1][j] = dp[i][j] + dp[i + 1][j - 1];
            if a[i] + 1 <= j {
                dp[i + 1][j] -= dp[i][j - a[i] - 1];
            };
            dp[j + 1][j] %= mm;
        }
    }
    dp[a.len()][m]
}

次章はデータ構造を扱う。