#1149469 (C++17) No.5023 Airlines Optimization

• 問題: No.5023 Airlines Optimization
• ユーザー: e869120
• 提出日時: 2026-02-26 01:23:23
• 言語: C++17 (gcc 15.2.0 + boost 1.89.0)
• 結果: AC
• 実行時間: 944 ms / 1,000 ms
• コード長: 12,084 bytes
• コンパイル時間: 9,373 ms
• コンパイル使用メモリ: 279,344 KB
• 実行使用メモリ: 8,092 KB
• スコア: 64,210,325
• 最終ジャッジ日時: 2026-02-26 01:25:07
• 合計ジャッジ時間: 103,071 ms
• ジャッジサーバーID （参考情報）: judge5 / judge4

ファイルパターン	結果
other	AC * 100

ソースコード

// Inspired by tetra4's solution
//   - https://yukicoder.me/submissions/1149387

#pragma GCC optimize("O3,unroll-loops")
#include <bits/stdc++.h>
using namespace std;

constexpr int N = 47, R = 1000, M = 400, K = 25;
constexpr int INF = 1e9;
constexpr int MAX_TIME_SLOTS = 182;
constexpr int MAX_FLIGHTS = 2500;

struct Flight {
    int u, s, v, t;
};

// スロット(0~180)を "HH:MM" 形式の文字列に変換
string to_time_str(int t) {
    int h = 6 + t / 12;
    int m = (t % 12) * 5;
    char buf[6];
    snprintf(buf, sizeof(buf), "%02d:%02d", h, m);
    return string(buf);
}

struct Input {
    array<int, N> x, y, w;
    array<int, M> a, s, b, t;
    bool valid_dist[N][N]{};
    long long weight[N][N]{};
    long long total_W = 0;

    Input() {
        int _;
        cin >> _ >> _;
        for (int i = 0; i < N; ++i) cin >> x[i] >> y[i] >> w[i];
        cin >> _;
        string s_str, t_str;
        for (int i = 0; i < M; ++i) {
            cin >> a[i] >> s_str >> b[i] >> t_str;
            --a[i], --b[i];
            s[i] = (stoi(s_str.substr(0, 2)) - 6) * 12 + stoi(s_str.substr(3, 2)) / 5;
            t[i] = (stoi(t_str.substr(0, 2)) - 6) * 12 + stoi(t_str.substr(3, 2)) / 5;
        }

        for (int i = 0; i < N; ++i) {
            for (int j = 0; j < N; ++j) {
                if (i == j) continue;
                long long dx = x[i] - x[j];
                long long dy = y[i] - y[j];
                if (dx * dx + dy * dy >= (R * R) / 16) {
                    valid_dist[i][j] = true;
                    weight[i][j] = (long long)w[i] * w[j];
                    total_W += weight[i][j] * 21;
                }
            }
        }
    }

    int get_duration(int u, int v) const {
        double dx = x[u] - x[v];
        double dy = y[u] - y[v];
        double d = sqrt(dx * dx + dy * dy);
        double minutes = 60.0 * d / 800.0 + 40.0;
        return (int)ceil(minutes / 5.0 - 1e-9);
    }
};

// 計算用のワークスペース (生配列で高速化)
struct ComputeSciWorkspace {
    int head[MAX_TIME_SLOTS + 1];
    int nxt[MAX_FLIGHTS];
    int f_u[MAX_FLIGHTS];
    int f_v[MAX_FLIGHTS];
    int f_t[MAX_FLIGHTS];
    int min_arr[MAX_TIME_SLOTS + 1][N];
};

struct MarginalGainWorkspace {
    int current_s_ci[N][N][21];
    long long edge_gain[N][N][MAX_TIME_SLOTS + 1];
    long long dp[N][MAX_TIME_SLOTS + 1];
    int nxt_v[N][MAX_TIME_SLOTS + 1];
    int dur[N][N];
};

void compute_s_ci(const Input& in, const vector<Flight>& flights,
                  int s_ci[N][N][21], ComputeSciWorkspace& ws) {
    for (int u = 0; u < N; ++u) {
        for (int dst = 0; dst < N; ++dst) {
            for (int t = 0; t < 21; ++t) {
                s_ci[u][dst][t] = -INF;
            }
        }
    }

    memset(ws.head, -1, sizeof(ws.head));
    int f_sz = flights.size();
    for (int i = 0; i < f_sz; ++i) {
        int s = flights[i].s;
        if (s >= 0 && s <= 180) {
            ws.nxt[i] = ws.head[s];
            ws.head[s] = i;
            ws.f_u[i] = flights[i].u;
            ws.f_v[i] = flights[i].v;
            ws.f_t[i] = flights[i].t;
        }
    }

    for (int dst = 0; dst < N; ++dst) {
        for (int t = 0; t <= 180; ++t) {
            for (int u = 0; u < N; ++u) ws.min_arr[t][u] = INF;
            ws.min_arr[t][dst] = t;
        }

        for (int t = 180; t >= 0; --t) {
            if (t + 1 <= 180) {
                memcpy(ws.min_arr[t], ws.min_arr[t + 1], sizeof(int) * N);
                ws.min_arr[t][dst] = t;
            }

            for (int e = ws.head[t]; e != -1; e = ws.nxt[e]) {
                int u = ws.f_u[e];
                int v = ws.f_v[e];
                int arr_time = ws.f_t[e];
                if (arr_time <= 180) {
                    if (ws.min_arr[arr_time][v] < ws.min_arr[t][u]) {
                        ws.min_arr[t][u] = ws.min_arr[arr_time][v];
                    }
                }
            }
        }

        for (int u = 0; u < N; ++u) {
            if (u == dst || !in.valid_dist[u][dst]) continue;
            int t = 180;
            for (int T_idx = 20; T_idx >= 0; --T_idx) {
                int T_target = 60 + T_idx * 6;
                t = min(t, T_target);
                while (t >= 0 && ws.min_arr[t][u] > T_target) {
                    t--;
                }
                if (t >= 0) {
                    s_ci[u][dst][T_idx] = t;
                } else {
                    break;
                }
            }
        }
    }
}

struct RivalSchedule {
    int s_sq[N][N][21];
    RivalSchedule(const Input& in, ComputeSciWorkspace& ws) {
        vector<Flight> r_flights;
        r_flights.reserve(M);
        for (int i = 0; i < M; ++i) r_flights.push_back({in.a[i], in.s[i], in.b[i], in.t[i]});
        compute_s_ci(in, r_flights, s_sq, ws);
    }
};

// ★ ビームサーチ用の状態構造体
struct State {
    vector<vector<Flight>> current_paths;
    vector<Flight> all_flights;
    long long score;

    bool operator<(const State& other) const {
        return score > other.score; // 降順ソート用
    }
};

// ★ ビームサーチ版のアルゴリズム
vector<vector<Flight>> solve_beam_search(const Input& in, const RivalSchedule& rival,
                                         ComputeSciWorkspace& c_ws,
                                         MarginalGainWorkspace& mg_ws) {
    const int BEAM_WIDTH = 3;
    const int BRANCHES = 15;

    for (int u = 0; u < N; ++u) {
        for (int v = 0; v < N; ++v) {
            mg_ws.dur[u][v] = (u != v) ? in.get_duration(u, v) : 0;
        }
    }

    vector<State> beam;
    State init_state;
    init_state.score = 0;
    beam.push_back(init_state);

    for (int k = 0; k < K; ++k) {
        vector<State> next_beam;
        
        for (const auto& state : beam) {
            // 現在の状態に対する s_ci を計算
            compute_s_ci(in, state.all_flights, mg_ws.current_s_ci, c_ws);

            memset(mg_ws.edge_gain, 0, sizeof(mg_ws.edge_gain));

            // 高速な枝刈り付き限界利益計算
            for (int u = 0; u < N; ++u) {
                for (int dst = 0; dst < N; ++dst) {
                    if (!in.valid_dist[u][dst]) continue;
                    long long w = in.weight[u][dst];
                    for (int T_idx = 0; T_idx < 21; ++T_idx) {
                        int s_sq_val = rival.s_sq[u][dst][T_idx];
                        int s_ci_val = mg_ws.current_s_ci[u][dst][T_idx];
                        
                        if (s_ci_val > s_sq_val) continue;

                        int min_s = max(0, s_sq_val + 1);
                        int T_target = 60 + T_idx * 6;

                        for (int v = 0; v < N; ++v) {
                            if (u == v) continue;
                            int max_t = (v == dst) ? T_target : mg_ws.current_s_ci[v][dst][T_idx];
                            if (max_t < 0) continue;
                            
                            int max_s = min(180 - mg_ws.dur[u][v], max_t - mg_ws.dur[u][v]);

                            if (min_s <= max_s) {
                                mg_ws.edge_gain[u][v][min_s] += w;
                                if (max_s + 1 <= 180) {
                                    mg_ws.edge_gain[u][v][max_s + 1] -= w;
                                }
                            }
                        }
                    }
                }
            }

            // いもす法展開
            for (int u = 0; u < N; ++u) {
                for (int v = 0; v < N; ++v) {
                    if (u == v) continue;
                    long long current_val = 0;
                    long long* gain_uv = mg_ws.edge_gain[u][v];
                    for (int s = 0; s <= 180; ++s) {
                        current_val += gain_uv[s];
                        gain_uv[s] = current_val;
                    }
                }
            }

            memset(mg_ws.dp, 0, sizeof(mg_ws.dp));
            memset(mg_ws.nxt_v, -1, sizeof(mg_ws.nxt_v));

            for (int t = 180; t >= 0; --t) {
                for (int u = 0; u < N; ++u) {
                    if (t + 1 <= 180) {
                        mg_ws.dp[u][t] = mg_ws.dp[u][t + 1];
                        mg_ws.nxt_v[u][t] = -1;
                    }
                    for (int v = 0; v < N; ++v) {
                        if (u == v) continue;
                        int t_next = t + mg_ws.dur[u][v];
                        if (t_next <= 180) {
                            long long gain = mg_ws.edge_gain[u][v][t];
                            if (gain + mg_ws.dp[v][t_next] > mg_ws.dp[u][t]) {
                                mg_ws.dp[u][t] = gain + mg_ws.dp[v][t_next];
                                mg_ws.nxt_v[u][t] = v;
                            }
                        }
                    }
                }
            }

            // ★ ビームサーチの分岐: 利益が高い順に開始都市をソート
            vector<pair<long long, int>> start_cands;
            for (int u = 0; u < N; ++u) {
                start_cands.emplace_back(mg_ws.dp[u][0], u);
            }
            sort(start_cands.rbegin(), start_cands.rend());

            // 上位 BRANCHES 個からパスを生成
            for (int i = 0; i < min(N, BRANCHES); ++i) {
                int start_u = start_cands[i].second;

                vector<Flight> path;
                int curr_u = start_u, curr_t = 0;
                while (curr_t <= 180) {
                    int next_city = mg_ws.nxt_v[curr_u][curr_t];
                    if (next_city != -1) {
                        int d = mg_ws.dur[curr_u][next_city];
                        path.push_back({curr_u, curr_t, next_city, curr_t + d});
                        curr_u = next_city;
                        curr_t += d;
                    } else {
                        curr_t += 1;
                    }
                }

                State nstate = state;
                nstate.current_paths.push_back(path);
                for (const auto& f : path) {
                    nstate.all_flights.push_back(f);
                }

                // 分岐先のスコアを正確に評価
                compute_s_ci(in, nstate.all_flights, mg_ws.current_s_ci, c_ws);
                long long v_ci = 0;
                for (int u = 0; u < N; ++u) {
                    for (int dst = 0; dst < N; ++dst) {
                        if (!in.valid_dist[u][dst]) continue;
                        long long w = in.weight[u][dst];
                        for (int T_idx = 0; T_idx < 21; ++T_idx) {
                            if (mg_ws.current_s_ci[u][dst][T_idx] > rival.s_sq[u][dst][T_idx]) {
                                v_ci += w;
                            }
                        }
                    }
                }
                nstate.score = (long long)(1000000.0 * (double)v_ci / in.total_W);
                next_beam.push_back(nstate);
            }
        }

        // ビームをスコア降順でソート
        sort(next_beam.begin(), next_beam.end());
        
        beam.clear();
        for (const auto& s : next_beam) {
            // スコアが全く同じ（実質的に同じルート）ものは除外して多様性を確保
            if (beam.empty() || beam.back().score != s.score) {
                beam.push_back(s);
                if (beam.size() == BEAM_WIDTH) break;
            }
        }
    }

    return beam[0].current_paths;
}

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

    Input in;
    auto c_ws = make_unique<ComputeSciWorkspace>();
    auto mg_ws = make_unique<MarginalGainWorkspace>();

    RivalSchedule rival(in, *c_ws);
    auto circle_flights = solve_beam_search(in, rival, *c_ws, *mg_ws);

    for (int k = 0; k < K; ++k) {
        cout << circle_flights[k].size() << "\n";
        for (const auto& f : circle_flights[k]) {
            cout << f.u + 1 << " " << to_time_str(f.s) << " " << f.v + 1 << " " << to_time_str(f.t) << "\n";
        }
    }

    return 0;
}