Merge pull request #6 from oscarhiggott/doxygen-docstrings

C++ docstrings added

src/pymatching/bindings.cpp

@@ -30,7 +30,7 @@ PYBIND11_MODULE(_cpp_mwpm, m) {
      .def("space_time_shortest_path", &IStabiliserGraph::SpaceTimeShortestPath, "node1"_a, "node2"_a)
      .def("qubit_ids", &IStabiliserGraph::QubitIDs, "node1"_a, "node2"_a)
      .def("get_num_qubits", &IStabiliserGraph::GetNumQubits)
-     .def("get_num_stabilisers", &IStabiliserGraph::GetNumStabilisers)
+     .def("get_num_nodes", &IStabiliserGraph::GetNumNodes)
      .def("get_num_edges", &IStabiliserGraph::GetNumEdges)
      .def("compute_all_pairs_shortest_paths", &IStabiliserGraph::ComputeAllPairsShortestPaths)
      .def("has_computed_all_pairs_shortest_paths", &IStabiliserGraph::HasComputedAllPairsShortestPaths)

src/pymatching/lemon_mwpm.h

@@ -17,6 +17,32 @@
 #include <pybind11/pybind11.h>
 #include <pybind11/numpy.h>
 
-
+/**
+ * @brief Given a stabiliser graph sg and a vector `defects` of indices of nodes that have a -1 syndrome, 
+ * find the find the minimum weight perfect matching in the complete graph with nodes in the defects 
+ * list, and where the edge between node i and j is given by the distance between i and j in sg. The 
+ * distances and shortest paths between nodes in the stabiliser graph sg are all precomputed and this 
+ * method returns the exact minimum-weight perfect matching. As a result it is suitable for matching graphs 
+ * with a few thousand nodes or less, but will be very memory and compute intensive for larger matching graphs.
+ * Returns a noise vector N for which N[i]=1 if qubit_id appeared an odd number of times in the minimum weight 
+ * perfect matching and N[i]=0 otherwise.
+ * 
+ * @param sg A stabiliser graph
+ * @param defects The indices of nodes that are associated with a -1 syndrome
+ * @return py::array_t<std::uint8_t> The noise vector for the minimum-weight perfect matching
+ */
 py::array_t<std::uint8_t> LemonDecode(IStabiliserGraph& sg, const py::array_t<int>& defects);
+/**
+ * @brief Given a stabiliser graph `sg`, a vector `defects` of indices of nodes that have a -1 syndrome and 
+ * a chosen `num_neighbours`, find the minimum weight perfect matching in the graph V where each defect node 
+ * is connected by an edge to each of the `num_neighbours` nearest other defect nodes in sg, and where the 
+ * weight of each edge is the distance between the two defect nodes in `sg`.
+ * Returns a noise vector N for which N[i]=1 if qubit_id appeared an odd number of times in the minimum weight 
+ * perfect matching and N[i]=0 otherwise.
+ * 
+ * @param sg A stabiliser graph
+ * @param defects The indices of nodes that are associated with a -1 syndrome
+ * @param num_neighbours The number of closest defects to connect each defect to in the matching graph
+ * @return py::array_t<std::uint8_t> The noise vector for the minimum-weight perfect matching
+ */
 py::array_t<std::uint8_t> LemonDecodeMatchNeighbourhood(WeightedStabiliserGraph& sg, const py::array_t<int>& defects, int num_neighbours);

src/pymatching/matching.py

@@ -270,7 +270,7 @@ def decode(self, z, num_neighbours=20):
                 <= self.num_stabilisers+len(self.boundary)):
             defects = z.nonzero()[0]
         elif len(z.shape) == 2 and z.shape[0] == self.num_stabilisers:
-            num_stabs = self.stabiliser_graph.get_num_stabilisers()
+            num_stabs = self.stabiliser_graph.get_num_nodes()
             max_num_defects = z.shape[0]*z.shape[1]
             if max_num_defects > num_stabs:
                 raise ValueError(f"Syndrome size {z.shape[0]}x{z.shape[1]} exceeds" \
@@ -305,8 +305,8 @@ def add_noise(self):
             Noise vector (binary numpy int array of length self.num_qubits)
         numpy.ndarray of dtype int
             Syndrome vector (binary numpy int array of length 
-            self.num_stabilisers if there is no boundary, or self.num_stabilisers+1
-            if there is a boundary)
+            self.num_stabilisers if there is no boundary, or self.num_stabilisers+len(self.boundary)
+            if there are boundary nodes)
         """
         if not self.stabiliser_graph.all_edges_have_error_probabilities:
             return None

src/pymatching/rand_gen.h

@@ -20,6 +20,18 @@ std::mt19937 & global_urng();
 
 void randomize();
 
+/**
+ * @brief Set the seed for the mt19937 random number generator
+ * 
+ * @param s 
+ */
 void set_seed(unsigned s);
 
+/**
+ * @brief A random double chosen uniformly at random between `from` and `to`
+ * 
+ * @param from
+ * @param to 
+ * @return double 
+ */
 double rand_float(double from, double to);

src/pymatching/stabiliser_graph.cpp

@@ -17,7 +17,7 @@
 
 
 double IStabiliserGraph::SpaceTimeDistance(int node1, int node2) {
-    int num_stab = GetNumStabilisers();
+    int num_stab = GetNumNodes();
     if ((node1 < num_stab) && (node2 < num_stab)){
         return Distance(node1, node2);
     }
@@ -29,7 +29,7 @@ double IStabiliserGraph::SpaceTimeDistance(int node1, int node2) {
 }
 
 std::vector<int> IStabiliserGraph::SpaceTimeShortestPath(int node1, int node2) {
-    int num_stab = GetNumStabilisers();
+    int num_stab = GetNumNodes();
     int r1 = node1 % num_stab;
     int r2 = node2 % num_stab;
     return ShortestPath(r1, r2);

src/pymatching/stabiliser_graph.h

@@ -26,7 +26,7 @@ class IStabiliserGraph{
         virtual std::set<int> QubitIDs(int node1, int node2) const = 0;
         virtual int GetNumEdges() const = 0;
         virtual int GetNumQubits() const = 0;
-        virtual int GetNumStabilisers() const = 0;
+        virtual int GetNumNodes() const = 0;
         virtual int GetNumConnectedComponents() const = 0;
         virtual std::vector<int> GetBoundary() const = 0;
         virtual void SetBoundary(std::vector<int>& boundary) = 0;

src/pymatching/weighted_stabiliser_graph.cpp

@@ -284,7 +284,7 @@ int WeightedStabiliserGraph::GetNumQubits() const {
     return num_qubits;
 }
 
-int WeightedStabiliserGraph::GetNumStabilisers() const {
+int WeightedStabiliserGraph::GetNumNodes() const {
     return boost::num_vertices(stabiliser_graph);
 };
 
@@ -299,7 +299,7 @@ std::set<int> WeightedStabiliserGraph::QubitIDs(int node1, int node2) const {
 }
 
 std::pair<py::array_t<std::uint8_t>,py::array_t<std::uint8_t>> WeightedStabiliserGraph::AddNoise() const {
-    auto syndrome = new std::vector<int>(GetNumStabilisers(), 0);
+    auto syndrome = new std::vector<int>(GetNumNodes(), 0);
     auto error = new std::vector<int>(GetNumQubits(), 0);
     double p;
     std::set<int> qids;

src/pymatching/weighted_stabiliser_graph.h

@@ -39,47 +39,194 @@ typedef boost::graph_traits < wgraph_t >::vertex_descriptor vertex_descriptor;
 typedef boost::graph_traits < wgraph_t >::edge_descriptor edge_descriptor;
 
 
+/**
+ * @brief 
+ * 
+ */
 class WeightedStabiliserGraph : public IStabiliserGraph{
     public:
         wgraph_t stabiliser_graph;
+        /**
+         * @brief Construct a new Weighted Stabiliser Graph object
+         * 
+         * @param num_stabilisers Number of stabiliser nodes (this excludes boundary nodes)
+         * @param boundary Indices of the boundary nodes
+         */
         WeightedStabiliserGraph(
             int num_stabilisers,
             std::vector<int>& boundary
             );
-        WeightedStabiliserGraph(
-            const py::array_t<int>& indices, 
-            const py::array_t<double>& weights,
-            std::vector<int>& boundary
-        );
-        WeightedStabiliserGraph(
-            const py::array_t<int>& indices, 
-            const py::array_t<double>& weights,
-            const py::array_t<double>& error_probabilies,
-            std::vector<int>& boundary
-        );
+        /**
+         * @brief Add an edge to the Weighted Stabiliser Graph object
+         * 
+         * @param node1 Index of the first node
+         * @param node2 Index of the second node
+         * @param qubit_ids Indices of the qubits associated with this edge
+         * @param weight Weight of the edge
+         * @param error_probability The error probability associated with an edge (optional, set to -1 if not needed)
+         * @param has_error_probability Flag whether a valid error probability has been supplied
+         */
         void AddEdge(int node1, int node2, std::set<int> qubit_ids,
                      double weight, double error_probability,
                      bool has_error_probability);
+        /**
+         * @brief Compute and store the shortest path between every pair of nodes in the matching graph.
+         * This is only used or needed if exact matching is used (rather than the default local matching).
+         * Note that this is only useful for relatively small matching graphs, and is very memory and 
+         * compute intensive matching graphs with many thousands of nodes.
+         * 
+         */
         virtual void ComputeAllPairsShortestPaths();
+        /**
+         * @brief Distance between two nodes in the matching graph. This method is used only for exact matching,
+         * since it uses the pre-computed all-pairs shortest paths data, and computes the all-pairs shortest paths
+         * if this has not yet been done.
+         * 
+         * @param node1 The index of the first node
+         * @param node2 The index of the second node
+         * @return double The distance between the first node and the second node in the matching graph
+         */
         virtual double Distance(int node1, int node2);
+        /**
+         * @brief Shortest path between two nodes in the matching graph. This method is used only for exact matching,
+         * since it uses the pre-computed all-pairs shortest paths data, and computes the all-pairs shortest paths
+         * if this has not yet been done.
+         * 
+         * @param node1 The index of the first node
+         * @param node2 The index of the second node
+         * @return std::vector<int> The vertex indices along the path between the first node and the second node in the matching graph
+         */
         virtual std::vector<int> ShortestPath(int node1, int node2);
+        /**
+         * @brief Get the qubit ids associated with the edge between node1 and node2
+         * 
+         * @param node1 The index of the first node
+         * @param node2 The index of the second node
+         * @return std::set<int> The set of qubit ids associated with the edge between node1 and node2
+         */
         virtual std::set<int> QubitIDs(int node1, int node2) const;
+        /**
+         * @brief Get the number of qubits associated with edges of the matching graph
+         * 
+         * @return int The number of qubits
+         */
         virtual int GetNumQubits() const;
-        virtual int GetNumStabilisers() const;
+        /**
+         * @brief Get the number of nodes in the stabiliser graph (this includes both boundaries and stabilisers)
+         * 
+         * @return int Number of nodes in stabiliser_graph
+         */
+        virtual int GetNumNodes() const;
+        /**
+         * @brief Get the number of edges in the stabiliser graph
+         * 
+         * @return int Number of edges in stabiliser_graph
+         */
         virtual int GetNumEdges() const;
+        /**
+         * @brief If an error_probability is assigned to every edge, flip each edge 
+         * with its corresponding error_probability. If an edge is flipped, add (mod 2) 
+         * an error to the associated qubits (specified by the qubit_ids edge data).
+         * The qubit errors are returned as a binary numpy array, as well as a syndrome 
+         * vector, also as a binary numpy array. The length of the syndrome vector is 
+         * is the number of nodes in the stabiliser graph (there is an element for each 
+         * stabiliser as well as for each boundary node). The syndromes of the boundary 
+         * nodes are all set to zero unless the parity of the stabiliser syndromes is odd, 
+         * in which case the first boundary node's syndrome is flipped.
+         * 
+         * @return std::pair<py::array_t<std::uint8_t>,py::array_t<std::uint8_t>> The first item in 
+         * the pair is the noise vector, and the second item is the syndrome vector.
+         */
         std::pair<py::array_t<std::uint8_t>,py::array_t<std::uint8_t>> AddNoise() const;
+        /**
+         * @brief The distance between every pair of nodes in the stabiliser graph, if 
+         * ComputeAllPairsShortestPaths has been run. all_distances[i][j] is the distance 
+         * between node i and node j in the stabiliser graph. Note that this is only used 
+         * if exact matching is used (the function LemonDecode in lemon_mwpm.cpp), and not if 
+         * (the Python default) local matching is used (the function LemonDecodeMatchNeighbourhood 
+         * in lemon_mwpm.cpp).
+         * 
+         */
         std::vector<std::vector<double>> all_distances;
+        /**
+         * @brief all_predecessors[i] is the PredecessorMap found by Boost's dijkstra_shortest_paths 
+         * from node i in stabiliser_graph to all other nodes. In other words, all_predecessors[i][j] 
+         * is the parent of node j in the shortest path from node i to node j in stabiliser_graph.
+         * Note that this is only used if exact matching is used (the function LemonDecode in lemon_mwpm.cpp), 
+         * and not if (the Python default) local matching is used (the function LemonDecodeMatchNeighbourhood 
+         * in lemon_mwpm.cpp).
+         * 
+         */
         std::vector<std::vector<vertex_descriptor>> all_predecessors;
-        bool all_edges_have_error_probabilities;
+        bool all_edges_have_error_probabilities; // true if every edge in stabiliser_graph has an error_probability in its edge data
+        /**
+         * @brief Get the indices of the boundary nodes
+         * 
+         * @return std::vector<int> The indices of the boundary nodes
+         */
         virtual std::vector<int> GetBoundary() const;
+        /**
+         * @brief Set the indices of the boundary nodes
+         * 
+         * @param boundary The indices of the boundary nodes
+         */
         virtual void SetBoundary(std::vector<int>& boundary);
+        /**
+         * @brief Flag whether or not the all-pairs shortest paths have been computed. 
+         * The all-pairs shortest paths are only needed for exact matching, which is not 
+         * a default option in the Python bindings.
+         * 
+         * @return true 
+         * @return false 
+         */
         virtual bool HasComputedAllPairsShortestPaths() const;
+        /**
+         * @brief Reset the _predecessors and _distances attributes, which 
+         * are used by WeightedStabiliserGraph::GetPath and WeightedStabiliserGraph::GetNearestNeighbours. 
+         * This just preallocates both of these vectors appropriately for use by the local dijkstra search 
+         * and will only do so if these vectors are not already the correct size. Once these attributes 
+         * are correctly preallocated the first time, the GetPath and GetNearestNeighbours methods always reset 
+         * only the elements in these vectors that have been used for efficiency (since these methods only use 
+         * a local search).
+         * 
+         */
         void ResetDijkstraNeighbours();
+        /**
+         * @brief Get the num_neighbours nearest neighbours i of a source node in the stabiliser_graph for which 
+         * defect_id[i] > -1. This is a local Dijkstra search that halts once num_neighbours nodes i have been found 
+         * that satisfy defect_id[i] > -1. The function returns a vector of pairs, where the first item in each 
+         * pair is the distance from source to one of the nearest nodes i, and the second item in the pair 
+         * is defect_id[i].
+         * This is used by LemonDecodeMatchNeighbourhood to find the num_neigbours nearest defects to be included 
+         * in the matching graph.
+         * 
+         * @param source 
+         * @param num_neighbours 
+         * @param defect_id 
+         * @return std::vector<std::pair<int, double>> 
+         */
         std::vector<std::pair<int, double>> GetNearestNeighbours(
             int source, int num_neighbours, std::vector<int>& defect_id);
+        /**
+         * @brief Get the shortest path from source to target using a version of Dijkstra 
+         * that terminates once target is found.
+         * 
+         * @param source 
+         * @param target 
+         * @return std::vector<int> 
+         */
         std::vector<int> GetPath(
             int source, int target);
+        /**
+         * @brief Get the number of connected components in the stabiliser graph
+         * 
+         * @return int The number of components
+         */
         virtual int GetNumConnectedComponents() const;
+        /**
+         * @brief The indices of the boundary nodes
+         * 
+         */
         std::vector<int> boundary;
         std::vector<double> _distances;
         std::vector<int> _predecessors;

tests/test_weighted_stabiliser_graph.py

@@ -35,7 +35,7 @@ def test_weighted_spacetime_shortest_path():
     assert(w.shortest_path(3, 5) == [3, 2, 5])
     assert(w.shortest_path(4, 2) == [4, 5, 2])
     assert(w.get_num_qubits() == 9)
-    assert(w.get_num_stabilisers() == 6)
+    assert(w.get_num_nodes() == 6)
 
 
 def test_weighted_num_qubits_and_stabilisers():
@@ -51,7 +51,7 @@ def test_weighted_num_qubits_and_stabilisers():
     w.add_edge(4, 5, {6}, 9.0)
     w.compute_all_pairs_shortest_paths()
     assert(w.get_num_qubits() == 7)
-    assert(w.get_num_stabilisers() == 6)
+    assert(w.get_num_nodes() == 6)
 
 
 @pytest.mark.parametrize("num_loops", range(1, 10))