[WIP] New algorithm for RB using transpiled Cliffords

docs/tutorials/randomized_benchmarking.rst

@@ -102,6 +102,13 @@ interleaved RB experiment will always give you accurate error value :math:`e_i`.
 
     # Run an RB experiment on qubit 0
     exp1 = StandardRB(qubits, lengths, num_samples=num_samples, seed=seed)
+
+    # basis_gates must be set for randomized benchmarking
+    transpiler_options = {
+        "basis_gates":  ["rz", "sx", "cx"],
+        "optimization_level": 1,
+        }
+    exp1.set_transpile_options(**transpiler_options)
     expdata1 = exp1.run(backend).block_for_results()
     results1 = expdata1.analysis_results()
 
@@ -179,14 +186,15 @@ The EPGs of two-qubit RB are analyzed with the corrected EPC if available.
         ],
         flatten_results=True,
     )
+    single_exps.set_transpile_options(**transpiler_options)
     expdata_1q = single_exps.run(backend).block_for_results()
 
 
 .. jupyter-execute::
 
     # Run an RB experiment on qubits 1, 4
     exp_2q = StandardRB(qubits, lengths_2_qubit, num_samples=num_samples, seed=seed)
-
+    exp_2q.set_transpile_options(**transpiler_options)
     # Use the EPG data of the 1-qubit runs to ensure correct 2-qubit EPG computation
     exp_2q.analysis.set_options(epg_1_qubit=expdata_1q.analysis_results())
 
@@ -213,6 +221,7 @@ Generating an example RB circuit:
 
     # Run an RB experiment on qubit 0
     exp = StandardRB(qubits=[0], lengths=[10], num_samples=1, seed=seed)
+    exp.set_transpile_options(**transpiler_options)
     c = exp.circuits()[0]
 
 We transpile the circuit into the backend’s basis gate set:
@@ -264,7 +273,8 @@ Running a 1-qubit interleaved RB experiment
     # The interleaved gate is the x gate
     int_exp1 = InterleavedRB(
         circuits.XGate(), qubits, lengths, num_samples=num_samples, seed=seed)
-
+    int_exp1.set_transpile_options(**transpiler_options)
+
     # Run
     int_expdata1 = int_exp1.run(backend).block_for_results()
     int_results1 = int_expdata1.analysis_results()
@@ -291,7 +301,8 @@ Running a 2-qubit interleaved RB experiment
     # The interleaved gate is the cx gate
     int_exp2 = InterleavedRB(
         circuits.CXGate(), qubits, lengths, num_samples=num_samples, seed=seed)
-
+    int_exp2.set_transpile_options(**transpiler_options)
+
     # Run
     int_expdata2 = int_exp2.run(backend).block_for_results()
     int_results2 = int_expdata2.analysis_results()
@@ -322,6 +333,7 @@ different qubits (see Ref. [5])
     exps = [StandardRB([i], lengths, num_samples=num_samples, seed=seed + i)
             for i in qubits]
     par_exp = ParallelExperiment(exps)
+    par_exp.set_transpile_options(**transpiler_options)
     par_expdata = par_exp.run(backend).block_for_results()
     par_results = par_expdata.analysis_results()
 

qiskit_experiments/library/randomized_benchmarking/clifford_utils.py

@@ -13,13 +13,26 @@
 Utilities for using the Clifford group in randomized benchmarking
 """
 
-from typing import Optional, Union
+import os
+from typing import List
 from functools import lru_cache
-from numpy.random import Generator, default_rng
+from math import isclose
+import numpy as np
+
 from qiskit import QuantumCircuit, QuantumRegister
 from qiskit.circuit import Gate
 from qiskit.circuit.library import SdgGate, HGate, SGate, SXdgGate
-from qiskit.quantum_info import Clifford, random_clifford
+from qiskit.exceptions import QiskitError
+from qiskit.quantum_info import Clifford
+
+from .clifford_data import (
+    CLIFF_SINGLE_GATE_MAP_1Q,
+    CLIFF_SINGLE_GATE_MAP_2Q,
+    CLIFF_COMPOSE_DATA_1Q,
+    CLIFF_COMPOSE_DATA_2Q,
+    CLIFF_INVERSE_DATA_1Q,
+    CLIFF_INVERSE_DATA_2Q,
+)
 
 
 class VGate(Gate):
@@ -64,67 +77,35 @@ class CliffordUtils:
         (2, 2, 3, 3, 3, 3, 4, 4),
         (2, 2, 3, 3, 4, 4),
     ]
+    GENERAL_CLIFF_LIST = ["id", "h", "sdg", "s", "x", "sx", "sxdg", "y", "z", "cx"]
+    TRANSPILED_CLIFF_LIST = ["sx", "rz", "cx"]
+    CLIFF_SINGLE_GATE_MAP = {1: CLIFF_SINGLE_GATE_MAP_1Q, 2: CLIFF_SINGLE_GATE_MAP_2Q}
+    CLIFF_COMPOSE_DATA = {1: CLIFF_COMPOSE_DATA_1Q, 2: CLIFF_COMPOSE_DATA_2Q}
+    CLIFF_INVERSE_DATA = {1: CLIFF_INVERSE_DATA_1Q, 2: CLIFF_INVERSE_DATA_2Q}
 
-    def clifford_1_qubit(self, num):
+    @classmethod
+    def clifford_1_qubit(cls, num):
         """Return the 1-qubit clifford element corresponding to `num`
         where `num` is between 0 and 23.
         """
-        return Clifford(self.clifford_1_qubit_circuit(num), validate=False)
+        return Clifford(cls.clifford_1_qubit_circuit(num), validate=False)
 
-    def clifford_2_qubit(self, num):
+    @classmethod
+    def clifford_2_qubit(cls, num):
         """Return the 2-qubit clifford element corresponding to `num`
         where `num` is between 0 and 11519.
         """
-        return Clifford(self.clifford_2_qubit_circuit(num), validate=False)
-
-    def random_cliffords(
-        self, num_qubits: int, size: int = 1, rng: Optional[Union[int, Generator]] = None
-    ):
-        """Generate a list of random clifford elements"""
-        if num_qubits > 2:
-            return random_clifford(num_qubits, seed=rng)
-
-        if rng is None:
-            rng = default_rng()
-
-        if isinstance(rng, int):
-            rng = default_rng(rng)
-
-        if num_qubits == 1:
-            samples = rng.integers(24, size=size)
-            return [Clifford(self.clifford_1_qubit_circuit(i), validate=False) for i in samples]
-        else:
-            samples = rng.integers(11520, size=size)
-            return [Clifford(self.clifford_2_qubit_circuit(i), validate=False) for i in samples]
-
-    def random_clifford_circuits(
-        self, num_qubits: int, size: int = 1, rng: Optional[Union[int, Generator]] = None
-    ):
-        """Generate a list of random clifford circuits"""
-        if num_qubits > 2:
-            return [random_clifford(num_qubits, seed=rng).to_circuit() for _ in range(size)]
-
-        if rng is None:
-            rng = default_rng()
-
-        if isinstance(rng, int):
-            rng = default_rng(rng)
-
-        if num_qubits == 1:
-            samples = rng.integers(24, size=size)
-            return [self.clifford_1_qubit_circuit(i) for i in samples]
-        else:
-            samples = rng.integers(11520, size=size)
-            return [self.clifford_2_qubit_circuit(i) for i in samples]
+        return Clifford(cls.clifford_2_qubit_circuit(num), validate=False)
 
+    @classmethod
     @lru_cache(maxsize=24)
-    def clifford_1_qubit_circuit(self, num):
+    def clifford_1_qubit_circuit(cls, num):
         """Return the 1-qubit clifford circuit corresponding to `num`
         where `num` is between 0 and 23.
         """
         # pylint: disable=unbalanced-tuple-unpacking
         # This is safe since `_unpack_num` returns list the size of the sig
-        (i, j, p) = self._unpack_num(num, self.CLIFFORD_1_QUBIT_SIG)
+        (i, j, p) = cls._unpack_num(num, cls.CLIFFORD_1_QUBIT_SIG)
         qr = QuantumRegister(1)
         qc = QuantumCircuit(qr)
         if i == 1:
@@ -141,12 +122,13 @@ def clifford_1_qubit_circuit(self, num):
             qc.z(0)
         return qc
 
+    @classmethod
     @lru_cache(maxsize=11520)
-    def clifford_2_qubit_circuit(self, num):
+    def clifford_2_qubit_circuit(cls, num):
         """Return the 2-qubit clifford circuit corresponding to `num`
         where `num` is between 0 and 11519.
         """
-        vals = self._unpack_num_multi_sigs(num, self.CLIFFORD_2_QUBIT_SIGS)
+        vals = cls._unpack_num_multi_sigs(num, cls.CLIFFORD_2_QUBIT_SIGS)
         qr = QuantumRegister(2)
         qc = QuantumCircuit(qr)
         if vals[0] == 0 or vals[0] == 3:
@@ -195,7 +177,8 @@ def clifford_2_qubit_circuit(self, num):
             qc.z(1)
         return qc
 
-    def _unpack_num(self, num, sig):
+    @staticmethod
+    def _unpack_num(num, sig):
         r"""Returns a tuple :math:`(a_1, \ldots, a_n)` where
         :math:`0 \le a_i \le \sigma_i` where
         sig=:math:`(\sigma_1, \ldots, \sigma_n)` and num is the sequential
@@ -207,7 +190,8 @@ def _unpack_num(self, num, sig):
             num //= k
         return res
 
-    def _unpack_num_multi_sigs(self, num, sigs):
+    @classmethod
+    def _unpack_num_multi_sigs(cls, num, sigs):
         """Returns the result of `_unpack_num` on one of the
         signatures in `sigs`
         """
@@ -216,6 +200,104 @@ def _unpack_num_multi_sigs(self, num, sigs):
             for k in sig:
                 sig_size *= k
             if num < sig_size:
-                return [i] + self._unpack_num(num, sig)
+                return [i] + cls._unpack_num(num, sig)
             num -= sig_size
         return None
+
+    @classmethod
+    def num_from_clifford_single_gate(cls, inst, qubits, rb_num_qubits, basis_gates):
+        """
+        This method does the reverse of clifford_1_qubit_circuit and clifford_2_qubit_circuit -
+        given a clifford, it returns the corresponding integer, with the mapping
+        defined in the above method.
+        The mapping is in the context of the basis_gates. Therefore, we define here
+        the possible supersets of basis gates, and verify that the given inst belong to
+        one of these sets.
+        """
+        name = inst.name
+        gates_with_delay = basis_gates.copy()
+        gates_with_delay.append("delay")
+
+        if not name in gates_with_delay:
+            raise QiskitError("Instruction {} is not in the basis gates".format(inst.name))
+        if set(basis_gates).issubset(set(cls.GENERAL_CLIFF_LIST)):
+            if name == "delay":
+                return 0
+            map_index = name
+
+        if set(basis_gates).issubset(set(cls.TRANSPILED_CLIFF_LIST)):
+            if name in {"sx", "cx"}:
+                map_index = name
+            elif name == "delay":
+                return 0
+            elif name == "rz":
+                # The next two are identical up to a phase, which makes no difference
+                # for the associated Cliffords
+                if isclose(inst.params[0], np.pi) or isclose(inst.params[0], -np.pi):
+                    map_index = "z"
+                elif isclose(inst.params[0], np.pi / 2):
+                    map_index = "s"
+                elif isclose(inst.params[0], -np.pi / 2):
+                    map_index = "sdg"
+                else:
+                    raise QiskitError("wrong param {} for rz in clifford".format(inst.params[0]))
+            else:
+                raise QiskitError(
+                    "Instruction {} could not be converted to Clifford gate".format(name)
+                )
+
+        return cls.CLIFF_SINGLE_GATE_MAP[rb_num_qubits][(map_index, str(qubits))]
+
+    @classmethod
+    def compose_num_with_clifford(
+        cls, num_qubits: int, composed_num: int, qc: QuantumCircuit, basis_gates: List[str]
+    ) -> int:
+        """Compose a number that represents a Clifford, with a single-gate Clifford, and return the
+        number that represents the resulting Clifford."""
+
+        # The numbers corresponding to single gate Cliffords are not in sequence -
+        # see num_from_1q_clifford_single_gate. To compute the index in
+        # the array CLIFF_COMPOSE_DATA_1Q, we map the numbers to [0, 8].
+        map_clifford_num_to_array_index = {}
+        num_single_gate_cliffs = len(cls.CLIFF_SINGLE_GATE_MAP[num_qubits])
+        for k in list(cls.CLIFF_SINGLE_GATE_MAP[num_qubits]):
+            map_clifford_num_to_array_index[cls.CLIFF_SINGLE_GATE_MAP[num_qubits][k]] = list(
+                cls.CLIFF_SINGLE_GATE_MAP[num_qubits].keys()
+            ).index(k)
+        if num_qubits == 1:
+            for inst, qargs, _ in qc:
+                num = cls.num_from_clifford_single_gate(
+                    inst=inst, qubits=[0], rb_num_qubits=1, basis_gates=basis_gates
+                )
+                index = num_single_gate_cliffs * composed_num + map_clifford_num_to_array_index[num]
+                composed_num = cls.CLIFF_COMPOSE_DATA[num_qubits][index]
+        else:
+            for inst, qargs, _ in qc:
+                if inst.num_qubits == 2:
+                    qubits = [qc.find_bit(qargs[0]).index, qc.find_bit(qargs[1]).index]
+                else:
+                    qubits = [qc.find_bit(qargs[0]).index]
+                num = cls.num_from_clifford_single_gate(
+                    inst=inst, qubits=qubits, rb_num_qubits=2, basis_gates=basis_gates
+                )
+                index = num_single_gate_cliffs * composed_num + map_clifford_num_to_array_index[num]
+                composed_num = cls.CLIFF_COMPOSE_DATA[num_qubits][index]
+        return composed_num
+
+    @classmethod
+    def clifford_inverse_by_num(cls, num: int, num_qubits: int):
+        """Return the number of the inverse Clifford to the input num"""
+        return cls.CLIFF_INVERSE_DATA[num_qubits][num]
+
+    @staticmethod
+    def file_name(num_qubits, basis_gates):
+        """Return the name of the file containing the transpiled Cliffords"""
+        suffix = ""
+        for n in basis_gates[0:-1]:
+            suffix += "_" + n
+        if num_qubits == 2:
+            suffix += "_" + basis_gates[-1]
+        circs_file_name = "/transpiled_circs_" + str(num_qubits) + "q" + suffix + ".qpy"
+        root_dir = os.path.dirname(os.path.abspath(__file__))
+        transpiled_circs_file = root_dir + circs_file_name
+        return transpiled_circs_file