#592522 (C++17(clang)) No.1316 Maximum Minimum Spanning Tree

• 問題: No.1316 Maximum Minimum Spanning Tree
• ユーザー: hitonanode
• 提出日時: 2020-12-12 21:29:41
• 言語: C++17(clang) (17.0.6 + boost 1.89.0)
• 結果: WA
• 実行時間: -
• コード長: 12,744 bytes
• コンパイル時間: 2,160 ms
• コンパイル使用メモリ: 143,672 KB
• 実行使用メモリ: 6,948 KB
• 最終ジャッジ日時: 2024-09-19 22:22:20
• 合計ジャッジ時間: 146,411 ms
• ジャッジサーバーID （参考情報）: judge5 / judge4

ファイルパターン	結果
sample	AC * 4
other	AC * 59 WA * 19

コンパイルメッセージ

main.cpp:403:25: warning: unqualified call to 'std::move' [-Wunqualified-std-cast-call]
  403 |                 depth = move(depth_);
      |                         ^
      |                         std::
main.cpp:404:26: warning: unqualified call to 'std::move' [-Wunqualified-std-cast-call]
  404 |                 parent = move(parent_);
      |                          ^
      |                          std::
main.cpp:405:24: warning: unqualified call to 'std::move' [-Wunqualified-std-cast-call]
  405 |                 edge = move(edge_);
      |                        ^
      |                        std::
3 warnings generated.

ソースコード

// https://yukicoder.me/submissions/592136 改変 山登り法
// #pragma GCC optimize("Ofast")
// #pragma GCC optimize("unroll-loops")
#include <algorithm>
#include <array>
#include <cassert>
#include <chrono>
#include <iostream>
#include <numeric>
#include <queue>
#include <tuple>
#include <vector>

//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// MinCostFlow by https://github.com/yosupo06/library-checker-problems/blob/master/graph/min_cost_b_flow/sol/ssp.cpp
// Copyright 2020 yosupo06/library-checker-problems
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

enum Objective {
    MINIMIZE = 1,
    MAXIMIZE = -1,
};
enum class Status {
    OPTIMAL,
    INFEASIBLE,
};

template <class Flow, class Cost, Objective objective = Objective::MINIMIZE>
class MinCostFlow {
    using V_id = uint32_t;
    using E_id = uint32_t;

    class Edge {
        friend class MinCostFlow;

        V_id src, dst;
        Flow flow, cap;
        Cost cost;
        E_id rev;

        public:
        Edge() = default;

        Edge(const V_id src, const V_id dst, const Flow cap, const Cost cost,
             const E_id rev)
            : src(src), dst(dst), flow(0), cap(cap), cost(cost), rev(rev) {}

        [[nodiscard]] Flow residual_cap() const { return cap - flow; }
    };

    public:
    class EdgePtr {
        friend class MinCostFlow;

        const MinCostFlow *instance;
        V_id v;
        E_id e;

        EdgePtr(const MinCostFlow * const instance, const V_id v, const E_id e)
            : instance(instance), v(v), e(e) {}

        [[nodiscard]] const Edge &edge() const { return instance->g[v][e]; }

        [[nodiscard]] const Edge &rev() const {
            const Edge &e = edge();
            return instance->g[e.dst][e.rev];
        }

        public:
        EdgePtr() = default;

        [[nodiscard]] V_id src() const { return v; }

        [[nodiscard]] V_id dst() const { return edge().dst; }

        [[nodiscard]] Flow flow() const { return edge().flow; }

        [[nodiscard]] Flow lower() const { return -rev().cap; }

        [[nodiscard]] Flow upper() const { return edge().cap; }

        [[nodiscard]] Cost cost() const { return edge().cost; }

        [[nodiscard]] Cost gain() const { return -edge().cost; }
    };

    private:
    V_id n;
    std::vector<std::vector<Edge>> g;
    std::vector<Flow> b;

    public:
    MinCostFlow() : n(0) {}

    V_id add_vertex() {
        ++n;
        g.resize(n);
        b.resize(n);
        return n-1;
    }

    std::vector<V_id> add_vertices(const size_t size) {
        std::vector<V_id> ret(size);
        std::iota(std::begin(ret), std::end(ret), n);
        n += size;
        g.resize(n);
        b.resize(n);
        return ret;
    }

    EdgePtr add_edge(const V_id src, const V_id dst, const Flow lower,
                     const Flow upper, const Cost cost) {
        const E_id e = g[src].size(), re = src == dst ? e + 1 : g[dst].size();
        assert(lower <= upper);
        g[src].emplace_back(Edge{src, dst, upper, cost * objective, re});
        g[dst].emplace_back(Edge{dst, src, -lower, -cost * objective, e});
        return EdgePtr{this, src, e};
    }

    void add_supply(const V_id v, const Flow amount) { b[v] += amount; }

    void add_demand(const V_id v, const Flow amount) { b[v] -= amount; }

    private:
    // Variables used in calculation
    const Cost unreachable = std::numeric_limits<Cost>::max();
    Cost farthest;
    std::vector<Cost> potential;
    std::vector<Cost> dist;
    std::vector<Edge *> parent; // out-forrest.
    std::priority_queue<std::pair<Cost, int>, std::vector<std::pair<Cost, int>>,
        std::greater<>>
            pq; // should be empty outside of dual()
    std::vector<V_id> excess_vs, deficit_vs;

    Edge &rev(const Edge &e) { return g[e.dst][e.rev]; }

    void push(Edge &e, const Flow amount) {
        e.flow += amount;
        g[e.dst][e.rev].flow -= amount;
    }

    Cost residual_cost(const V_id src, const V_id dst, const Edge &e) {
        return e.cost + potential[src] - potential[dst];
    }

    bool dual(const Flow delta) {
        dist.assign(n, unreachable);
        parent.assign(n, nullptr);
        excess_vs.erase(std::remove_if(std::begin(excess_vs), std::end(excess_vs),
                                       [&](const V_id v) { return b[v] < delta; }),
                        std::end(excess_vs));
        deficit_vs.erase(std::remove_if(std::begin(deficit_vs),
                                        std::end(deficit_vs),
                                        [&](const V_id v) { return b[v] > -delta; }),
                         std::end(deficit_vs));
        for (const auto v : excess_vs) pq.emplace(dist[v] = 0, v);
        farthest = 0;
        std::size_t deficit_count = 0;
        while (!pq.empty()) {
            Cost d;
            std::size_t u;
            std::tie(d, u) = pq.top();
            // const auto [d, u] = pq.top();
            pq.pop();
            if (dist[u] < d) continue;
            farthest = d;
            if (b[u] <= -delta) ++deficit_count;
            if (deficit_count >= deficit_vs.size()) break;
            for (auto &e : g[u]) {
                if (e.residual_cap() < delta) continue;
                const auto v = e.dst;
                const auto new_dist = d + residual_cost(u, v, e);
                if (new_dist >= dist[v]) continue;
                pq.emplace(dist[v] = new_dist, v);
                parent[v] = &e;
            }
        }
        pq = decltype(pq)(); // pq.clear() doesn't exist.
        for (V_id v = 0; v < n; ++v) {
            potential[v] += std::min(dist[v], farthest);
        }
        return deficit_count > 0;
    }

    void primal(const Flow delta) {
        for (const auto t : deficit_vs) {
            if (dist[t] > farthest) continue;
            Flow f = -b[t];
            V_id v;
            for (v = t; parent[v] != nullptr; v = parent[v]->src) {
                f = std::min(f, parent[v]->residual_cap());
            }
            f = std::min(f, b[v]);
            f -= f % delta;
            if (f <= 0) continue;
            for (v = t; parent[v] != nullptr;) {
                auto &e = *parent[v];
                push(e, f);
                int u = parent[v]->src;
                if (e.residual_cap() <= 0) parent[v] = nullptr;
                v = u;
            }
            b[t] += f;
            b[v] -= f;
        }
    }

    void saturate_negative(const Flow delta) {
        for (auto &es : g) for (auto &e : es) {
            Flow rcap = e.residual_cap();
            rcap -= rcap % delta;
            const Cost rcost = residual_cost(e.src, e.dst, e);
            if (rcost < 0 || rcap < 0) {
                push(e, rcap);
                b[e.src] -= rcap;
                b[e.dst] += rcap;
            }
        }
        for (V_id v = 0; v < n; ++v) if (b[v] != 0) {
            (b[v] > 0 ? excess_vs : deficit_vs).emplace_back(v);
        }
    }

    public:
    std::pair<Status, Cost> solve() {
        potential.resize(n);

        saturate_negative(1);
        while (dual(1)) primal(1);

        Cost value = 0;
        for (const auto &es : g) for (const auto &e : es) {
            value += e.flow * e.cost;
        }
        value /= 2;

        if (excess_vs.empty() && deficit_vs.empty()) {
            return { Status::OPTIMAL, value / objective };
        } else {
            return { Status::INFEASIBLE, value / objective };
        }
    }

    std::vector<Cost> get_potential() {
        return potential;
    }
    template<class T>
    T get_result_value() {
        T value = 0;
        for (const auto &es : g) for (const auto &e : es) {
            value += (T)(e.flow) * (T)(e.cost);
        }
        value /= (T)2;
        return value;
    }
};
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// MinCostFlow END
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////


using namespace std;
constexpr long long INFll = (1ll << 60) - 1;

int n, m, k;
vector<array<int, 4>> E;
vector<int> depth, parent, edge;
vector<int> depth_, parent_, edge_;
void build(const vector<int> &v) {
    depth.assign(n, 0);
    parent.assign(n, 0);
    edge.assign(n, 0);
    // depth[0] = 0;
    vector<vector<pair<int, int>>> g(n);
    for(int i : v){
        int a = E[i][0], b = E[i][1];
        g[a].emplace_back(b, i);
        g[b].emplace_back(a, i);
    }
    auto dfs = [&](int from, int at, auto dfs) -> void {
        const int d2 = depth[at] + 1;
        for(const auto& [to, id] : g[at]) if(to != from){
            depth[to] = d2;
            parent[to] = at;
            edge[to] = id;
            dfs(at, to, dfs);
        }
    };
    dfs(-1, 0, dfs);
}
vector<int> exchangable_edges(int e){
    vector<int> ans;
    int a = E[e][0], b = E[e][1];
    while(a != b){
        if(depth[a] < depth[b]) swap(a, b);
        ans.push_back(edge[a]);
        a = parent[a];
    }
    return ans;
}

int64_t flow_us;
// 「G の最小全域木として vl が採用されうる」という条件の下で x を動かしたときの最大値を求める．
long long subsolve(const vector<int>& vl, const vector<int>& vr){
    long long res = 0;
    for (auto i : vl) res += E[i][2];
    res *= k;
    
    build(vl);
    
    MinCostFlow<long long, long long> mcf;
    const auto vs = mcf.add_vertices(m + 2);
    int s = m, t = m + 1;
    const long long BIG = 1e10;
    
    for (auto j : vr) {
        // vl に含まれない辺 j を追加したとき，代わりに取り除ける辺 i を列挙
        for (auto i : exchangable_edges(j)) mcf.add_edge(i, j, 0, BIG, E[j][2] - E[i][2]);
    }
    for (auto i : vl) mcf.add_edge(s, i, max(k - E[i][3], 0), BIG, 0);
    for (auto j : vr) mcf.add_edge(j, t, 0, E[j][3], 0);
    mcf.add_edge(t, s, 0, BIG, 0);

    auto START = std::chrono::system_clock::now();
    auto [status, f] = mcf.solve();
    flow_us += std::chrono::duration_cast<std::chrono::microseconds>(std::chrono::system_clock::now() - START).count();

    if (status == Status::INFEASIBLE) return INFll;
    res += f;
    return res;
}

struct UnionFind{
    vector<int> data;
    UnionFind(int n): data(n, -1){}
    bool unite(int a, int b){
        a = root(a); b = root(b);
        if(a == b) return 0;
        if(data[a] > data[b]) swap(a, b);
        data[a] += data[b];
        data[b] = a;
        return 1;
    }
    bool find(int a, int b){ return root(a) == root(b); }
    int root(int a){ return data[a] < 0 ? a : data[a] = root(data[a]); }
    int size(int a){ return -data[root(a)]; }
    int operator[](int a){ return root(a); }
};

struct Xorshift64{
    using result_type = uint32_t;
    static constexpr result_type min(){ return 0; }
    static constexpr result_type max(){ return -1; }
    uint64_t x = 1;
    result_type operator()(){
        x ^= (x << 13);
        x ^= (x >> 7);
        x ^= (x << 17);
        return static_cast<uint32_t>(x);
    }
}rnd;


int main(){
    cin.tie(nullptr);
    ios::sync_with_stdio(false);
    auto start = chrono::system_clock::now();

    cin >> n >> m >> k;
    E.resize(m);
    for (auto& [a, b, c, d] : E){
        cin >> a >> b >> c >> d;
        a--; b--;
    }
    sort(E.begin(), E.end(), [](const auto& a, const auto& b){ return a[2] < b[2]; });
    UnionFind uf(n);
    vector<int> vl, vr;
    for (int i = 0; i < m; i++) {
        int a = E[i][0], b = E[i][1];
        if(uf.unite(a, b)) vl.push_back(i);
        else vr.push_back(i);
    }

    long long ans = subsolve(vl, vr);

    int nbupd = 0, nb_loop = 0;
    if(vr.size()){
        while (chrono::duration_cast<chrono::milliseconds>(chrono::system_clock::now() - start).count() < 1950) {
            nb_loop++;
            auto p = vr.begin() + rnd() % vr.size();
            auto v = exchangable_edges(*p);
            auto q = find(vl.begin(), vl.end(), v[rnd() % v.size()]);
            iter_swap(p, q);
            swap(depth, depth_);
            swap(parent, parent_);
            swap(edge, edge_);
            const long long ans2 = subsolve(vl, vr);
            if (ans < ans2) nbupd++;
            ans = max(ans, ans2);
            if(ans != ans2){
                iter_swap(p, q);
                depth = move(depth_);
                parent = move(parent_);
                edge = move(edge_);
            }
        }
        
    }
    cerr << nb_loop << ' ' << nbupd << ' ' << flow_us / 1000 << '\n';

    if (ans == INFll) ans = -1;
    cout << ans << '\n';
}