バイトニック経路版の平面巡回セールスマン問題

巡回セールスマン問題(TSP)は有名なNP完全問題だが、ユークリッド平面におけるTSPは経路がbitonic tourであるという制約をつければ動的計画法で簡単に解けるようになる。以下そのメモ。
巡回セールスマン問題 - Wikipedia

bitonic tourとは

bitonic tourとは、左端の点から出発し、必ず左から右に向かって右端の点まで行き、つぎに右から左の方向に出発点まで戻ってくるような経路。wikipediaの図がわかりやすい。
Bitonic tour - Wikipedia, the free encyclopedia

解法

巡回するユークリッド座標上のN個の点を、x座標の昇順にソートした点列p1..pNを考える。1<=i2 , i=j-1) \end{array} \right]

3つ目の式は図で表現するとこんな感じ。下図では、3つのうち経路最短となる真ん中が部分bitonic tourの最適解となる。

で、これはL[N][N]をボトムアップに埋めていく動的計画法に簡単に落とし込める。計算量はO(n^2)。

#include <iostream>
#include <cmath>
#include <limits>
#include <vector>
#include <utility>

using namespace std;
typedef pair<int, int> Point;

double dist(const Point& p1, const Point& p2) {
    return sqrt(pow(p1.first - p2.first, 2.0) + pow(p1.second - p2.second, 2.0));
}

double bitonic_tsp_distance(const vector<Point>& p) {
    const int N = p.size();
    vector<vector<double> > L(N, vector<double>(N, 0.0));

    for (int j = 1; j < N; j++) {
        for (int i = 0; i < j; i++) {
            if (i == 0 && j == 1) {
                L[i][j] = dist(p[i], p[j]);
            } else if (i < j - 1) {
                L[i][j] = L[i][j-1] + dist(p[j-1], p[j]);
            } else {
                L[i][j] = numeric_limits<double>::infinity();
                for (int k = 0; k < i; k++)
                    L[i][j] = min(L[i][j], L[k][i] + dist(p[k], p[j]));
            }
        }
    }

    double ans = numeric_limits<double>::infinity();

    for (int k = 0; k < N-1; k++)
        ans = min(ans, L[k][N-1] + dist(p[k], p[N-1]));

    return ans;
}

int main(void) {
    vector<pair<int,int> > points;
 
    points.push_back(make_pair(0, 0));
    points.push_back(make_pair(1, 6));
    points.push_back(make_pair(2, 3));
    points.push_back(make_pair(5, 2));
    points.push_back(make_pair(6, 5));
    points.push_back(make_pair(7, 1));
    points.push_back(make_pair(8, 4));

    cout << "answer:" << bitonic_tsp_distance(points) << endl;
}

どういう時に使えそうか

普通に考えると、x座標が重複するような点群に対してはうまくいかない。また、(bitonic tourに限定されない)最適解に対して、精度の保証も取れない。普通のケースなら、O(n^2)で精度も最適解の2倍以下になる最近追加法でいい。DPが綺麗に適用できて面白いけど、実際の使いどころは無いのかなーという残念な結論。