けんちょんの競プロ精進記録

競プロの精進記録や小ネタを書いていきます

競プロ典型 90 問 077 - Planes on a 2D Plane(★7)

二部マッチングの練習問題

問題概要

二次元平面上に  N 体の飛行機がある。飛行機  i は座標  (A_{x, i}, A_{y, i}) にいる。各飛行機に対して適切に 8 種類の向き付けをしたい (上・下・左・右・左上・左下・右上・右下)。

具体的には、 T 秒後において  N 体の飛行機が  (B_{x, 0}, B_{y, 0}), (B_{x, 1}, B_{y, 1}), \dots, (B_{x, N-1}, B_{y, N-1}) にいる状態にしたい。

そのようなことが可能となる向き付けを構成せよ。不可能である場合はその旨を報告せよ。

制約

  •  1 \le N \le 20000
  •  T, 座標値 \le 10^{9}

解法:二部マッチング

この問題は、 N 個の座標  (A_{x, 0}, A_{y, 0}), (A_{x, 1}, A_{y, 1}), \dots, (A_{x, N-1}, A_{y, N-1}) と、 N 個の座標  (B_{x, 0}, B_{y, 0}), (B_{x, 1}, B_{y, 1}), \dots, (B_{x, N-1}, B_{y, N-1}) との間のマッチング問題だと言えます。

たとえば、 T = 1 A = (0, 0), (0, 2), (2, 0) B = (1, 1), (0, 1), (1, 0) のとき、下図のような二部グラフにおいて、3 個のペアが作れるかを問う問題となります。この場合は、

  • (0, 0) を上に 1 秒だけ動かして (0, 1)
  • (0, 2) を右下に 1 秒だけ動かして (1, 1)
  • (2, 0) を左に 1 秒だけ動かして (1, 0)

という解があります (他にも 3 ペア取れる解があります)。

このように、2 つのカテゴリ間でできるだけ多くのペアを作る問題を最大二部マッチング問題と言います。今回の問題の場合、 N 個の座標と  N 個の座標の間で最大マッチングを求めてあげて、その値が  N であれば "Yes"、 N 未満であれば "No" ということになる。

最大二部マッチング問題の解き方

最大二部マッチング問題は最大フロー問題に帰着させて解くことができる。その詳細については、次の記事に詳しく書きました。

qiita.com

なお、最大フロー問題を解くアルゴリズムとして Dinic 法を採用したとき、最大二部マッチング問題を解くときの計算量は、グラフの頂点数を  V、辺数を  E として  O(E\sqrt{V}) となることが知られています (詳しいことは Hopcroft-Karp 法の計算量解析より)。

今回は頂点数は  O(N) であり、辺数も  O(N) に限られる (各初期座標につき、多くとも 8 個までしか対応できる頂点が存在しないため) ことに注意しましょう。よって計算量は  O(N \sqrt{N}) であって十分間に合います。

 

コード

ACL を用いた実装

まずは ACL の最大流アルゴリズムを用いて通します。使い方はドキュメントが参考になります。

マッチングに使われたペアを復元する際には「各辺を流れているフローの値が 1 になっているかどうか」を判断しています。

#include <bits/stdc++.h>
#include <atcoder/maxflow>
using namespace std;
using namespace atcoder;

int main() {
    // 8 方向への移動
    const vector<int> dx = {1, 1, 0, -1, -1, -1, 0, 1};
    const vector<int> dy = {0, 1, 1, 1, 0, -1, -1, -1};
    
    // 入力
    int N, T;
    cin >> N >> T;
    vector<int> AX(N), AY(N), BX(N), BY(N);
    for (int i = 0; i < N; ++i) cin >> AX[i] >> AY[i];
    for (int i = 0; i < N; ++i) cin >> BX[i] >> BY[i];
    map<pair<int,int>, int> B;
    for (int i = 0; i < N; ++i) B[{BX[i], BY[i]}] = i;
    
    // フローネットワークを作る
    // 左側頂点:0, 1, ..., N-1、右側頂点:N, N+2, ..., 2N-1 とする
    mf_graph<int> G(N*2+2);  // ソースとシンクも込みで頂点は 2N+2 個
    int s = N*2, t = N*2 + 1;  // ソースとシンク
    
    // ソースと左側頂点、右側頂点とシンクを繋ぐ
    for (int i = 0; i < N; ++i) {
        G.add_edge(s, i, 1);
        G.add_edge(i+N, t, 1);
    }
    
    // 左側頂点と右側頂点のうち、時刻 T で移り合う部分を繋ぐ
    for (int i = 0; i < N; ++i) {
        for (int dir = 0; dir < 8; ++dir) {
            int bx = AX[i] + dx[dir] * T, by = AY[i] + dy[dir] * T;
            if (B.count({bx, by})) {
                G.add_edge(i, B[{bx, by}] + N, 1);
            }
        }
    }
    
    // 最大流を流す
    int max_flow = G.flow(s, t);
        
    // 復元する
    if (max_flow < N) cout << "No" << endl;
    else {
        cout << "Yes" << endl;
        vector<int> res(N);
        const auto &edges = G.edges();
        for (const auto &e : edges) {
            if (e.from == s || e.to == t || e.flow == 0) continue;
            
            // e.from から e.to への移動ベクトルを特定する
            int dir = 0;
            for (; dir < 8; ++dir) {
                if (BX[e.to-N] == AX[e.from] + dx[dir] * T &&
                    BY[e.to-N] == AY[e.from] + dy[dir] * T) break;
            }
            res[e.from] = dir+1;
        }
        for (auto val : res) cout << val << " ";
        cout << endl;
    }
}

自前の Dinic 法でも AC

自前実装の Dinic 法でも通します。ACL の最大流ライブラリとインターフェースを概ね一緒にしています。

#include <bits/stdc++.h>
using namespace std;

// edge class (for network-flow)
template<class FLOWTYPE> struct FlowEdge {
    // core members
    int rev, from, to;
    FLOWTYPE cap, icap, flow;
    
    // constructor
    FlowEdge(int r, int f, int t, FLOWTYPE c)
    : rev(r), from(f), to(t), cap(c), icap(c), flow(0) {}
    void reset() { cap = icap, flow = 0; }
    
    // debug
    friend ostream& operator << (ostream& s, const FlowEdge& E) {
        return s << E.from << "->" << E.to << '(' << E.flow << '/' << E.icap << ')';
    }
};

// graph class (for network-flow)
template<class FLOWTYPE> struct FlowGraph {
    // core members
    vector<vector<FlowEdge<FLOWTYPE>>> list;
    vector<pair<int,int>> pos;  // pos[i] := {vertex, order of list[vertex]} of i-th edge
    
    // constructor
    FlowGraph(int n = 0) : list(n) { }
    void init(int n = 0) {
        list.assign(n, FlowEdge<FLOWTYPE>());
        pos.clear();
    }
    
    // getter
    vector<FlowEdge<FLOWTYPE>> &operator [] (int i) {
        return list[i];
    }
    const vector<FlowEdge<FLOWTYPE>> &operator [] (int i) const {
        return list[i];
    }
    size_t size() const {
        return list.size();
    }
    FlowEdge<FLOWTYPE> &get_rev_edge(const FlowEdge<FLOWTYPE> &e) {
        if (e.from != e.to) return list[e.to][e.rev];
        else return list[e.to][e.rev + 1];
    }
    FlowEdge<FLOWTYPE> &get_edge(int i) {
        return list[pos[i].first][pos[i].second];
    }
    const FlowEdge<FLOWTYPE> &get_edge(int i) const {
        return list[pos[i].first][pos[i].second];
    }
    vector<FlowEdge<FLOWTYPE>> get_edges() const {
        vector<FlowEdge<FLOWTYPE>> edges;
        for (int i = 0; i < (int)pos.size(); ++i) {
            edges.push_back(get_edge(i));
        }
        return edges;
    }
    
    // change edges
    void reset() {
        for (int i = 0; i < (int)list.size(); ++i) {
            for (FlowEdge<FLOWTYPE> &e : list[i]) e.reset();
        }
    }
    void change_edge(FlowEdge<FLOWTYPE> &e, FLOWTYPE new_cap, FLOWTYPE new_flow) {
        FlowEdge<FLOWTYPE> &re = get_rev_edge(e);
        e.cap = new_cap - new_flow, e.icap = new_cap, e.flow = new_flow;
        re.cap = new_flow;
    }
    
    // add_edge
    void add_edge(int from, int to, FLOWTYPE cap) {
        pos.emplace_back(from, (int)list[from].size());
        list[from].push_back(FlowEdge<FLOWTYPE>((int)list[to].size(), from, to, cap));
        list[to].push_back(FlowEdge<FLOWTYPE>((int)list[from].size() - 1, to, from, 0));
    }

    // debug
    friend ostream& operator << (ostream& s, const FlowGraph &G) {
        const auto &edges = G.get_edges();
        for (const auto &e : edges) s << e << endl;
        return s;
    }
};

// Dinic
template<class FLOWTYPE> FLOWTYPE Dinic
 (FlowGraph<FLOWTYPE> &G, int s, int t, FLOWTYPE limit_flow)
{
    FLOWTYPE current_flow = 0;
    vector<int> level((int)G.size(), -1), iter((int)G.size(), 0);
    
    // Dinic BFS
    auto bfs = [&]() -> void {
        level.assign((int)G.size(), -1);
        level[s] = 0;
        queue<int> que;
        que.push(s);
        while (!que.empty()) {
            int v = que.front();
            que.pop();
            for (const FlowEdge<FLOWTYPE> &e : G[v]) {
                if (level[e.to] < 0 && e.cap > 0) {
                    level[e.to] = level[v] + 1;
                    if (e.to == t) return;
                    que.push(e.to);
                }
            }
        }
    };
    
    // Dinic DFS
    auto dfs = [&](auto self, int v, FLOWTYPE up_flow) {
        if (v == t) return up_flow;
        FLOWTYPE res_flow = 0;
        for (int &i = iter[v]; i < (int)G[v].size(); ++i) {
            FlowEdge<FLOWTYPE> &e = G[v][i], &re = G.get_rev_edge(e);
            if (level[v] >= level[e.to] || e.cap == 0) continue;
            FLOWTYPE flow = self(self, e.to, min(up_flow - res_flow, e.cap));
            if (flow <= 0) continue;
            res_flow += flow;
            e.cap -= flow, e.flow += flow;
            re.cap += flow, re.flow -= flow;
            if (res_flow == up_flow) break;
        }
        return res_flow;
    };
    
    // flow
    while (current_flow < limit_flow) {
        bfs();
        if (level[t] < 0) break;
        iter.assign((int)iter.size(), 0);
        while (current_flow < limit_flow) {
            FLOWTYPE flow = dfs(dfs, s, limit_flow - current_flow);
            if (!flow) break;
            current_flow += flow;
        }
    }
    return current_flow;
};

template<class FLOWTYPE> FLOWTYPE Dinic(FlowGraph<FLOWTYPE> &G, int s, int t) {
    return Dinic(G, s, t, numeric_limits<FLOWTYPE>::max());
}


int main() {
    // 8 方向への移動
    const vector<int> dx = {1, 1, 0, -1, -1, -1, 0, 1};
    const vector<int> dy = {0, 1, 1, 1, 0, -1, -1, -1};
    
    // 入力
    int N, T;
    cin >> N >> T;
    vector<int> AX(N), AY(N), BX(N), BY(N);
    for (int i = 0; i < N; ++i) cin >> AX[i] >> AY[i];
    for (int i = 0; i < N; ++i) cin >> BX[i] >> BY[i];
    map<pair<int,int>, int> B;
    for (int i = 0; i < N; ++i) B[{BX[i], BY[i]}] = i;
    
    // フローネットワークを作る
    // 左側頂点:0, 1, ..., N-1、右側頂点:N, N+2, ..., 2N-1 とする
    FlowGraph<int> G(N*2+2);  // ソースとシンクも込みで頂点は 2N+2 個
    int s = N*2, t = N*2 + 1;  // ソースとシンク
    
    // ソースと左側頂点、右側頂点とシンクを繋ぐ
    for (int i = 0; i < N; ++i) {
        G.add_edge(s, i, 1);
        G.add_edge(i+N, t, 1);
    }
    
    // 左側頂点と右側頂点のうち、時刻 T で移り合う部分を繋ぐ
    for (int i = 0; i < N; ++i) {
        for (int dir = 0; dir < 8; ++dir) {
            int bx = AX[i] + dx[dir] * T, by = AY[i] + dy[dir] * T;
            if (B.count({bx, by})) {
                G.add_edge(i, B[{bx, by}] + N, 1);
            }
        }
    }
    
    // 最大流を流す
    int max_flow = Dinic(G, s, t);
        
    // 復元する
    if (max_flow < N) cout << "No" << endl;
    else {
        cout << "Yes" << endl;
        vector<int> res(N);
        const auto &edges = G.get_edges();
        for (const auto &e : edges) {
            if (e.from == s || e.to == t || e.flow == 0) continue;
            
            // e.from から e.to への移動ベクトルを特定する
            int dir = 0;
            for (; dir < 8; ++dir) {
                if (BX[e.to-N] == AX[e.from] + dx[dir] * T &&
                    BY[e.to-N] == AY[e.from] + dy[dir] * T) break;
            }
            res[e.from] = dir+1;
        }
        for (auto val : res) cout << val << " ";
        cout << endl;
    }
}