#include #include #include #include class Edge { public: int u, v; bool operator==(const Edge &e) const { return u == e.u && v == e.v; } }; // Its assumed that u,v are distinct. // Check if all four vertices are distinct bool edgeConflicts(const Edge &e1, const Edge &e2) { return (e1.u == e2.u || e1.u == e2.v || e1.v == e2.u || e1.v == e2.v); } std::ostream &operator<<(std::ostream &s, const Edge &e) { s << '{' << e.u << ',' << e.v << '}'; return s; } class DiDegree { public: int in; int out; }; typedef std::vector DegreeSequence; typedef std::vector DiDegreeSequence; class Graph { public: Graph() : edgeCount(0) {} Graph(int n) : edgeCount(0) { adj.resize(n); } ~Graph() {} bool createFromSequence(const DegreeSequence &d) { edgeCount = std::accumulate(d.begin(), d.end(), 0); // // TODO // return false; } bool hasEdge(const Edge &e) const { for (int v : adj[e.u]) { if (v == e.v) return true; } return false; } // There are two possible edge exchanges // switchType indicates which one is desired // Returns false if the switch is not possible bool exchangeEdges(const Edge &e1, const Edge &e2, bool switchType) { // The new edges configuration is one of these two // A) e1.u - e2.u and e1.v - e2.v // B) e1.u - e2.v and e1.v - e2.u // First check if the move is possible if (switchType) { if (hasEdge({e1.u, e2.u}) || hasEdge({e1.v, e2.v})) return false; // conflicting edges } else { if (hasEdge({e1.u, e2.v}) || hasEdge({e1.v, e2.u})) return false; // conflicting edges } // Find the edges in the adjacency lists int i1, j1, i2, j2; for (i1 = 0; i1 < (int)adj[e1.u].size(); ++i1) { if (adj[e1.u][i1] == e1.v) break; } for (j1 = 0; j1 < (int)adj[e1.v].size(); ++j1) { if (adj[e1.v][j1] == e1.u) break; } for (i2 = 0; i2 < (int)adj[e2.u].size(); ++i2) { if (adj[e2.u][i2] == e2.v) break; } for (j2 = 0; j2 < (int)adj[e2.v].size(); ++j2) { if (adj[e2.v][j2] == e2.u) break; } // Remove the old edges bool removedOne = false; for (auto iter = edges.begin(); iter != edges.end();) { if (*iter == e1) { iter = edges.erase(iter); if (removedOne) break; removedOne = true; } else if (*iter == e2) { iter = edges.erase(iter); if (removedOne) break; removedOne = true; } else { ++iter; } } // Add the new edges if (switchType) { adj[e1.u][i1] = e2.u; adj[e1.v][j1] = e2.v; adj[e2.u][i2] = e1.u; adj[e2.v][j2] = e1.v; edges.push_back({e1.u, e2.u}); edges.push_back({e1.v, e2.v}); } else { adj[e1.u][i1] = e2.v; adj[e1.v][j1] = e2.u; adj[e2.u][i2] = e1.v; adj[e2.v][j2] = e1.u; edges.push_back({e1.u, e2.v}); edges.push_back({e1.v, e2.u}); } return true; } // Graph is saved in two formats for speed // The two should be kept consistent at all times std::vector> adj; std::vector edges; int edgeCount; }; class SwitchChain { public: SwitchChain() : permutationDistribution(0, 2) { // random_device uses hardware entropy if available std::random_device rd; mt.seed(rd()); } ~SwitchChain() {} bool initialize(const DegreeSequence &d) { if (!g.createFromSequence(d)) return false; edgeDistribution.param( std::uniform_int_distribution<>::param_type(0, g.edgeCount - 1)); return true; } bool doMove() { Edge e1 = g.edges[edgeDistribution(mt)]; Edge e2 = g.edges[edgeDistribution(mt)]; // Keep regenerating while conflicting edges int timeout = 0; while (edgeConflicts(e1, e2)) { e1 = g.edges[edgeDistribution(mt)]; e2 = g.edges[edgeDistribution(mt)]; ++timeout; if (timeout % 100 == 0) { std::cerr << "Warning: sampled " << timeout << " random edges but they keep conflicting.\n"; } } // Consider one of the three possible permutations // 1) e1.u - e1.v and e2.u - e2.v (original) // 2) e1.u - e2.u and e1.v - e2.v // 3) e1.u - e2.v and e1.v - e2.u // Note that it might be that these new edges already exist // in which case we also reject the move // This is checked in exchangeEdges int perm = permutationDistribution(mt); if (perm == 0) // Original permutation return false; return g.exchangeEdges(e1, e2, perm == 1); } Graph g; std::mt19937 mt; std::uniform_int_distribution<> edgeDistribution; std::uniform_int_distribution<> permutationDistribution; }; int main() { SwitchChain chain; if (!chain.initialize({3, 2, 4, 2, 2, 1})) { std::cerr << "Could not initialize Markov chain.\n"; return 1; } std::cout << "Starting switch Markov chain" << std::endl; int movesDone = 0; for (int i = 0; i < 100; ++i) if (chain.doMove()) ++movesDone; std::cout << movesDone << "/100 moves succeeded." << std::endl; return 0; }