graph/minimumsteinertree.hpp

• View this file on GitHub
• Last update: 2025-05-09 08:40:27+00:00
• Include: #include "graph/minimumsteinertree.hpp"

Depends on

• graph/graphtemplate.hpp

Verified with

• verify/yosupo-minimum_steiner_tree.test.cpp

Code

#pragma once
#include "graphtemplate"
#include <bits/stdc++.h>
using namespace std;
template <class T, bool directed = false, bool weighted = true>
graph<T, false, true> minimumsteinertree(graph<T, directed, weighted> &g, vector<int> &v) {
    vector<vector<T>> dp(1<<v.size(), vector<T>(g.size(), numeric_limits<T>::max()));
    vector<vector<T>> d(g.size(), vector<T>(g.size(), numeric_limits<T>::max()));
    vector<vector<int>> id(g.size(), vector<int>(g.size(), -1));
    vector<vector<pair<int, int>>> par(1<<v.size(), vector<pair<int, int>>(g.size(), {-1, -1}));
    for (auto& _e : g._edges) {
        if (_e.cost < d[_e.from][_e.to]) {
            d[_e.from][_e.to] = _e.cost;
            d[_e.to][_e.from] = _e.cost;
            id[_e.from][_e.to] = _e.id;
            id[_e.to][_e.from] = _e.id;
        }
    }
    for (int i=0; i<g.size(); i++) {
        d[i][i] = 0;
    }
    for (int i=0; i<v.size(); i++) {
        dp[1<<i][v[i]] = 0;
    }
    for (int i=0; i<(1<<v.size()); i++) {
        for (int j=i; 0<j; j=(j-1)&i) {
            for (int k=0; k<g.size(); k++) {
                if (dp[j][k] == numeric_limits<T>::max() || dp[i^j][k] == numeric_limits<T>::max()) continue;
                if (dp[j][k] + dp[i^j][k] < dp[i][k]) {
                    dp[i][k] = dp[j][k] + dp[i^j][k];
                    par[i][k] = {0, j};
                }
            }
        }
        priority_queue<pair<T, int>, vector<pair<T, int>>, greater<pair<T, int>>> q;
        for (int j=0; j<g.size(); j++) {
            if (dp[i][j] != numeric_limits<T>::max()) q.push({dp[i][j], j});
        }
        while (!q.empty()) {
            auto [x, y] = q.top(); q.pop();
            if (dp[i][y] < x) continue;
            for (auto& _e : g[y]) {
                if (x + _e.cost < dp[i][_e.to]) {
                    dp[i][_e.to] = x + _e.cost;
                    q.push({x + _e.cost, _e.to});
                    par[i][_e.to] = {1, _e.from};
                }
            }
        }
    }
    int c = -1;
    T ans = numeric_limits<T>::max();
    for (int i=0; i<g.size(); i++) {
        if (dp.back()[i] < ans) {
            ans = dp.back()[i];
            c = i;
        }
    }
    graph<T, false, true> res(g.size());
    vector<int> used(g._edges.size());
    if (c == -1) return res;
    stack<pair<int, int>> s;
    s.push({(1<<v.size())-1, c});
    while (!s.empty()) {
        auto [x, y] = s.top(); s.pop();
        auto [X, Y] = par[x][y];
        if (X == -1) continue;
        else if (X == 0) {
            s.push({Y, y});
            s.push({x^Y, y});
        } else if (X == 1) {
            s.push({x, Y});
            int z = id[y][Y];
            if (!used[z]) {
                used[z] = 1;
                res.add_edge(g._edges[z]);
            }
        }
    }
    return res;
}

#line 2 "graph/graphtemplate.hpp"
#include<bits/stdc++.h>
using namespace std;
// 辺の構造体 edge(from, to, cost, id)
template<class T = int>
struct edge {
    int from, to;
    T cost;
    int id;
};
// 頂点の構造体 vector<edge<T>>
template<class T = int>
using edges = vector<edge<T>>;
// グラフの構造体 graph<T, directed, weighted> 
template <class T = int, bool directed = false, bool weighted = false>
struct graph {
    bool isdirected, isweighted;
    edges<T> _edges;
    vector<edges<T>> data;
    T sumcost;
    graph() = default;
    // 頂点数 n のグラフを作成する
    graph(int n) : isdirected(directed), isweighted(weighted), data(n), sumcost(T{}) {}
    // from から to へ辺を追加する
    void add_edge(int from, int to, T cost = 1, int id = -1) {
        if (id == -1) id = _edges.size();
        data[from].push_back(edge<T>(from, to, cost, id));
        _edges.push_back(edge<T>(from, to, cost, id));
        if (!isdirected) {
            data[to].push_back(edge<T>(to, from, cost, id));
        }
        sumcost += cost;
    }
    // 辺を追加する
    void add_edge(edge<T> _e) {
        add_edge(_e.from, _e.to, _e.cost, _e.id);
    }
    // 標準入力から辺を読み込む
    void read(int m, int indexed = 1) {
        for (int i=0; i<m; i++) {
            int from, to;
            T cost = 1;
            cin >> from >> to;
            if (isweighted) cin >> cost;
            add_edge(from - indexed, to - indexed, cost);
        }
    }
    // 頂点数を返す
    int size() {
        return data.size();
    }
    // 頂点を返す
    edges<T> operator[](int k) {
        return data[k];
    }
    vector<int> path_to_vertex(edges<T>& _e) {
        vector<int> re;
        if (_e.size() == 0) {
            return re;
        }
        if (_e.size() == 1) {
            re.push_back(_e[0].from);
            re.push_back(_e[0].to);
            return re;
        }
        int x=_e[0].from,y=_e[0].to;
        if (x==_e[1].to || x == _e[1].from) swap(x, y);
        re.push_back(x);
        for (int i=1; i<_e.size(); i++) {
            re.push_back(y);
            x = _e[i].to;
            if (x == y) x = _e[i].from;
            swap(x, y);
        }
        return re;
    }
    edges<T> vetex_to_path (vector<int>& v){
        edges<T> re;
        for (int i=0; i+1<v.size(); i++) {
            for (auto& _e : this[v[i]]) {
                if (_e.to == v[i+1]) {
                    re.push_back(_e);
                    break;
                }
            }
        }
        return re;
    }
};
#line 4 "graph/minimumsteinertree.hpp"
using namespace std;
template <class T, bool directed = false, bool weighted = true>
graph<T, false, true> minimumsteinertree(graph<T, directed, weighted> &g, vector<int> &v) {
    vector<vector<T>> dp(1<<v.size(), vector<T>(g.size(), numeric_limits<T>::max()));
    vector<vector<T>> d(g.size(), vector<T>(g.size(), numeric_limits<T>::max()));
    vector<vector<int>> id(g.size(), vector<int>(g.size(), -1));
    vector<vector<pair<int, int>>> par(1<<v.size(), vector<pair<int, int>>(g.size(), {-1, -1}));
    for (auto& _e : g._edges) {
        if (_e.cost < d[_e.from][_e.to]) {
            d[_e.from][_e.to] = _e.cost;
            d[_e.to][_e.from] = _e.cost;
            id[_e.from][_e.to] = _e.id;
            id[_e.to][_e.from] = _e.id;
        }
    }
    for (int i=0; i<g.size(); i++) {
        d[i][i] = 0;
    }
    for (int i=0; i<v.size(); i++) {
        dp[1<<i][v[i]] = 0;
    }
    for (int i=0; i<(1<<v.size()); i++) {
        for (int j=i; 0<j; j=(j-1)&i) {
            for (int k=0; k<g.size(); k++) {
                if (dp[j][k] == numeric_limits<T>::max() || dp[i^j][k] == numeric_limits<T>::max()) continue;
                if (dp[j][k] + dp[i^j][k] < dp[i][k]) {
                    dp[i][k] = dp[j][k] + dp[i^j][k];
                    par[i][k] = {0, j};
                }
            }
        }
        priority_queue<pair<T, int>, vector<pair<T, int>>, greater<pair<T, int>>> q;
        for (int j=0; j<g.size(); j++) {
            if (dp[i][j] != numeric_limits<T>::max()) q.push({dp[i][j], j});
        }
        while (!q.empty()) {
            auto [x, y] = q.top(); q.pop();
            if (dp[i][y] < x) continue;
            for (auto& _e : g[y]) {
                if (x + _e.cost < dp[i][_e.to]) {
                    dp[i][_e.to] = x + _e.cost;
                    q.push({x + _e.cost, _e.to});
                    par[i][_e.to] = {1, _e.from};
                }
            }
        }
    }
    int c = -1;
    T ans = numeric_limits<T>::max();
    for (int i=0; i<g.size(); i++) {
        if (dp.back()[i] < ans) {
            ans = dp.back()[i];
            c = i;
        }
    }
    graph<T, false, true> res(g.size());
    vector<int> used(g._edges.size());
    if (c == -1) return res;
    stack<pair<int, int>> s;
    s.push({(1<<v.size())-1, c});
    while (!s.empty()) {
        auto [x, y] = s.top(); s.pop();
        auto [X, Y] = par[x][y];
        if (X == -1) continue;
        else if (X == 0) {
            s.push({Y, y});
            s.push({x^Y, y});
        } else if (X == 1) {
            s.push({x, Y});
            int z = id[y][Y];
            if (!used[z]) {
                used[z] = 1;
                res.add_edge(g._edges[z]);
            }
        }
    }
    return res;
}

Back to top page