Improve performance of multiplication of cliffords

qiskit_experiments/library/randomized_benchmarking/clifford_utils.py

@@ -17,10 +17,9 @@
 import os
 from functools import lru_cache
 from numbers import Integral
-from typing import Optional, Union, Tuple, Sequence
+from typing import Optional, Union, Tuple, Sequence, Iterable
 
 import numpy as np
-import scipy.sparse
 from numpy.random import Generator, default_rng
 
 from qiskit.circuit import CircuitInstruction, Qubit
@@ -36,11 +35,13 @@
 
 _CLIFFORD_COMPOSE_1Q = np.load(f"{_DATA_FOLDER}/clifford_compose_1q.npz")["table"]
 _CLIFFORD_INVERSE_1Q = np.load(f"{_DATA_FOLDER}/clifford_inverse_1q.npz")["table"]
-_CLIFFORD_COMPOSE_2Q = scipy.sparse.lil_matrix(
-    scipy.sparse.load_npz(f"{_DATA_FOLDER}/clifford_compose_2q_sparse.npz")
-)
 _CLIFFORD_INVERSE_2Q = np.load(f"{_DATA_FOLDER}/clifford_inverse_2q.npz")["table"]
-
+_clifford_compose_2q_data = np.load(f"{_DATA_FOLDER}/clifford_compose_2q_dense_selected.npz")
+_CLIFFORD_COMPOSE_2Q_DENSE = _clifford_compose_2q_data["table"]
+# valid indices for the columns of the _CLIFFORD_COMPOSE_2Q_DENSE table
+_valid_sparse_indices = _clifford_compose_2q_data["valid_sparse_indices"]
+# map a clifford number to the index of _CLIFFORD_COMPOSE_2Q_DENSE
+_clifford_num_to_dense_index = {idx: ii for ii, idx in enumerate(_valid_sparse_indices)}
 
 # Transpilation utilities
 def _transpile_clifford_circuit(
@@ -422,8 +423,9 @@ def compose_2q(lhs: Integral, rhs: Integral) -> Integral:
     """Return the composition of 2-qubit clifford integers."""
     num = lhs
     for layer, idx in enumerate(_layer_indices_from_num(rhs)):
-        circ = _CLIFFORD_LAYER[layer][idx]
-        num = _compose_num_with_circuit_2q(num, circ)
+        gate_numbers = _CLIFFORD_LAYER_NUMS[layer][idx]
+        for n in gate_numbers:
+            num = _CLIFFORD_COMPOSE_2Q_DENSE[num, _clifford_num_to_dense_index[n]]
     return num
 
 
@@ -434,17 +436,9 @@ def inverse_2q(num: Integral) -> Integral:
 
 def num_from_2q_circuit(qc: QuantumCircuit) -> Integral:
     """Convert a given 2-qubit Clifford circuit to the corresponding integer."""
-    return _compose_num_with_circuit_2q(0, qc)
-
-
-def _compose_num_with_circuit_2q(num: Integral, qc: QuantumCircuit) -> Integral:
-    """Compose a number that represents a Clifford, with a Clifford circuit, and return the
-    number that represents the resulting Clifford."""
-    lhs = num
-    for inst in qc:
-        qubits = tuple(qc.find_bit(q).index for q in inst.qubits)
-        rhs = _num_from_2q_gate(op=inst.operation, qubits=qubits)
-        lhs = _CLIFFORD_COMPOSE_2Q[lhs, rhs]
+    lhs = 0
+    for rhs in _clifford_2q_nums_from_2q_circuit(qc):
+        lhs = _CLIFFORD_COMPOSE_2Q_DENSE[lhs, _clifford_num_to_dense_index[rhs]]
     return lhs
 
 
@@ -570,6 +564,20 @@ def _create_cliff_2q_layer_2():
 _NUM_LAYER_2 = 16
 
 
+def _clifford_2q_nums_from_2q_circuit(qc: QuantumCircuit) -> Iterable[Integral]:
+    """Yield Clifford numbers that represents the 2Q Clifford circuit."""
+    for inst in qc:
+        qubits = tuple(qc.find_bit(q).index for q in inst.qubits)
+        yield _num_from_2q_gate(op=inst.operation, qubits=qubits)
+
+
+# Construct mapping from Clifford layers to series of Clifford numbers
+_CLIFFORD_LAYER_NUMS = [
+    [tuple(_clifford_2q_nums_from_2q_circuit(qc)) for qc in _CLIFFORD_LAYER[layer]]
+    for layer in [0, 1, 2]
+]
+
+
 @lru_cache(maxsize=None)
 def _transformed_clifford_layer(
     layer: int, index: Integral, basis_gates: Tuple[str, ...]

qiskit_experiments/library/randomized_benchmarking/data/generate_clifford_data.py

@@ -109,7 +109,27 @@ def gen_clifford_compose_2q_gate():
         for rhs in _CLIFF_SINGLE_GATE_MAP_2Q.values():
             composed = cliff_lhs.compose(_CLIFF_2Q[rhs])
             products[lhs, rhs] = _TO_INT_2Q[_hash_cliff(composed)]
-    return products.tocsr()
+    return products.tocsc()
+
+
+def gen_clifford_compose_2q_dense() -> tuple[np.typing.NDArray[int], list[int]]:
+    """Generate a dense multiplication table for 1-qubit Clifford numbers
+
+    The multiplication table is generated from the sparse table generated by :meth:`gen_clifford_compose_2q_gate`.
+
+    Returns:
+        Tuple with a dense multiplication table and the valid indices for the columns
+    """
+    _CLIFFORD_COMPOSE_2Q = gen_clifford_compose_2q_gate()
+    number_of_cliffords = _CLIFFORD_COMPOSE_2Q.shape[0]
+    valid_sparse_indices = [
+        num
+        for num in range(number_of_cliffords)
+        if _CLIFFORD_COMPOSE_2Q[:, num].nnz == number_of_cliffords - 1
+    ]
+    _CLIFFORD_COMPOSE_2Q_DENSE = _CLIFFORD_COMPOSE_2Q[:, valid_sparse_indices].toarray()
+
+    return _CLIFFORD_COMPOSE_2Q_DENSE, valid_sparse_indices
 
 
 _GATE_LIST_1Q = [
@@ -184,4 +204,11 @@ def gen_cliff_single_2q_gate_map():
             "_CLIFF_SINGLE_GATE_MAP_2Q must be generated by gen_cliff_single_2q_gate_map()"
         )
     np.savez_compressed("clifford_inverse_2q.npz", table=gen_clifford_inverse_2q())
-    scipy.sparse.save_npz("clifford_compose_2q_sparse.npz", gen_clifford_compose_2q_gate())
+
+    _CLIFFORD_COMPOSE_2Q_DENSE, valid_sparse_indices = gen_clifford_compose_2q_dense()
+
+    np.savez_compressed(
+        "clifford_compose_2q_dense_selected.npz",
+        table=_CLIFFORD_COMPOSE_2Q_DENSE,
+        valid_sparse_indices=valid_sparse_indices,
+    )

test/library/randomized_benchmarking/test_clifford_utils.py

@@ -20,6 +20,7 @@
 from numpy.random import default_rng
 
 from qiskit import QuantumCircuit
+from qiskit.exceptions import QiskitError
 from qiskit.circuit.library import (
     IGate,
     XGate,
@@ -43,6 +44,7 @@
     _num_from_layer_indices,
     _layer_indices_from_num,
     _CLIFFORD_LAYER,
+    _CLIFFORD_INVERSE_2Q,
 )
 
 
@@ -195,3 +197,27 @@ def test_num_from_layer(self):
                 circ.compose(_CLIFFORD_LAYER[layer][idx], inplace=True)
             layered = Clifford(circ)
             self.assertEqual(standard, layered)
+
+    def test_num_from_2q_circuit(self):
+        """Check conversion of circuits to integers with num_from_2q_circuit"""
+        qc = QuantumCircuit(2)
+        qc.h(0)
+        num = num_from_2q_circuit(qc)
+        self.assertEqual(num, 5760)
+        qc = QuantumCircuit(2)
+        qc.u(0, 0, np.pi, 0)
+        with self.assertRaises(QiskitError):
+            # raising an error is ok, num_from_2q_circuit does not support all 2-qubit gates
+            num_from_2q_circuit(qc)
+
+        # regression test for using the dense multiplication table
+        qc = QuantumCircuit(2)
+        qc.cz(1, 0)
+        num = num_from_2q_circuit(qc)
+        self.assertEqual(num, 368)
+
+    def test_clifford_inverse_table(self):
+        """Check correctness of the Clifford inversion table"""
+        for lhs, rhs in enumerate(_CLIFFORD_INVERSE_2Q):
+            c = compose_2q(lhs, rhs)
+            self.assertEqual(c, 0)