library
This documentation is automatically generated by online-judge-tools/verification-helper
Matroid
(Algorithm/matroid.hpp)
- View this file on GitHub
- Last update: 2024-04-26 03:18:17+09:00
- Include:
#include "Algorithm/matroid.hpp"
Depends on
Bipartite Matching
(Graph/bimatching.hpp)
- DM decomposition
(Graph/dmdecomp.hpp)
- Strongly Connected Components
(Graph/scc.hpp)
- Random
(Utility/random.hpp)
Code
#pragma once
#include "Graph/bimatching.hpp"
#include "Graph/dmdecomp.hpp"
struct GraphicMatroid {
using P = pair<int, int>;
int n;
vector<P> E;
vector<vector<P>> G;
vector<int> allow;
GraphicMatroid() {}
GraphicMatroid(int n, vector<P> &E) : n(n), E(E), G(n) {
rep(i, 0, SZ(E)) {
auto [u, v] = E[i];
G[u].push_back({v, i});
G[v].push_back({u, i});
}
}
void set(vector<int> &I) {
allow = I;
}
vector<int> circuit(int e) {
auto [x, y] = E[e];
vector<int> ret;
ret.push_back(e);
auto dfs = [&](auto &dfs, int v, int p) -> bool {
if (v == y)
return true;
for (auto &[to, i] : G[v])
if (allow[i] and to != p) {
ret.push_back(i);
if (dfs(dfs, to, v))
return true;
ret.pop_back();
}
return false;
};
if (dfs(dfs, x, -1))
return ret;
else
return {};
}
};
struct PartitionMatroid {
vector<int> grp; // -1:not assign
vector<int> lim;
vector<vector<int>> cnt;
PartitionMatroid() {}
PartitionMatroid(vector<int> &grp, vector<int> lim) : grp(grp), lim(lim) {}
void set(vector<int> &I) {
cnt.assign(SZ(lim), {});
rep(i, 0, SZ(grp)) if (I[i] != 0 and grp[i] != -1) {
cnt[grp[i]].push_back(i);
}
}
vector<int> circuit(int e) {
if (grp[e] == -1)
return {};
if (SZ(cnt[grp[e]]) + 1 > lim[grp[e]]) {
vector<int> ret = cnt[grp[e]];
ret.push_back(e);
return ret;
} else
return {};
}
};
struct TransversalMatroid {
using P = pair<int, int>;
int n, m;
vector<vector<int>> G, aug;
vector<int> match, fixed;
TransversalMatroid() {}
TransversalMatroid(int n, int m, vector<P> &E)
: n(n), m(m), G(n), aug(n), fixed(n) {
for (auto &[u, v] : E)
G[u].push_back(v);
}
void set(vector<int> &I) {
vector g(n, vector<int>());
rep(e, 0, n) if (I[e]) {
for (auto &to : G[e])
g[e].push_back(to);
}
auto match = BiMatching(n, m, g);
auto dm = DMdecomposition(n, m, g, match);
fixed.assign(m, 1);
for (auto &v : dm.front())
if (v >= n)
fixed[v - n] = 0;
aug.assign(n + m, {});
rep(e, 0, n) {
for (auto &to : G[e])
aug[e].push_back(to + n);
}
rep(i, 0, n) if (match[i] != -1) aug[match[i] + n].push_back(i);
}
vector<int> circuit(int e) {
for (auto &to : G[e])
if (!fixed[to]) {
return {};
}
vector<int> used(n + m);
queue<int> que;
used[e] = 1;
que.push(e);
while (!que.empty()) {
int v = que.front();
que.pop();
for (auto &to : aug[v])
if (!used[to]) {
used[to] = 1;
que.push(to);
}
}
vector<int> ret;
rep(i, 0, n) if (used[i]) ret.push_back(i);
return ret;
}
};
template <typename M1, typename M2, typename Val>
pair<vector<Val>, vector<vector<int>>>
MinimumMatroidIntersection(M1 &m1, M2 &m2, vector<Val> &ws) {
using P = pair<Val, int>;
int n = SZ(ws);
vector<Val> ret1;
vector<vector<int>> ret2;
Val profit = 0;
vector<int> I(n);
ret1.push_back(profit);
ret2.push_back(I);
for (;;) {
vector G(n + 2, vector<P>());
int S = n, T = n + 1;
m1.set(I);
m2.set(I);
rep(e, 0, n) {
if (I[e])
continue;
auto c1 = m1.circuit(e);
if (c1.empty())
G[S].push_back({ws[e] * (n + 1) + 1, e});
else
for (auto &to : c1)
if (e != to) {
G[to].push_back({ws[e] * (n + 1) + 1, e});
}
auto c2 = m2.circuit(e);
if (c2.empty())
G[e].push_back({Val(0), T});
else
for (auto &to : c2)
if (e != to) {
G[e].push_back({-ws[to] * (n + 1) + 1, to});
}
}
vector<ll> dist(n + 2, INF);
dist[S] = 0;
vector<int> pre(n + 2, -1);
for (;;) {
bool upd = 0;
rep(v, 0, n + 2) if (dist[v] != INF) {
for (auto &[cost, to] : G[v]) {
if (chmin(dist[to], dist[v] + cost)) {
pre[to] = v;
upd = 1;
}
}
}
if (!upd)
break;
}
if (dist[T] == INF)
break;
int cur = T;
while (cur != S) {
cur = pre[cur];
if (cur != S) {
I[cur] ^= 1;
if (I[cur])
profit += ws[cur];
else
profit -= ws[cur];
}
}
ret1.push_back(profit);
ret2.push_back(I);
}
return {ret1, ret2};
}
/**
* @brief Matroid
*/#line 2 "Utility/random.hpp"
namespace Random {
mt19937_64 randgen(chrono::steady_clock::now().time_since_epoch().count());
using u64 = unsigned long long;
u64 get() {
return randgen();
}
template <typename T> T get(T L) { // [0,L]
return get() % (L + 1);
}
template <typename T> T get(T L, T R) { // [L,R]
return get(R - L) + L;
}
double uniform() {
return double(get(1000000000)) / 1000000000;
}
string str(int n) {
string ret;
rep(i, 0, n) ret += get('a', 'z');
return ret;
}
template <typename Iter> void shuffle(Iter first, Iter last) {
if (first == last)
return;
int len = 1;
for (auto it = first + 1; it != last; it++) {
len++;
int j = get(0, len - 1);
if (j != len - 1)
iter_swap(it, first + j);
}
}
template <typename T> vector<T> select(int n, T L, T R) { // [L,R]
if (n * 2 >= R - L + 1) {
vector<T> ret(R - L + 1);
iota(ALL(ret), L);
shuffle(ALL(ret));
ret.resize(n);
return ret;
} else {
unordered_set<T> used;
vector<T> ret;
while (SZ(used) < n) {
T x = get(L, R);
if (!used.count(x)) {
used.insert(x);
ret.push_back(x);
}
}
return ret;
}
}
void relabel(int n, vector<pair<int, int>> &es) {
shuffle(ALL(es));
vector<int> ord(n);
iota(ALL(ord), 0);
shuffle(ALL(ord));
for (auto &[u, v] : es)
u = ord[u], v = ord[v];
}
template <bool directed, bool simple> vector<pair<int, int>> genGraph(int n) {
vector<pair<int, int>> cand, es;
rep(u, 0, n) rep(v, 0, n) {
if (simple and u == v)
continue;
if (!directed and u > v)
continue;
cand.push_back({u, v});
}
int m = get(SZ(cand));
vector<int> ord;
if (simple)
ord = select(m, 0, SZ(cand) - 1);
else {
rep(_, 0, m) ord.push_back(get(SZ(cand) - 1));
}
for (auto &i : ord)
es.push_back(cand[i]);
relabel(n, es);
return es;
}
vector<pair<int, int>> genTree(int n) {
vector<pair<int, int>> es;
rep(i, 1, n) es.push_back({get(i - 1), i});
relabel(n, es);
return es;
}
}; // namespace Random
/**
* @brief Random
*/
#line 3 "Graph/bimatching.hpp"
vector<int> BiMatching(int n, int m, vector<vector<int>> &g) {
rep(v, 0, n) Random::shuffle(ALL(g[v]));
vector<int> L(n, -1), R(m, -1);
for (;;) {
vector<int> par(n, -1), root(n, -1);
queue<int> que;
rep(i, 0, n) if (L[i] == -1) {
que.push(i);
root[i] = i;
}
bool upd = 0;
while (!que.empty()) {
int v = que.front();
que.pop();
if (L[root[v]] != -1)
continue;
for (auto u : g[v]) {
if (R[u] == -1) {
while (u != -1) {
R[u] = v;
swap(L[v], u);
v = par[v];
}
upd = 1;
break;
}
int to = R[u];
if (par[to] == -1) {
root[to] = root[v];
par[to] = v;
que.push(to);
}
}
}
if (!upd)
break;
}
return L;
}
/**
* @brief Bipartite Matching
*/
#line 2 "Graph/scc.hpp"
struct SCC{
int n,m,cur;
vector<vector<int>> g;
vector<int> low,ord,id;
SCC(int _n=0):n(_n),m(0),cur(0),g(_n),low(_n),ord(_n,-1),id(_n){}
void resize(int _n){
n=_n;
g.resize(n);
low.resize(n);
ord.resize(n,-1);
id.resize(n);
}
void add_edge(int u,int v){g[u].emplace_back(v);}
void dfs(int v,vector<int>& used){
ord[v]=low[v]=cur++;
used.emplace_back(v);
for(auto& nxt:g[v]){
if(ord[nxt]==-1){
dfs(nxt,used); chmin(low[v],low[nxt]);
}
else{
chmin(low[v],ord[nxt]);
}
}
if(ord[v]==low[v]){
while(1){
int add=used.back(); used.pop_back();
ord[add]=n; id[add]=m;
if(v==add)break;
}
m++;
}
}
void run(){
vector<int> used;
rep(v,0,n)if(ord[v]==-1)dfs(v,used);
for(auto& x:id)x=m-1-x;
}
};
/**
* @brief Strongly Connected Components
*/
#line 4 "Graph/dmdecomp.hpp"
vector<vector<int>> DMdecomposition(int n, int m, vector<vector<int>> &g,
vector<int> match) {
if (match.empty())
match = BiMatching(n, m, g);
vector G(n + m, vector<int>()), REV(n + m, vector<int>());
rep(i, 0, n) for (auto &j : g[i]) {
G[i].push_back(j + n);
REV[j + n].push_back(i);
}
vector<int> R(m, -1);
rep(i, 0, n) if (match[i] != -1) {
G[match[i] + n].push_back(i);
REV[i].push_back(match[i] + n);
R[match[i]] = i;
}
vector<int> V0, Vinf;
queue<int> que;
vector<int> used(n + m);
rep(i, 0, n) if (match[i] == -1) {
used[i] = 1;
que.push(i);
}
while (!que.empty()) {
int v = que.front();
que.pop();
Vinf.push_back(v);
for (auto &to : G[v])
if (!used[to]) {
used[to] = 1;
que.push(to);
}
}
rep(i, 0, m) if (R[i] == -1) {
used[i + n] = 1;
que.push(i + n);
}
while (!que.empty()) {
int v = que.front();
que.pop();
V0.push_back(v);
for (auto &to : REV[v])
if (!used[to]) {
used[to] = 1;
que.push(to);
}
}
SCC scc(n + m);
rep(i, 0, n + m) for (auto &to : G[i]) {
if (!used[i] and !used[to])
scc.add_edge(i, to);
}
scc.run();
vector group(scc.m, vector<int>());
rep(i, 0, n + m) if (!used[i]) group[scc.id[i]].push_back(i);
vector<vector<int>> ret;
ret.push_back(V0);
for (auto &v : group)
ret.push_back(v);
ret.push_back(Vinf);
return ret;
}
/**
* @brief DM decomposition
*/
#line 4 "Algorithm/matroid.hpp"
struct GraphicMatroid {
using P = pair<int, int>;
int n;
vector<P> E;
vector<vector<P>> G;
vector<int> allow;
GraphicMatroid() {}
GraphicMatroid(int n, vector<P> &E) : n(n), E(E), G(n) {
rep(i, 0, SZ(E)) {
auto [u, v] = E[i];
G[u].push_back({v, i});
G[v].push_back({u, i});
}
}
void set(vector<int> &I) {
allow = I;
}
vector<int> circuit(int e) {
auto [x, y] = E[e];
vector<int> ret;
ret.push_back(e);
auto dfs = [&](auto &dfs, int v, int p) -> bool {
if (v == y)
return true;
for (auto &[to, i] : G[v])
if (allow[i] and to != p) {
ret.push_back(i);
if (dfs(dfs, to, v))
return true;
ret.pop_back();
}
return false;
};
if (dfs(dfs, x, -1))
return ret;
else
return {};
}
};
struct PartitionMatroid {
vector<int> grp; // -1:not assign
vector<int> lim;
vector<vector<int>> cnt;
PartitionMatroid() {}
PartitionMatroid(vector<int> &grp, vector<int> lim) : grp(grp), lim(lim) {}
void set(vector<int> &I) {
cnt.assign(SZ(lim), {});
rep(i, 0, SZ(grp)) if (I[i] != 0 and grp[i] != -1) {
cnt[grp[i]].push_back(i);
}
}
vector<int> circuit(int e) {
if (grp[e] == -1)
return {};
if (SZ(cnt[grp[e]]) + 1 > lim[grp[e]]) {
vector<int> ret = cnt[grp[e]];
ret.push_back(e);
return ret;
} else
return {};
}
};
struct TransversalMatroid {
using P = pair<int, int>;
int n, m;
vector<vector<int>> G, aug;
vector<int> match, fixed;
TransversalMatroid() {}
TransversalMatroid(int n, int m, vector<P> &E)
: n(n), m(m), G(n), aug(n), fixed(n) {
for (auto &[u, v] : E)
G[u].push_back(v);
}
void set(vector<int> &I) {
vector g(n, vector<int>());
rep(e, 0, n) if (I[e]) {
for (auto &to : G[e])
g[e].push_back(to);
}
auto match = BiMatching(n, m, g);
auto dm = DMdecomposition(n, m, g, match);
fixed.assign(m, 1);
for (auto &v : dm.front())
if (v >= n)
fixed[v - n] = 0;
aug.assign(n + m, {});
rep(e, 0, n) {
for (auto &to : G[e])
aug[e].push_back(to + n);
}
rep(i, 0, n) if (match[i] != -1) aug[match[i] + n].push_back(i);
}
vector<int> circuit(int e) {
for (auto &to : G[e])
if (!fixed[to]) {
return {};
}
vector<int> used(n + m);
queue<int> que;
used[e] = 1;
que.push(e);
while (!que.empty()) {
int v = que.front();
que.pop();
for (auto &to : aug[v])
if (!used[to]) {
used[to] = 1;
que.push(to);
}
}
vector<int> ret;
rep(i, 0, n) if (used[i]) ret.push_back(i);
return ret;
}
};
template <typename M1, typename M2, typename Val>
pair<vector<Val>, vector<vector<int>>>
MinimumMatroidIntersection(M1 &m1, M2 &m2, vector<Val> &ws) {
using P = pair<Val, int>;
int n = SZ(ws);
vector<Val> ret1;
vector<vector<int>> ret2;
Val profit = 0;
vector<int> I(n);
ret1.push_back(profit);
ret2.push_back(I);
for (;;) {
vector G(n + 2, vector<P>());
int S = n, T = n + 1;
m1.set(I);
m2.set(I);
rep(e, 0, n) {
if (I[e])
continue;
auto c1 = m1.circuit(e);
if (c1.empty())
G[S].push_back({ws[e] * (n + 1) + 1, e});
else
for (auto &to : c1)
if (e != to) {
G[to].push_back({ws[e] * (n + 1) + 1, e});
}
auto c2 = m2.circuit(e);
if (c2.empty())
G[e].push_back({Val(0), T});
else
for (auto &to : c2)
if (e != to) {
G[e].push_back({-ws[to] * (n + 1) + 1, to});
}
}
vector<ll> dist(n + 2, INF);
dist[S] = 0;
vector<int> pre(n + 2, -1);
for (;;) {
bool upd = 0;
rep(v, 0, n + 2) if (dist[v] != INF) {
for (auto &[cost, to] : G[v]) {
if (chmin(dist[to], dist[v] + cost)) {
pre[to] = v;
upd = 1;
}
}
}
if (!upd)
break;
}
if (dist[T] == INF)
break;
int cur = T;
while (cur != S) {
cur = pre[cur];
if (cur != S) {
I[cur] ^= 1;
if (I[cur])
profit += ws[cur];
else
profit -= ws[cur];
}
}
ret1.push_back(profit);
ret2.push_back(I);
}
return {ret1, ret2};
}
/**
* @brief Matroid
*/