bit DP:TSP(ABC180E) - わりとよくある備忘録

Traveling Salesman Problem

巡回セールスマン問題（じゅんかいセールスマンもんだい、英: traveling salesman problem、TSP）は、都市の集合と各2都市間の移動コスト（たとえば距離）が与えられたとき、全ての都市をちょうど一度ずつ巡り出発地に戻る巡回路のうちで総移動コストが最小のものを求める（セールスマンが所定の複数の都市を1回だけ巡回する場合の最短経路を求める）組合せ最適化問題である。
by 巡回セールスマン問題 - Wikipedia

全探索を行えば勿論最適解は求まるが、あり得る訪問経路の数は $(N-1)!$ 通りになるのですぐに組合せ爆発する。(例えば $N-1=16$ なら $16!≒2e+13$ )
wikipediaにあるようにTSPはいわゆるNP困難な問題で、多項式時間で解けるアルゴリズムは見つかっていないので現実的な時間で最適解を求められるかどうかは都市の数( $=N$ )に依存する。

AtCoderの時間制約的には $N$ が10以上程度になると既に全探索では厳しくなってくるが、 $N≦17$ 程度ならばDPによる高速化で十分に解ける場合が多い。

ちなみに全探索版で解ける問題はABC183Cで出てましたね。
atcoder.jp

どう解くか

というわけでDPによる高速化をしていく。

重要なポイントは、現状態から終端までの最適な訪問順は「これまでどの都市を訪問してきたか」ということと「今いる都市」という情報のみに依存し、「これまで訪問した都市の訪問順序」には依らないということ。つまり残りの未訪問の都市集合と、スタート地点が同じなら、最適な経路は当然一意に定まる。
※このエントリの図が非常に分かりやすく参考になりました。
dalekspritner.hatenablog.com

よってあり得る状態は、「これまでどの都市を訪問してきたか」という状態 $S$ と、「今いる都市」である $v$ で表現でき、 $dp[S][v]$ を状態 $S$ - $v$ に至るまでの最小コストとすると $dp[S_{all}][0]$ を求めれば良いということになる。
遷移は、 $S$ から今いる都市を除いた状態(現状態の前状態としてあり得る状態)を $S'$ と表すと、 $dp[S][v]$ は $dp[S'][v']+cost[v'][v]$ ,( $v' \in S'$ )の中から最小のものを選べば良いので下記のようになる。

$dp [S][v]=min(dp[S][v],dp[S'][v']+cost[v'][v])$

bitによる状態Sの表現

上記までで遷移は分かったので、あとは $S$ , $S'$ を配列のindexで取り扱えるような値で表せればok。ここで $S$ がなんだったかを振り返ると、ある都市を訪問したか否かを表現できれば良いのでbitで取り扱えば良さそう。

例えば、都市1,3,4を訪問した場合だと下記のように表せば良い。

$S＝\{1,3,4\} => 0b11010$
これで各状態を整数値で取り扱えるようになったので、前述の遷移を実装すれば終了。

計算量

あり得る状態数は $2^n$ であり、各状態に対して遷移 $n^2$ を集約するので計算量は $O(2^n*n^2)$ になる。

実装例

ここまでのロジックを忠実に実装すると下記のようになるが、実際最内ループの条件分岐は無くても動く。（初期値をINFにしていることでどちらにしても更新されないので。）

#include <algorithm>
#include <cmath>
#include <cstring>
#include <deque>
#include <iomanip>
#include <iostream>
#include <map>
#include <queue>
#include <set>
#include <stack>
#include <string>
#include <vector>

#define rep(i, n) for (int i = 0; i < (n); i++)
#define reps(i, n, s) for (int i = (s); i < (n); i++)

using ll = long long;
using namespace std;
constexpr long long INF = 1LL << 60;
constexpr int MAX_N = 17;

ll dp[(1 << MAX_N)][MAX_N];

int main() {
    cin.tie(0);
    ios::sync_with_stdio(false);

    int n;
    cin >> n;

    vector<ll> x(n), y(n), z(n);
    rep(i, n) cin >> x[i] >> y[i] >> z[i];
    // calc cost
    vector<vector<ll>> cost(n, vector<ll>(n, 0));
    rep(i, n) {
        reps(j, n, i + 1) {
            ll c = abs(x[i] - x[j]) + abs(y[i] - y[j]);
            cost[i][j] = c + max(0LL, z[j] - z[i]);
            cost[j][i] = c + max(0LL, z[i] - z[j]);
        }
    }

    // TSP
    fill(dp[0], dp[(1 << MAX_N)], INF);
    dp[0][0] = 0;
    reps(S, 1 << n, 1) {
        rep(v, n) {
            if ((S >> v) & 1) {
                rep(j, n) {
                    ll S_hat = S - (1 << v);
                    if (S_hat >> j & 1 || S_hat == 0) { dp[S][v] = min(dp[S][v], dp[S_hat][j] + cost[j][v]); }
                }
            }
        }
    }

    cout << dp[(1 << n) - 1][0] << endl;
    return 0;
}

bit DP楽しい。