SAT mapper fix (#28)

demos/qiskit_patterns/demo_src/graph.py

@@ -2,20 +2,22 @@
 import numpy as np
 import matplotlib.pyplot as plt
 
+
 def generate_demo_graph():
     # Generating a graph of 5 nodes
     n = 5  # Number of nodes in graph
     G = nx.Graph()
     G.add_nodes_from(np.arange(0, n, 1))
-    elist = [(0, 1, 1.0), (0, 2, 1.0), (0, 4, 1.0), (1, 2, 1.0), (2, 3, 1.0),(3, 4, 1.0)]
+    elist = [(0, 1, 1.0), (0, 2, 1.0), (0, 4, 1.0), (1, 2, 1.0), (2, 3, 1.0), (3, 4, 1.0)]
     # tuple is (i,j,weight) where (i,j) is the edge
     G.add_weighted_edges_from(elist)
     return G
 
+
 def draw_graph(G):
-    colors = ['tab:grey' for node in G.nodes()]
+    colors = ["tab:grey" for node in G.nodes()]
     pos = nx.spring_layout(G)
     default_axes = plt.axes(frameon=True)
     nx.draw_networkx(G, node_color=colors, node_size=600, alpha=0.8, ax=default_axes, pos=pos)
     edge_labels = nx.get_edge_attributes(G, "weight")
-    nx.draw_networkx_edge_labels(G, pos=pos, edge_labels=edge_labels)
+    nx.draw_networkx_edge_labels(G, pos=pos, edge_labels=edge_labels)

demos/qiskit_patterns/demo_src/map.py

@@ -1,15 +1,18 @@
 from qiskit_optimization.applications import Maxcut
 from qiskit.circuit.library import QAOAAnsatz
 
+
 def map_graph_to_qubo(graph):
     max_cut = Maxcut(graph)
     qp = max_cut.to_quadratic_program()
     return qp
 
+
 def map_qubo_to_ising(qubo):
     qubitOp, offset = qubo.to_ising()
     return qubitOp, offset
 
+
 def map_ising_to_circuit(hamiltonian, num_layers):
     ansatz = QAOAAnsatz(cost_operator=hamiltonian, reps=num_layers)
     ansatz.measure_all()

demos/qiskit_patterns/demo_src/post.py

@@ -23,30 +23,31 @@ def sample_most_likely(state_vector, num_bits):
     most_likely_bitstring.reverse()
     return np.asarray(most_likely_bitstring)
 
+
 # auxiliary function to plot graphs
 def plot_result(G, x):
-    colors = ['tab:grey' if i == 0 else 'tab:purple' for i in x]
+    colors = ["tab:grey" if i == 0 else "tab:purple" for i in x]
     pos, default_axes = nx.spring_layout(G), plt.axes(frameon=True)
     nx.draw_networkx(G, node_color=colors, node_size=600, alpha=0.8, pos=pos)
 
 
 def plot_distribution(final_distribution):
-    matplotlib.rcParams.update({'font.size': 10})
+    matplotlib.rcParams.update({"font.size": 10})
     final_bits = final_distribution.binary_probabilities()
     values = np.abs(list(final_bits.values()))
-    top_4_values = sorted(values, reverse = True)[:4]
+    top_4_values = sorted(values, reverse=True)[:4]
     positions = []
     for value in top_4_values:
         positions.append(np.where(values == value)[0])
-    fig = plt.figure(figsize = (11,6))
-    ax=fig.add_subplot(1,1,1)
+    fig = plt.figure(figsize=(11, 6))
+    ax = fig.add_subplot(1, 1, 1)
     plt.xticks(rotation=45)
     plt.title("Result Distribution")
     plt.xlabel("Bitstrings (reversed)")
     plt.ylabel("Probability")
-    ax.bar(list(final_bits.keys()), list(final_bits.values()), color='tab:grey')
+    ax.bar(list(final_bits.keys()), list(final_bits.values()), color="tab:grey")
     for p in positions:
-        ax.get_children()[int(p)].set_color('tab:purple')
+        ax.get_children()[int(p)].set_color("tab:purple")
     plt.show()
 
 
@@ -55,7 +56,6 @@ def samples_to_objective_values(samples, qp):
     """Convert the samples to values of the objective function."""
     objective_values = defaultdict(float)
     for bit_str, prob in samples.items():
-
         # Qiskit use little endian hence the [::-1]
         candidate_sol = [int(bit) for bit in bit_str[::-1]]
         fval = qp.objective.evaluate(candidate_sol)
@@ -88,7 +88,6 @@ def load_data(qp):
 
 # auxiliary function to load the QP from the saved Paulis
 def load_qp():
-
     # First, load the Paulis that encode the MaxCut problem.
     with open("data/125node_example_ising.txt", "r") as fin:
         paulis, lines = [], fin.read()
@@ -102,7 +101,10 @@ def load_qp():
     # Next, convert the Paulis to a weighted graph.
     wedges = []
     for pauli_str, coefficient in paulis:
-        wedges.append([idx for idx, char in enumerate(pauli_str[::-1]) if char == "Z"] + [{"weight": coefficient}])
+        wedges.append(
+            [idx for idx, char in enumerate(pauli_str[::-1]) if char == "Z"]
+            + [{"weight": coefficient}]
+        )
 
     weighted_graph = nx.DiGraph(wedges)
 
@@ -136,7 +138,8 @@ def _plot_cdf(objective_values: dict, ax, label):
     x_vals = sorted(objective_values.keys())
     y_vals = np.cumsum([objective_values[x] for x in x_vals])
     ax.plot(x_vals, y_vals, label=label)
-
+
+
 def plot_cdf(depth_zero, depth_one, max_cut, min_cut, ax, title):
     _plot_cdf(depth_zero, ax, "Depth-zero")
     _plot_cdf(depth_one, ax, "Depth-one")
@@ -146,4 +149,4 @@ def plot_cdf(depth_zero, depth_one, max_cut, min_cut, ax, title):
     ax.set_title(title + f" {approx}%")
     ax.set_xlabel("Objective function value")
     ax.set_ylabel("Cumulative distribution function")
-    ax.legend(loc=2);
+    ax.legend(loc=2)

demos/qiskit_patterns/demo_src/transpile.py

@@ -1,4 +1,5 @@
 from qiskit import transpile
 
+
 def optimize_circuit(circuit, backend):
-    return transpile(circuit, backend=backend)
+    return transpile(circuit, backend=backend)

qopt_best_practices/sat_mapping/sat_mapper.py

@@ -119,8 +119,12 @@ def interrupt(solver):
         # number of swap layers that satisfies the subgraph isomorphism problem.
         while min_layers < max_layers:
             num_layers = (min_layers + max_layers) // 2
-            distance_matrix = swap_strategy.distance_matrix
-            connectivity_matrix = (distance_matrix <= num_layers).astype(int)
+
+            # Create the connectivity matrix. Note that if the swap strategy cannot reach
+            # full connectivity then its distance matrix will have entries with -1. These
+            # entries must be treated as False.
+            d_matrix = swap_strategy.distance_matrix
+            connectivity_matrix = ((-1 < d_matrix) & (d_matrix <= num_layers)).astype(int)
             # Make a cnf for the adjacency constraint
             cnf2 = []
             for e_0, e_1 in program_graph.edges:

test/test_sat_mapping.py

@@ -5,6 +5,7 @@
 import os
 import networkx as nx
 
+from qiskit.transpiler import CouplingMap
 from qiskit.transpiler.passes.routing.commuting_2q_gate_routing import SwapStrategy
 
 from qopt_best_practices.utils import build_max_cut_graph, build_max_cut_paulis
@@ -59,3 +60,36 @@ def test_remap_graph_with_sat(self):
         )
 
         self.assertTrue(nx.is_isomorphic(remapped_g, self.mapped_graph))
+
+    def test_deficient_strategy(self):
+        """Test that the SAT mapper works when the SWAP strategy is deficient.
+
+        Note: a deficient strategy does not result in full connectivity but
+        may still be useful.
+        """
+        cmap = CouplingMap([(idx, idx + 1) for idx in range(10)])
+
+        # This swap strategy is deficient but can route the graph below.
+        swaps = (
+            ((0, 1), (2, 3), (4, 5), (6, 7), (8, 9)),
+            ((1, 2), (3, 4), (5, 6), (7, 8)),
+            ((2, 3), (4, 5), (6, 7), (8, 9)),
+            (),
+            (),
+            (),
+            (),
+            (),
+        )
+        swap_strategy = SwapStrategy(cmap, swaps)
+        graph = nx.random_regular_graph(3, 10, seed=2)
+
+        mapper = SATMapper()
+
+        _, permutation, min_layer = mapper.remap_graph_with_sat(graph, swap_strategy)
+
+        # Spot check a few permutations.
+        self.assertEqual(permutation[0], 9)
+        self.assertEqual(permutation[8], 1)
+
+        # Crucially, if the `connectivity_matrix` in `find_initial_mappings` we get a wrong result.
+        self.assertEqual(min_layer, 3)