New QEM benchmarking method: "rotated" RB circuits

docs/source/apidoc.md

@@ -23,6 +23,13 @@
    :members:
 ```
 
+### Rotated Randomized Benchmarking Circuits
+
+```{eval-rst}
+.. automodule:: mitiq.benchmarks.rotated_randomized_benchmarking
+   :members:
+```
+
 ### GHZ Circuits
 ```{eval-rst}
 .. automodule:: mitiq.benchmarks.ghz_circuits

docs/source/conf.py

@@ -161,6 +161,9 @@
     r"https://dl.acm.org/doi/10.1145/3466752.3480059",
     r"https://doi.org/10.1145/3466752.3480059",
     r"https://doi.org/10.7566/jpsj.90.032001",
+    r"https://arxiv.org/abs/1612.02058",
+    r"https://arxiv.org/abs/1805.04492",
+    r"https://arxiv.org/abs/1807.05572",
 ]
 
 

mitiq/benchmarks/__init__.py

@@ -4,6 +4,9 @@
 # LICENSE file in the root directory of this source tree.
 
 from mitiq.benchmarks.randomized_benchmarking import generate_rb_circuits
+from mitiq.benchmarks.rotated_randomized_benchmarking import (
+    generate_rotated_rb_circuits,
+)
 from mitiq.benchmarks.mirror_circuits import generate_mirror_circuit
 from mitiq.benchmarks.ghz_circuits import generate_ghz_circuit
 from mitiq.benchmarks.quantum_volume_circuits import (

mitiq/benchmarks/rotated_randomized_benchmarking.py

@@ -0,0 +1,90 @@
+# Copyright (C) Unitary Fund
+#
+# This source code is licensed under the GPL license (v3) found in the
+# LICENSE file in the root directory of this source tree.
+
+"""Functions for generating rotated randomized benchmarking circuits."""
+from typing import List, Optional, cast
+
+import cirq
+import numpy as np
+
+from mitiq import QPROGRAM
+from mitiq.benchmarks import generate_rb_circuits
+from mitiq.interface import convert_from_mitiq
+
+
+def generate_rotated_rb_circuits(
+    n_qubits: int,
+    num_cliffords: int,
+    theta: Optional[float] = None,
+    trials: int = 1,
+    return_type: Optional[str] = None,
+) -> List[QPROGRAM]:
+    r"""
+    Generates a list of "rotated" randomized benchmarking circuits.
+    This benchmarking method enables testing QEM techniques in more general
+    scenarios, closer to real-world applications in which expectation values
+    can take arbitrary values.
+
+    Rotated randomized bencmarking circuits are randomized benchmarking
+    circuits with an :math:`R_z(\theta)` rotation inserted in the middle, such
+    that:
+
+    .. math::
+        C(\theta) =    G_n \dots G_{n/2 +1} R_z(\theta)G_{n/2} \dots G_2 G_1
+
+    where :math:`G_j` are Clifford elements or Clifford gates.
+
+    The circuits generate expectation values which are sinusoidal functions of
+    :math:`\theta`, which in the ideal (noiseless) case vary in a continuous
+    interval of :math:`E_{\rm ideal} \in  [-1, 1]`.
+
+    Since (up to factors of 2) we have
+    :math:`R_z(\theta) =cos(\theta) I +  i \ sin(\theta) Z`, the rotated
+    Clifford circuit :math:`C(\theta)` can be written as a linear combination
+    of just two Clifford circuits, and therefore it is still easy to
+    classically simulate.
+
+    Args:
+        n_qubits: The number of qubits. Can be either 1 or 2.
+        num_cliffords: The number of Clifford group elements in the
+            random circuits. This is proportional to the depth per circuit.
+        theta: The rotation angle about the :math:`Z` axis.
+        trials: The number of random circuits to return.
+        return_type: String which specifies the type of the
+            returned circuits. See the keys of
+            ``mitiq.SUPPORTED_PROGRAM_TYPES`` for options. If ``None``, the
+            returned circuits have type ``cirq.Circuit``.
+
+
+    Returns:
+        A list of rotated randomized benchmarking circuits.
+    """
+
+    circuits = cast(
+        List[cirq.Circuit],
+        generate_rb_circuits(n_qubits, num_cliffords, trials),
+    )
+
+    if theta is None:
+        theta = 2 * np.pi * np.random.rand()
+
+    for circ in circuits:
+        qubits = list(circ.all_qubits())
+        circ.insert(len(circ) // 2, cirq.Rz(rads=theta).on(qubits[0]))
+        if (
+            theta != 0
+            and cirq.DensityMatrixSimulator()
+            .simulate(circ)
+            .final_density_matrix[0, 0]
+            .real
+            == 1
+        ):
+            circ = cast(
+                cirq.Circuit,
+                generate_rotated_rb_circuits(n_qubits, num_cliffords)[0],
+            )
+
+    return_type = "cirq" if not return_type else return_type
+    return [convert_from_mitiq(circuit, return_type) for circuit in circuits]

mitiq/benchmarks/tests/test_rotated_randomized_benchmarking.py

@@ -0,0 +1,57 @@
+# Copyright (C) Unitary Fund
+#
+# This source code is licensed under the GPL license (v3) found in the
+# LICENSE file in the root directory of this source tree.
+
+"""Tests for rotated randomized benchmarking circuits."""
+
+import cirq
+import numpy as np
+import pytest
+
+from mitiq import SUPPORTED_PROGRAM_TYPES
+from mitiq.benchmarks.rotated_randomized_benchmarking import (
+    generate_rotated_rb_circuits,
+)
+
+
+@pytest.mark.parametrize("n_qubits", (1, 2))
+@pytest.mark.parametrize("theta", [0, np.pi / 2, np.pi])
+def test_rotated_rb_circuits(n_qubits, theta):
+    depth = 10
+    for trials in [2, 3]:
+        circuits = generate_rotated_rb_circuits(
+            n_qubits=n_qubits, num_cliffords=depth, theta=theta, trials=trials
+        )
+        for circ in circuits:
+            zero_prob = (
+                cirq.DensityMatrixSimulator()
+                .simulate(circ)
+                .final_density_matrix[0, 0]
+                .real
+            )
+            assert -1.0001 <= zero_prob <= 1.0001
+
+
+@pytest.mark.parametrize("n_qubits", (1, 2))
+@pytest.mark.parametrize("theta", [0, np.pi / 2, np.pi])
+@pytest.mark.parametrize("return_type", SUPPORTED_PROGRAM_TYPES.keys())
+def test_rotated_rb_conversion(n_qubits, theta, return_type):
+    depth = 10
+    for trials in [2, 3]:
+        circuits = generate_rotated_rb_circuits(
+            n_qubits=n_qubits,
+            num_cliffords=depth,
+            theta=theta,
+            trials=trials,
+            return_type=return_type,
+        )
+        for circ in circuits:
+            assert return_type in circ.__module__
+
+
+def test_rotated_rb_circuit_no_theta():
+    circuit = generate_rotated_rb_circuits(n_qubits=1, num_cliffords=5)[0]
+    assert (
+        len(list(circuit.findall_operations_with_gate_type(cirq.ops.Rz))) > 0
+    )