make SingleQubitCliffordGate immutable singletons and use it in qubit_characterizations for a 37% speedup

cirq-core/cirq/experiments/qubit_characterizations.py

@@ -14,6 +14,7 @@
 
 import dataclasses
 import itertools
+import functools
 
 from typing import (
     Any,
@@ -58,11 +59,11 @@ class Cliffords:
         s1_y
     """
 
-    c1_in_xy: List[List[ops.Gate]]
-    c1_in_xz: List[List[ops.Gate]]
-    s1: List[List[ops.Gate]]
-    s1_x: List[List[ops.Gate]]
-    s1_y: List[List[ops.Gate]]
+    c1_in_xy: List[List[ops.SingleQubitCliffordGate]]
+    c1_in_xz: List[List[ops.SingleQubitCliffordGate]]
+    s1: List[List[ops.SingleQubitCliffordGate]]
+    s1_x: List[List[ops.SingleQubitCliffordGate]]
+    s1_y: List[List[ops.SingleQubitCliffordGate]]
 
 
 class RandomizedBenchMarkResult:
@@ -299,17 +300,14 @@ def parallel_single_qubit_randomized_benchmarking(
         A dictionary from qubits to RandomizedBenchMarkResult objects.
     """
 
-    cliffords = _single_qubit_cliffords()
-    c1 = cliffords.c1_in_xy if use_xy_basis else cliffords.c1_in_xz
-    clifford_mats = np.array([_gate_seq_to_mats(gates) for gates in c1])
+    clifford_group = _single_qubit_cliffords()
+    c1 = clifford_group.c1_in_xy if use_xy_basis else clifford_group.c1_in_xz
 
     # create circuits
     circuits_all: List['cirq.AbstractCircuit'] = []
     for num_cliffords in num_clifford_range:
         for _ in range(num_circuits):
-            circuits_all.append(
-                _create_parallel_rb_circuit(qubits, num_cliffords, c1, clifford_mats)
-            )
+            circuits_all.append(_create_parallel_rb_circuit(qubits, num_cliffords, c1))
 
     # run circuits
     results = sampler.run_batch(circuits_all, repetitions=repetitions)
@@ -562,11 +560,9 @@ def _measurement(two_qubit_circuit: circuits.Circuit) -> np.ndarray:
 
 
 def _create_parallel_rb_circuit(
-    qubits: Iterator['cirq.Qid'], num_cliffords: int, c1: list, clifford_mats: np.ndarray
+    qubits: Iterator['cirq.Qid'], num_cliffords: int, c1: list
 ) -> 'cirq.Circuit':
-    circuits_to_zip = [
-        _random_single_q_clifford(qubit, num_cliffords, c1, clifford_mats) for qubit in qubits
-    ]
+    circuits_to_zip = [_random_single_q_clifford(qubit, num_cliffords, c1) for qubit in qubits]
     circuit = circuits.Circuit.zip(*circuits_to_zip)
     return circuits.Circuit.from_moments(*circuit, ops.measure_each(*qubits))
 
@@ -612,16 +608,12 @@ def _two_qubit_clifford_matrices(
 
 
 def _random_single_q_clifford(
-    qubit: 'cirq.Qid',
-    num_cfds: int,
-    cfds: Sequence[Sequence['cirq.Gate']],
-    cfd_matrices: np.ndarray,
+    qubit: 'cirq.Qid', num_cfds: int, cfds: Sequence[Sequence['cirq.Gate']]
 ) -> 'cirq.Circuit':
     clifford_group_size = 24
     gate_ids = list(np.random.choice(clifford_group_size, num_cfds))
     gate_sequence = [gate for gate_id in gate_ids for gate in cfds[gate_id]]
-    idx = _find_inv_matrix(_gate_seq_to_mats(gate_sequence), cfd_matrices)
-    gate_sequence.extend(cfds[idx])
+    gate_sequence.append(_reduce_gate_seq(gate_sequence) ** -1)
     circuit = circuits.Circuit(gate(qubit) for gate in gate_sequence)
     return circuit
 
@@ -681,11 +673,13 @@ def _matrix_bar_plot(
         ax.set_title(title)
 
 
-def _gate_seq_to_mats(gate_seq: Sequence['cirq.Gate']) -> np.ndarray:
-    mat_rep = protocols.unitary(gate_seq[0])
+def _reduce_gate_seq(
+    gate_seq: Sequence[ops.SingleQubitCliffordGate],
+) -> ops.SingleQubitCliffordGate:
+    cur = gate_seq[0]
     for gate in gate_seq[1:]:
-        mat_rep = np.dot(protocols.unitary(gate), mat_rep)
-    return mat_rep
+        cur = cur.merged_with(gate)
+    return cur
 
 
 def _two_qubit_clifford(
@@ -793,11 +787,16 @@ def _single_qubit_gates(
         yield gate(qubit)
 
 
+@functools.cache
 def _single_qubit_cliffords() -> Cliffords:
-    X, Y, Z = ops.X, ops.Y, ops.Z
+    X, Y, Z = (
+        ops.SingleQubitCliffordGate.X,
+        ops.SingleQubitCliffordGate.Y,
+        ops.SingleQubitCliffordGate.Z,
+    )
 
-    c1_in_xy: List[List['cirq.Gate']] = []
-    c1_in_xz: List[List['cirq.Gate']] = []
+    c1_in_xy: List[List[ops.SingleQubitCliffordGate]] = []
+    c1_in_xz: List[List[ops.SingleQubitCliffordGate]] = []
 
     for phi_0, phi_1 in itertools.product([1.0, 0.5, -0.5], [0.0, 0.5, -0.5]):
         c1_in_xy.append([X**phi_0, Y**phi_1])
@@ -820,8 +819,20 @@ def _single_qubit_cliffords() -> Cliffords:
     for z0, x, z1 in phi_xz:
         c1_in_xz.append([Z**z0, X**x, Z**z1])
 
-    s1: List[List['cirq.Gate']] = [[X**0.0], [Y**0.5, X**0.5], [X**-0.5, Y**-0.5]]
-    s1_x: List[List['cirq.Gate']] = [[X**0.5], [X**0.5, Y**0.5, X**0.5], [Y**-0.5]]
-    s1_y: List[List['cirq.Gate']] = [[Y**0.5], [X**-0.5, Y**-0.5, X**0.5], [Y, X**0.5]]
+    s1: List[List[ops.SingleQubitCliffordGate]] = [
+        [X**0.0],
+        [Y**0.5, X**0.5],
+        [X**-0.5, Y**-0.5],
+    ]
+    s1_x: List[List[ops.SingleQubitCliffordGate]] = [
+        [X**0.5],
+        [X**0.5, Y**0.5, X**0.5],
+        [Y**-0.5],
+    ]
+    s1_y: List[List[ops.SingleQubitCliffordGate]] = [
+        [Y**0.5],
+        [X**-0.5, Y**-0.5, X**0.5],
+        [Y, X**0.5],
+    ]
 
     return Cliffords(c1_in_xy, c1_in_xz, s1, s1_x, s1_y)

cirq-core/cirq/experiments/qubit_characterizations_test.py

@@ -75,9 +75,32 @@ def check_distinct(unitaries):
 
     # Check that XZ decomposition has at most one X gate per clifford.
     for gates in cliffords.c1_in_xz:
-        num_x = len([gate for gate in gates if isinstance(gate, cirq.XPowGate)])
-        num_z = len([gate for gate in gates if isinstance(gate, cirq.ZPowGate)])
-        assert num_x + num_z == len(gates)
+        num_i = len([gate for gate in gates if gate == cirq.ops.SingleQubitCliffordGate.I])
+        num_x = len(
+            [
+                gate
+                for gate in gates
+                if gate
+                in (
+                    cirq.ops.SingleQubitCliffordGate.X,
+                    cirq.ops.SingleQubitCliffordGate.X_sqrt,
+                    cirq.ops.SingleQubitCliffordGate.X_nsqrt,
+                )
+            ]
+        )
+        num_z = len(
+            [
+                gate
+                for gate in gates
+                if gate
+                in (
+                    cirq.ops.SingleQubitCliffordGate.Z,
+                    cirq.ops.SingleQubitCliffordGate.Z_sqrt,
+                    cirq.ops.SingleQubitCliffordGate.Z_nsqrt,
+                )
+            ]
+        )
+        assert num_x + num_z + num_i == len(gates)
         assert num_x <= 1
 
 

cirq-core/cirq/ops/clifford_gate.py

@@ -14,11 +14,13 @@
 
 from typing import Any, Dict, List, Optional, Sequence, Tuple, TYPE_CHECKING, Union
 
-
+import functools
+from dataclasses import dataclass
 import numpy as np
 
 from cirq import protocols, value, linalg, qis
 from cirq._import import LazyLoader
+from cirq._compat import cached_property, cached_method
 from cirq.ops import common_gates, named_qubit, raw_types, pauli_gates, phased_x_z_gate
 from cirq.ops.pauli_gates import Pauli
 from cirq.type_workarounds import NotImplementedType
@@ -356,6 +358,8 @@ def _get_sqrt_map(
 class CliffordGate(raw_types.Gate, CommonCliffordGates):
     """Clifford rotation for N-qubit."""
 
+    _clifford_tableau: qis.CliffordTableau
+
     def __init__(self, *, _clifford_tableau: qis.CliffordTableau) -> None:
         # We use the Clifford tableau to represent a Clifford gate.
         # It is crucial to note that the meaning of tableau here is different
@@ -376,7 +380,7 @@ def __init__(self, *, _clifford_tableau: qis.CliffordTableau) -> None:
         # more precisely the conjugate transformation of ZI by this gate, becomes -ZI.
         # (Note the real clifford tableau has to satify the Symplectic property.
         # here is just for illustration)
-        self._clifford_tableau = _clifford_tableau.copy()
+        object.__setattr__(self, '_clifford_tableau', _clifford_tableau.copy())
 
     @property
     def clifford_tableau(self):
@@ -399,6 +403,12 @@ def _has_stabilizer_effect_(self) -> Optional[bool]:
     def __pow__(self, exponent) -> 'CliffordGate':
         if exponent == -1:
             return CliffordGate.from_clifford_tableau(self.clifford_tableau.inverse())
+        if exponent == 0:
+            return CliffordGate.from_clifford_tableau(
+                qis.CliffordTableau(num_qubits=self._num_qubits_())
+            )
+        if exponent == 1:
+            return self
         if exponent > 0 and int(exponent) == exponent:
             base_tableau = self.clifford_tableau.copy()
             for _ in range(int(exponent) - 1):
@@ -457,6 +467,7 @@ def _act_on_(
         return NotImplemented
 
 
+@dataclass(frozen=True, init=False, eq=False, repr=False)
 @value.value_equality(manual_cls=True)
 class SingleQubitCliffordGate(CliffordGate):
     """Any single qubit Clifford rotation."""
@@ -468,6 +479,7 @@ def _num_qubits_(self):
         return 1
 
     @staticmethod
+    @functools.cache
     def from_clifford_tableau(tableau: qis.CliffordTableau) -> 'SingleQubitCliffordGate':
         if not isinstance(tableau, qis.CliffordTableau):
             raise ValueError('Input argument has to be a CliffordTableau instance.')
@@ -679,6 +691,10 @@ def to_phased_xz_gate(self) -> phased_x_z_gate.PhasedXZGate:
             * {middle point of xyz in 4 Quadrant} * 120 is [[0, 1], [1, 1]]
             * {middle point of xyz in 4 Quadrant} * 240 is [[1, 1], [1, 0]]
         """
+        return self._to_phased_xz_gate
+
+    @cached_property
+    def _to_phased_xz_gate(self) -> phased_x_z_gate.PhasedXZGate:
         x_to_flip, z_to_flip = self.clifford_tableau.rs
         flip_index = int(z_to_flip) * 2 + x_to_flip
         a, x, z = 0.0, 0.0, 0.0
@@ -716,7 +732,7 @@ def to_phased_xz_gate(self) -> phased_x_z_gate.PhasedXZGate:
             z = -0.5 if x_to_flip else 0.5
         return phased_x_z_gate.PhasedXZGate(x_exponent=x, z_exponent=z, axis_phase_exponent=a)
 
-    def __pow__(self, exponent) -> 'SingleQubitCliffordGate':
+    def __pow__(self, exponent: Union[float, int]) -> 'SingleQubitCliffordGate':
         # First to check if we can get the sqrt and negative sqrt Clifford.
         if self._get_sqrt_map().get(exponent, None):
             pow_gate = self._get_sqrt_map()[exponent].get(self, None)
@@ -761,6 +777,7 @@ def commutes_with_pauli(self, pauli: Pauli) -> bool:
         to, flip = self.pauli_tuple(pauli)
         return to == pauli and not flip
 
+    @cached_method
     def merged_with(self, second: 'SingleQubitCliffordGate') -> 'SingleQubitCliffordGate':
         """Returns a SingleQubitCliffordGate such that the circuits
             --output-- and --self--second--
@@ -773,6 +790,10 @@ def _has_unitary_(self) -> bool:
         return True
 
     def _unitary_(self) -> np.ndarray:
+        return self._unitary
+
+    @cached_property
+    def _unitary(self) -> np.ndarray:
         mat = np.eye(2)
         qubit = named_qubit.NamedQubit('arbitrary')
         for op in protocols.decompose_once_with_qubits(self, (qubit,)):
@@ -787,6 +808,10 @@ def decompose_gate(self) -> Sequence['cirq.Gate']:
             clifford gate if applied in order. This decomposition agrees with
             cirq.unitary(self), including global phase.
         """
+        return self._decompose_gate
+
+    @cached_property
+    def _decompose_gate(self) -> Sequence['cirq.Gate']:
         if self == SingleQubitCliffordGate.H:
             return [common_gates.H]
         rotations = self.decompose_rotation()
@@ -802,6 +827,10 @@ def decompose_rotation(self) -> Sequence[Tuple[Pauli, int]]:
         Note that the combined unitary effect of these rotations may
         differ from cirq.unitary(self) by a global phase.
         """
+        return self._decompose_rotation
+
+    @cached_property
+    def _decompose_rotation(self) -> Sequence[Tuple[Pauli, int]]:
         x_rot = self.pauli_tuple(pauli_gates.X)
         y_rot = self.pauli_tuple(pauli_gates.Y)
         z_rot = self.pauli_tuple(pauli_gates.Z)
@@ -895,6 +924,10 @@ def _circuit_diagram_info_(
         )
 
     def _value_equality_values_(self):
+        return self._value_equality_values
+
+    @cached_property
+    def _value_equality_values(self):
         return self._clifford_tableau.matrix().tobytes() + self._clifford_tableau.rs.tobytes()
 
     def _value_equality_values_cls_(self):

cirq-core/cirq/qis/clifford_tableau.py

@@ -17,7 +17,7 @@
 import numpy as np
 
 from cirq import protocols
-from cirq._compat import proper_repr
+from cirq._compat import proper_repr, cached_method
 from cirq.qis import quantum_state_representation
 from cirq.value import big_endian_int_to_digits, linear_dict, random_state
 
@@ -652,3 +652,7 @@ def measure(
         self, axes: Sequence[int], seed: 'cirq.RANDOM_STATE_OR_SEED_LIKE' = None
     ) -> List[int]:
         return [self._measure(axis, random_state.parse_random_state(seed)) for axis in axes]
+
+    @cached_method
+    def __hash__(self) -> int:
+        return hash(self.matrix().tobytes() + self.rs.tobytes())