#include <cmath>
#include <functional>
#include <random>
#include <vector>
#include <limits>

#include "randutils/randutils.hpp"
#include "pcg-cpp/include/pcg_random.hpp"

#include "blossom5-v2.05.src/PerfectMatching.h"

using Rng = randutils::random_generator<pcg32>;


class Graph {
public:
class Edge;
class HalfEdge;

class Coordinate {
  std::array<int, 2> r;
  public:
  void operator=(const std::array<int, 2>& rp) {
    r = rp;
  }
  void operator+=(const Coordinate& x) {
    r[0] += x.r[0];
    r[1] += x.r[1];
  }
  void operator-=(const Coordinate& x) {
    r[0] -= x.r[0];
    r[1] -= x.r[1];
  }
  int operator[](unsigned i) const {
    return r[i];
  }
};

class Vertex {
public:
  unsigned index;
  std::vector<std::reference_wrapper<HalfEdge>> neighbors;
  Coordinate coordinate;

  void addEdge(HalfEdge& e) {
    neighbors.push_back(e);
  }
};

class HalfEdge {
private:
  Vertex* tail;
  Vertex* head;

public:
  const Edge& edge;

  HalfEdge(const Edge& e) : edge(e) {}
  void setVertices(Vertex& vt, Vertex& vh) {
    tail = &vt;
    head = &vh;
  }
  const Vertex& getHead() const {
    return *head;
  }
  const Vertex& getTail() const {
    return *tail;
  }
};

class Edge {
public:
  std::array<HalfEdge, 2> halfedges;
  double weight;

  Edge() : halfedges{*this, *this} {};
  void setVertices(Vertex& red, Vertex& blue) {
    halfedges[0].setVertices(red, blue);
    halfedges[1].setVertices(blue, red);
    red.addEdge(halfedges[0]);
    blue.addEdge(halfedges[1]);
  }
};
  std::vector<Vertex> vertices;
  std::vector<Edge> edges;

  Graph(int n, Rng& r) : vertices(2 * n * (n + 1)), edges(pow(2 * n, 2)) {
    unsigned M = vertices.size() / 2;
    for (int i = 0; i < M; i++) {
      vertices[i].index = i;
      vertices[i].coordinate = {2 * (i % (n + 1)), 2 * (i / (n + 1)) + 1};
      vertices[M + i].index = M + i;
      vertices[M + i].coordinate = {2 * (i % n) + 1, 2 * (i / n)};
    }
    for (unsigned i = 0; i < edges.size(); i++) {
      Vertex& redVertex = vertices[(1 + (i % (2 * n))) / 2 + (n + 1) * ((i / 4) / n)];
      Vertex& blueVertex = vertices[M + (i % (2 * n)) / 2 + n * (((i + 2 * n) / 4) / n)];
      edges[i].setVertices(redVertex, blueVertex);
      edges[i].weight = r.variate<double, std::exponential_distribution>(1);
    }
  }
};

bool edgeMatched(PerfectMatching& pm, const Graph::Edge& e) {
  return e.halfedges[0].getTail().index == pm.GetMatch(e.halfedges[0].getHead().index);
}