Posted at

動的計画法をデータ構造で高速化するテクニック

More than 1 year has passed since last update.

本記事は U-TOKYO AP Advent Calendar 2017 の 7 日目の記事です。

こんにちは!計数工学科に所属する @kagamiz と言うものです。これまで Qiita には色んな場面でお世話になってきたので、自分が投稿する機会にめぐりあえて少し感動しています(?)。

今回は動的計画法をデータ構造で高速化する幾つかのテクニックについて、問題を解きながら紹介、そして自分が作った問題を宣伝していきたいと思います。


そもそも動的計画法とは?

動的計画法の説明として分かりやすい記事は大量にあります。幾つか挙げると

などの記事があります。基礎的な考え方はこれらの記事を参照してください。


データ構造を用いた高速化

動的計画法の遷移式に於いて、区間に対して操作をする必要がある場合、


  • 区間の最大値や最小値を取得する遷移式になった場合 → セグメント木

  • 固定長の区間の最大値や最小値が欲しい場合 → 両端キュー

を使うことで、愚直な DP の遷移式を高速化できることがあります。

セグメント木についての解説はぼくの過去のブログに載せています。


問題例①:Minimizing Maximizer

(問題引用元 : 1769 -- Minimizing Maximizer)


問題概要

$m$ 個の閉区間 $[a_i, b_i]$ が与えられる。 $(1 \leq a_i < b_i \leq n, a_i, b_i \in \mathbb{Z})$

このとき、次の条件をみたすように閉区間 $[1, n]$ を覆いたい。


  • 番号順に区間を使うか使わないかを決める。

  • 区間 $[a_i, b_i]$ を使うとき、区間 $[1, a_i]$ は既に覆われているようにしたい。ただし、区間が 1 つも選ばれていないときは $a_i = 1$ であれば使っても良い。

式で上記の条件を書くと

\begin{cases}

\displaystyle\bigcup_{k}[a_{i_k}, b_{i_k}] = [1, n] & \\
[1, b_{i_{j}}] \cap [a_{i_{j+1}}, b_{i_{j+1}}] \neq \emptyset
\end{cases}

となる。このときに使う区間の個数を最小化せよ。


制約


  • $2 \leq n \leq 50,000$

  • $1 \leq m \leq 500,000$

  • $1 \leq a_i < b_i \leq n \ (1 \leq i \leq m)$


解き方

先頭から区間を使っていくことを考えると、

$\mathrm{dp}[x] :=$ $[1, x]$ を覆うために必要な区間の個数の最小値

という dp テーブルを定義し、

$\mathrm{dp}[1] = 0$ ($[1, 1]$ は無条件で覆えているとみなす)

$\mathrm{dp}[k] = \infty$ ($k \geq 2$)

と初期化し、

$\mathrm{dp}[b_i] = \min \{ \mathrm{dp}[b_i], \displaystyle\min_{a_i \leq k < b_i} \mathrm{dp}[k] + 1 \} $

と dp テーブルを $i$ の昇順に更新すれば良いことがわかります。しかし、愚直に最小値の検索をすると時間計算量は $\mathrm{O}(nm)$ 時間となってしまいます。


minimizing_maximizer_slow.cpp

#include <cstdio>

#include <algorithm>

#define MAX_N (50000)

//dp[x] := [1, x] を覆うのに必要な区間数の最小値
int dp[MAX_N + 1];

int main()
{
int n, m;

scanf("%d %d", &n, &m);

std::fill(dp, dp + n + 1, MAX_N);

dp[1] = 0;

for (int i = 0; i < m; i++){
int a, b;
scanf("%d %d", &a, &b);

int tmp = dp[b];
for (int j = a; j < b; j++){
tmp = std::min(tmp, dp[j] + 1);
}

dp[b] = std::min(dp[b], tmp);
}

printf("%d\n", dp[n]);

return 0;
}


この dp のボトルネックは、最小値を求める操作に $\mathrm{O}(n)$ 時間かかってしまっている所です。なので、最小値を RMQ によって求めれば、最小値を求める操作が $\mathrm{O}(\log n)$ となり、この問題が $\mathrm{O}(m \log n)$ 時間で解けることがわかります。


minimizing_maximizer_fast.cpp

#include <cstdio>

#include <algorithm>

#define MAX_N (50000)
#define SEG_SIZE (65536)

int seg[2 * SEG_SIZE - 1];

void update(int k, int x)
{
k += SEG_SIZE - 1;
seg[k] = x;
while (k){
k = (k - 1) / 2;
seg[k] = std::min(seg[k * 2 + 1], seg[k * 2 + 2]);
}
}

int getMin(int a, int b, int k = 0, int l = 0, int r = SEG_SIZE)
{
if (b <= l || r <= a){
return MAX_N;
}
if (a <= l && r <= b){
return seg[k];
}
int lf = getMin(a, b, k * 2 + 1, l, (l + r) / 2),
rg = getMin(a, b, k * 2 + 2, (l + r) / 2, r);

return std::min(lf, rg);
}

//dp[x] := [1, x] を覆うのに必要な区間数の最小値
int dp[MAX_N + 1];

int main()
{
int n, m;

scanf("%d %d", &n, &m);

std::fill(dp, dp + n + 1, MAX_N);
std::fill(seg, seg + 2 * SEG_SIZE, MAX_N);

dp[1] = 0;
update(1, 0);

for (int i = 0; i < m; i++){
int a, b;
scanf("%d %d", &a, &b);

int tmp = std::min(dp[b], getMin(a, b) + 1);

dp[b] = tmp;
update(b, tmp);
}

printf("%d\n", dp[n]);

return 0;
}



問題例② : Watching Fireworks is Fun

(問題引用元 : Problem - 372C - Codeforces)

自作の問題なので引っ張ってきました。


問題概要

$m$ 回花火の打ち上げが数直線 $[1, n]$ 上で行われる。$i$ 番目の花火は時刻 $t_i$ に場所 $a_i$ で打ち上げられる。

$i$ 番目の花火を場所 $x$ で観測すると、幸福度 $b_i - |a_i - x|$ が得られる。

$1$ 秒間に $d$ だけ左右に動け、時刻 $1$ に任意の $[1, n]$ 内の整数点にいても良いときに、得られる幸福度の和を最大化せよ。ただし整数点にのみしか移動できないとする。


制約


  • $1 \leq m \leq 300$

  • $1 \leq n \leq 100,000$

  • $1 \leq d \leq n$

  • $1 \leq a_i \leq n$

  • $1 \leq b_i \leq 1,000,000,000$

  • $1 \leq t_i \leq 1,000,000,000$

  • 入力はすべて整数


解き方


愚直な解き方

$\mathrm{dp}[i][j] :=$ $i$ 番目の花火を場所 $j$ で見たときの幸福度の和の最大値、とすると、DP の遷移式は

\begin{cases}

\mathrm{dp}[1][x] &= b_1 - |a_1 - x|\\
\mathrm{dp}[k][x] &= \displaystyle\max_{x - (t_k - t_{k-1})d \leq i \leq x + (t_k - t_{k-1})d} b_k - |a_k - x| + \mathrm{dp}[k - 1][i] & (k > 1)
\end{cases}

となります。この遷移式をそのまま実装すると $\mathrm{O}(mn^2)$ 時間でこの問題を解くことが出来ます。


セグメント木を用いた高速化

先程の DP の遷移式のうち、$k > 1$ の部分は

$$\mathrm{dp}[k][x] = b_k - |a_k - x| + \displaystyle\max_{x - (t_k - t_{k-1})d \leq i \leq x + (t_k - t_{k-1})d} \mathrm{dp}[k - 1][i]$$

という形で、 max の中の部分を $k$ に依存しない形で書くことが出来ます。そこで、「Minimizing Maximizer」でやったように、max を求める処理を RMQ (Range Maximum Query のほう) を用いて行えば、この問題は $\mathrm{O}(mn \log n)$ 時間で解くことが出来ます。


両端キューを用いた高速化

更に遷移式に注目すると、同じ $k$ について、max を取る区間の長さが $2d(t_k - t_{k-1}) + 1$ で固定されているという特徴があります。この特徴に注目すると、 Sliding Window Maximium という方法を用いることで、ひとつの $k$ について $\mathrm{dp}[k][\cdot]$ を $\mathrm{O}(n)$ ですべて計算することが出来ます。 $k$ は $1$ から $m$ まで動くため、全体の時間計算量は $\mathrm{O}(mn)$ となります。


  • Sliding Window Maximum の説明

長さ $n$ の数列 $\{a_i\}$ が与えられたとき、長さ $n - k + 1$ の数列 $\{b_i\}$ を

$$b_i := \displaystyle\max_{i \leq j \leq i + k - 1} a_j$$

で定義します。数列 $\{b_i\}$ は、両端キューを用いると $\mathrm{O}(n)$ 時間で求めることが出来ます。

両端キューには数列 $\{a_i\}$ の添字を入れて管理します。

以下のように、現在考えている区間の最大値のうち、最も後ろに現れるものがキューの先頭に来るようにキューの状態を更新します。

キューに要素を追加するときは、キューの末尾に要素を追加します。このとき、$a_{キューの末尾} \leq a_{追加する要素の位置}$ が満たされるあいだ、キューの末尾の要素を削除します。キューが空になるか、$a_{キューの末尾} > a_{追加する要素の位置}$ となったときに、追加する要素の添字をキューの末尾に追加します。

この条件でキューに要素を追加していくことで、両端キューの中身は、先頭から見ていくと常に単調に減少しているようになります。そこで、$b_j$ を両端キューの先頭の値としてやれば良いです。

そして、キューの先頭の要素の添字が今見ている区間を外れている間、キューの先頭を削除する、ということを繰り返していけば良いです。

この操作において、数列の 1 つの要素はたかだか 1 回キューに追加され、たかだか 1 回キューから削除されます。よって、Sliding Window Maximum を用いると $\{b_i\}$ の構成は $\mathrm{O}(n)$ 時間で行なえます。

C++ では、 deque コンテナを用いることで両端キューを扱うことが出来ます。


fireworks_deque

#include <cstdio>

#include <climits>
#include <deque>
#include <algorithm>

using namespace std;

typedef long long lint;

int main()
{
int N, M, D;
int a[1000], b[1000], t[1000];
lint dp[2][300001];

scanf("%d %d %d", &N, &M, &D);

for (int i = 0; i < M; i++){
scanf("%d %d %d", a + i, b + i, t + i);
}

for (int i = 1; i <= N; i++){
dp[0][i] = 0;
}

int p = 1;
int pT = 1;
for (int i = 0; i < M; i++){
deque<int> dq;
lint interval = min((lint)N, (lint)(t[i] - pT) * D);

for (lint j = 1; j <= interval; j++){
while (dq.size() && dp[1 - p][dq[dq.size() - 1]] <= dp[1 - p][(int)j]) dq.pop_back();
dq.push_back((int)j);
}

for (lint j = 1; j <= (lint)N; j++){
lint k = j + interval;
if (k <= (lint)N){
while (dq.size() && dp[1 - p][dq[dq.size() - 1]] <= dp[1 - p][(int)k]) dq.pop_back();
dq.push_back((int)k);
}
dp[p][j] = dp[1 - p][dq[0]] + (lint)(b[i] - abs(a[i] - j));
while (dq.size() && dq[0] <= (int)(j - interval)) dq.pop_front();
}
pT = t[i];
p = 1 - p;
}

lint ans = LLONG_MIN;
for (int i = 1; i <= N; i++)
ans = max(ans, dp[1 - p][i]);

printf("%lld\n", ans);

return (0);
}



  • 補足

Sliding Window Maximum は、着目している区間の端点が単調非減少であれば用いることができるので、尺取り法との相性も良いです。


注意すべきこと

動的計画法では、いかにうまく状態を定めるかが肝になります。今回紹介したテクニックは当然必ず使えるわけではないので、あくまでアイディアの 1 つとして覚えていて貰えれば幸いです。


問題集

この節では、データ構造で DP を高速化する問題の紹介を行います。今後記事中で高速化法について詳しく述べるかもしれません。



  1. 個数制限つき Knapsack 問題


    • 問題概要


      • $n$ 種類の品物があります。 $i$ 種類目の品物は $m_i$ 個あり、価値は $v_i$ で、重さは $w_i$ です。
        重さの和が $W$ 以下である荷物の集合の価値の和を最大化してください。



    • 制約


      • $1 \leq n \leq 100$

      • $1 \leq m_i \leq 10,000$

      • $1 \leq w_i, v_i \leq 100$

      • $1 \leq W \leq 10,000$



    • 愚直な DP をすると時間計算量が $\mathrm{O}(nmW)$ となってしまいます。データ構造を用いて時間計算量を $\mathrm{O}(nW)$ に高速化することが出来ます。もしくは、荷物の参照をそれぞれの種類につき$\mathrm{O}(\log m)$回に減らし、全体の計算量を $\mathrm{O}(nW \log m)$時間にすることもできます。




  2. House Moving




  3. (おまけ) Easy Dynamic Programming


    • E : はじめての動的計画法


    • DP の高速化とはまったく関係ない。DP にちなんで自分の作った DP に関する問題の逆問題を紹介します :)




最後に

データ構造を用いる方法以外にも、動的計画法を高速化するテクニックは色々存在します(Monge 性など)。今後機会があれば、それらの方法について投稿していきたいと思います。


参考文献


  • プログラミングコンテストチャレンジブック (Minimizing Maximier, 個数制限つき Knapsack 問題)