Issue 344, Issue 436 and Issue 431 (#437)

* Updated location of QuantumErrorChannel in aer

* linting and Updates for qc.data accesses

* Updated requirements-dev.txt wrt black

* linter change

* Fix for issue #431

* Linting

* Adjusting one int type to long int + corrected doc

requirements-dev.txt

@@ -6,4 +6,4 @@ sphinx-autodoc-typehints
 nbsphinx
 ddt~=1.4.2
 matplotlib>=3.3.0
-black[jupyter]~=22.1
+black[jupyter]

src/qiskit_qec/analysis/distance.h

@@ -6,6 +6,8 @@
 #include <set>
 #include <exception>
 #include <iterator>
+#include <iostream>
+#include <cassert>
 
 #include "linear.h"
 #include "combinations.h"

src/qiskit_qec/analysis/distance.py

@@ -1,4 +1,5 @@
 """Compute code distance."""
+
 from typing import List
 from itertools import combinations, product
 import functools
@@ -28,6 +29,8 @@ def _distance_test(stab: np.ndarray, logic_op: np.ndarray, weight: int) -> bool:
     Pauli operators on n qubits are represented as integer numpy arrays
     of length 2n in the order [X-part, Z-part].
 
+    Restrictions: n < 63 to avoid overflows in numpy
+
     Returns True or False.
     """
     # number of qubits
@@ -94,6 +97,8 @@ def _minimum_distance_2_python(stabilizer: np.ndarray, gauge: np.ndarray, max_we
 
     Second method based on parititioning errors.
 
+    Restrictions: n < 63 to avoid overflows in numpy
+
     Returns the minimum distance of the code, or 0 if greater than max_weight.
     """
     if symplectic.is_stabilizer_group(gauge):
@@ -109,17 +114,21 @@ def _minimum_distance_2_python(stabilizer: np.ndarray, gauge: np.ndarray, max_we
         zl = np.setdiff1d(zp_rows, z_rows).view(zp.dtype).reshape(-1, zp.shape[1])
     if xl.shape[0] == 0:  # k = 0, fall back to first method
         return _minimum_distance_1_python(stabilizer, gauge, max_weight)
+
     weight = max_weight + 1
+
     for row in range(xl.shape[0]):
         for w in range(1, max_weight + 1):
-            if _distance_test(stabilizer.astype(int), xl[row].astype(int), w):
+            if _distance_test(stabilizer.astype(int), xl[row].astype(int), w) is True:
                 weight = min(weight, w)
                 break
+
     for row in range(zl.shape[0]):
         for w in range(1, max_weight + 1):
-            if _distance_test(stabilizer.astype(int), zl[row].astype(int), w):
+            if _distance_test(stabilizer.astype(int), zl[row].astype(int), w) is True:
                 weight = min(weight, w)
                 break
+
     if weight < max_weight + 1:
         return weight
     else:
@@ -223,7 +232,7 @@ def _minimum_distance_2_compiled(stabilizer: np.ndarray, gauge: np.ndarray, max_
 def minimum_distance(
     stabilizer_or_gauge: np.ndarray,
     max_weight: int = 10,
-    method: str = "enumerate",
+    method: str = "partition",
     try_compiled: bool = True,
 ) -> int:
     """Minimum distance of (subsystem) stabilizer code.
@@ -233,6 +242,12 @@ def minimum_distance(
 
     method and try_compiled select an algorithm and implementation.
 
+    Restrictions:
+
+    Partition: Current methods require n-k < 63
+
+    Enumerate: None, although very slow for larger n
+
     Returns the minimum distance of the code, or 0 if greater than max_weight.
     """
     method_enumerate: str = "enumerate"
@@ -253,11 +268,15 @@ def minimum_distance(
         if method == method_enumerate:
             distance = _c_minimum_distance(inputform1, inputform2, max_weight)
         elif method == method_partition:
+            if method == method_partition and not stabilizer.shape[0] < 63:
+                raise QiskitQECError(f"Compiled method {method} only supports n-k < 63.")
             distance = _minimum_distance_2_compiled(stabilizer, gauge, max_weight)
     else:
         logger.exception("from compiled extension was not loaded: switching to no compiled")
         if method == method_enumerate:
             distance = _minimum_distance_1_python(stabilizer, gauge, max_weight)
         elif method == method_partition:
+            if method == method_partition and not stabilizer.shape[0] < 63:
+                raise QiskitQECError(f"Python method {method} only supports n-k < 63.")
             distance = _minimum_distance_2_python(stabilizer, gauge, max_weight)
     return distance

src/qiskit_qec/analysis/faultsampler.py

@@ -1,4 +1,5 @@
 """Quantum circuit fault path sampler."""
+
 from typing import Tuple, List
 
 import numpy as np

src/qiskit_qec/analysis/intern/distance.cpp

@@ -62,6 +62,9 @@ int minimum_distance(std::vector<std::vector<int> > &symplectic_vectors,
  *
  * Pauli operators on n qubits are represented as integer arrays
  * of length 2n in the order [X-part, Z-part].
+ * 
+ * As the powers of 3 are stored in a unit64_t we must have 
+ * symplectic_vectors.size() < 64 to prevent overflow 
  *
  * Returns True or False.
  */
@@ -71,9 +74,10 @@ bool distance_test(std::vector<std::vector<int> > &symplectic_vectors,
 {
   int n = symplectic_vectors[0].size() / 2;  // num. qubits
   int m = symplectic_vectors.size();
-  std::vector<int> pow2;
+  assert(m<63 && "n-k must be < 63 for this version not to overflow.");
+  std::vector<uint64_t> pow2;
   for (int i = 0; i < m+1; i++)
-    pow2.push_back(1 << i);
+    pow2.push_back(static_cast<uint64_t>(1) << i);
   int w1 = 0;
   int w2 = 0;
   if (weight % 2 == 0) {
@@ -83,17 +87,18 @@ bool distance_test(std::vector<std::vector<int> > &symplectic_vectors,
     w1 = (weight + 1) / 2;
     w2 = (weight - 1) / 2;
   }
+
   // Compute syndromes of all single-qubit errors
-  std::vector<int> single_qubit_syndromes;
+  std::vector<uint64_t> single_qubit_syndromes;
   for (int q=0; q < n; q++) {
-    int syndX = 0;
+    uint64_t syndX = 0;
     for (int i=0; i < m; i++) {
       if (symplectic_vectors[i][n + q] == 1)
         syndX += pow2[i];
     }
     if (symplectic_logical_op[n + q] == 1)
       syndX += pow2[m];
-    int syndZ = 0;
+    uint64_t syndZ = 0;
     for (int i=0; i < m; i++) {
       if (symplectic_vectors[i][q] == 1)
         syndZ += pow2[i];
@@ -105,14 +110,14 @@ bool distance_test(std::vector<std::vector<int> > &symplectic_vectors,
     single_qubit_syndromes.push_back(syndX ^ syndZ);
   }
   // Examine all errors with weight w1
-  int mask1 = 1 << m;
-  std::set<int> T1c;
-  std::set<int> T1a;
+  uint64_t mask1 = static_cast<uint64_t>(1) << m;
+  std::set<uint64_t> T1c;
+  std::set<uint64_t> T1a;
   if (w1 > 0) {
     Combinations c(single_qubit_syndromes.size(), w1);
     std::vector<int> state;
     while (c.next_combination(state)) {
-      int synd = 0;
+      uint64_t synd = 0;
       for (int i=0; i < w1; i++)
         synd ^= single_qubit_syndromes[state[i]];
       if (synd & mask1)
@@ -122,14 +127,14 @@ bool distance_test(std::vector<std::vector<int> > &symplectic_vectors,
     }
   }
   // Examine all errors with weight w2
-  std::set<int> T2c;
-  std::set<int> T2a;
+  std::set<uint64_t> T2c;
+  std::set<uint64_t> T2a;
   if (w1 != w2) {
     if (w2 > 0) {
       Combinations c(single_qubit_syndromes.size(), w2);
       std::vector<int> state;
       while (c.next_combination(state)) {
-        int synd = 0;
+        uint64_t synd = 0;
         for (int i=0; i < w2; i++)
           synd ^= single_qubit_syndromes[state[i]];
         if (synd & mask1)
@@ -142,7 +147,7 @@ bool distance_test(std::vector<std::vector<int> > &symplectic_vectors,
     T2c = T1c;
     T2a = T1a;
   }
-  std::set<int> intersect;
+  std::set<uint64_t> intersect;
   std::set_intersection(T1c.begin(), T1c.end(), T2a.begin(), T2a.end(),
                         std::inserter(intersect, intersect.begin()));
   if (intersect.size() > 0) return true;
@@ -169,7 +174,9 @@ int minimum_distance_by_tests(std::vector<std::vector<int> > &symplectic_vectors
                               std::vector<std::vector<int> > &symplectic_xl,
                               std::vector<std::vector<int> > &symplectic_zl,
                               int max_weight) {
+
     int weight = max_weight + 1;
+
     for (int row=0; row < symplectic_xl.size(); row++) {
       for (int w=1; w < max_weight + 1; w++) {
         if (distance_test(symplectic_vectors, symplectic_xl[row], w)) {
@@ -178,6 +185,7 @@ int minimum_distance_by_tests(std::vector<std::vector<int> > &symplectic_vectors
         }
       }
     }
+
     for (int row=0; row < symplectic_zl.size(); row++) {
       for (int w=1; w < max_weight + 1; w++) {
         if (distance_test(symplectic_vectors, symplectic_zl[row], w)) {
@@ -186,6 +194,7 @@ int minimum_distance_by_tests(std::vector<std::vector<int> > &symplectic_vectors
         }
       }
     }
+
     if (weight < max_weight + 1)
       return weight;
     else

src/qiskit_qec/analysis/intern/faultenumerator.cpp

@@ -78,7 +78,7 @@ std::vector<FaultPath> FaultEnumerator::enumerate(int blocksize)
   // Reset if the state is empty
   if(_state.size() == 0) { _index = 0; _done = false; }
   Combinations c(_faulty_ops_indices.size(), _order);
-  long initial_index = _index;
+  long int initial_index = _index;
   while(c.next_combination(_state))
   {
     std::vector<int> indices;

src/qiskit_qec/analysis/intern/productiterator.cpp

@@ -2,7 +2,7 @@
 
 ProductIterator::ProductIterator(std::vector<int> &ni_)
 {
-  for(int i=0; i < ni_.size(); i++)
+  for(uint64_t i=0; i < ni_.size(); i++)
     if(ni_[i] <= 0) throw std::invalid_argument("each element must be positive");
   ni = ni_;
 }

src/qiskit_qec/analysis/productiterator.h

@@ -3,6 +3,8 @@
 
 #include <vector>
 #include <stdexcept>
+#include <cstdint> 
+
 
 // Iterate over a product of finite sets {0, 1, ..., n_i-1}
 // for i = 1, 2, ..., m.

src/qiskit_qec/circuits/css_code.py

@@ -144,7 +144,10 @@ def _get_code_properties(self):
             stabilizers = [[], []]
             logicals = [[], []]
 
-            for (raw_ops, ops,) in zip(
+            for (
+                raw_ops,
+                ops,
+            ) in zip(
                 [raw_gauges, raw_stabilizers, raw_logicals],
                 [gauges, stabilizers, logicals],
             ):

src/qiskit_qec/circuits/surface_code.py

@@ -23,7 +23,6 @@
 
 
 class SurfaceCodeCircuit(CodeCircuit):
-
     """Distance d rotated surface code with  T syndrome measurement rounds."""
 
     def __init__(self, d: int, T: int, basis: str = "z", resets=True):

src/qiskit_qec/decoders/decoding_graph.py

@@ -485,7 +485,7 @@ def get_edge_graph(self):
                             spares.add((n, n2))
                             spares.add((n2, n))
         # find bulk-boundary pairs not covered by the above
-        for (n0, n1) in self.graph.edge_list():
+        for n0, n1 in self.graph.edge_list():
             n2 = None
             for n in (n0, n1):
                 if nodes[n].is_logical:

src/qiskit_qec/decoders/hdrg_decoders.py

@@ -508,9 +508,11 @@ def _grow_clusters(self) -> List[FusionEntry]:
                             fully_grown_edges=set(),
                             edge_support=set(),
                             atypical_nodes=set(),
-                            boundary_nodes=set([edge.neighbour_vertex])
-                            if self.graph[edge.neighbour_vertex].is_logical
-                            else set([]),
+                            boundary_nodes=(
+                                set([edge.neighbour_vertex])
+                                if self.graph[edge.neighbour_vertex].is_logical
+                                else set([])
+                            ),
                             nodes=set([edge.neighbour_vertex]),
                             size=1,
                         )

src/qiskit_qec/decoders/temp_code_util.py

@@ -1,4 +1,5 @@
 """Temporary module with methods for codes."""
+
 from functools import partial
 from typing import List
 

src/qiskit_qec/decoders/temp_graph_util.py

@@ -1,4 +1,5 @@
 """Temporary module with methods for graphs."""
+
 import json
 import networkx as nx
 import rustworkx as rx

src/qiskit_qec/geometry/tiles/diagonalbartile.py

@@ -23,7 +23,6 @@
 # pylint: disable=anomalous-backslash-in-string)
 class DiagonalBarTile(Tile):
     """Diagonal Bar Tile
-
     The diagram is as follows::
 
         q0           q1     q1           q2

src/qiskit_qec/linear/symplectic.py

@@ -196,7 +196,7 @@ def _symplectic_product_vv_boolean(vec1: np.ndarray, vec2: np.ndarray, n: int) -
     assert vec1.dtype == bool
     r = False
     for i in range(n):
-        r += vec1[i] & vec2[n + i] ^ vec1[n + i] & vec2[i]
+        r ^= vec1[i] & vec2[n + i] ^ vec1[n + i] & vec2[i]
     return int(r)
 
 

src/qiskit_qec/noise/paulinoisemodel.py

@@ -1,4 +1,5 @@
 """Pauli circuit-level noise model."""
+
 import copy
 from typing import Union, Dict
 from qiskit_aer import noise

src/qiskit_qec/operators/pauli.py

@@ -704,7 +704,11 @@ def instrs2symplectic(instr: Union[Instruction, QuantumCircuit]):
                 output_encoding=pauli_rep.INTERNAL_PHASE_ENCODING,
             )
         # Recursively apply instructions
-        for dinstr, qregs, cregs in instr.data:
+        # for dinstr, qregs, cregs in instr.data:
+        for instruction in instr.data:
+            dinstr = instruction.operation
+            qregs = instruction.qubits
+            cregs = instruction.clbits
             if cregs:
                 raise QiskitError(
                     f"Cannot apply instruction with classical registers: {dinstr.name}"

src/qiskit_qec/utils/stim_tools.py

@@ -28,7 +28,7 @@
 import rustworkx as rx
 
 from qiskit import QuantumCircuit
-from qiskit_aer.noise.errors.quantum_error import QuantumChannelInstruction
+from qiskit_aer.noise.errors.base_quantum_error import QuantumChannelInstruction
 from qiskit_aer.noise import pauli_error
 from qiskit_qec.utils.decoding_graph_attributes import (
     DecodingGraphNode,
@@ -131,7 +131,11 @@ def get_stim_circuits(
         creg_offset = {}
         prevq_offset = 0
         prevc_offset = 0
-        for inst, qargs, cargs in circ.data:
+
+        for instruction in circ.data:
+            inst = instruction.operation
+            qargs = instruction.qubits
+            cargs = instruction.clbits
             for qubit in qargs:
                 if qubit._register.name not in qreg_offset:
                     qreg_offset[qubit._register.name] = prevq_offset

test/code_circuits/test_rep_codes.py

@@ -86,7 +86,7 @@ def single_error_test(
             circuit_name = {}
             depth = len(qc)
             for j in range(depth):
-                qubits = qc.data[j][1]
+                qubits = qc.data[j].qubits
                 for qubit in qubits:
                     for error in ["x", "y", "z"]:
                         temp_qc = blank_qc.copy()
@@ -101,7 +101,7 @@ def single_error_test(
             job = simulator.run(list(error_circuit.values()), shots=1)
 
             for j in range(depth):
-                qubits = qc.data[j][1]
+                qubits = qc.data[j].qubits
                 for qubit in qubits:
                     for error in ["x", "y", "z"]:
                         results = job.result().get_counts(str((j, qubit, error)))
@@ -352,7 +352,7 @@ def test_single_error_202s(self):
                         error_qc.add_register(creg)
                     barrier_num = 0
                     for gate in qc.data:
-                        if gate[0].name == "barrier":
+                        if gate.operation.name == "barrier":
                             barrier_num += 1
                             if barrier_num == 2 * t + 1:
                                 if q % 2 == 0:

test/distance/test_minimum_distance.py

@@ -1,4 +1,5 @@
 """Test minimum distance computations."""
+
 import unittest
 import numpy as np