Combine RB and XEB to compute inferred errors

cirq-core/cirq/experiments/__init__.py

@@ -65,4 +65,8 @@
 
 from cirq.experiments.xeb_fitting import XEBPhasedFSimCharacterizationOptions
 
-from cirq.experiments.two_qubit_xeb import TwoQubitXEBResult, parallel_two_qubit_xeb
+from cirq.experiments.two_qubit_xeb import (
+    TwoQubitXEBResult,
+    CombinedXEBRBResult,
+    parallel_two_qubit_xeb,
+)

cirq-core/cirq/experiments/two_qubit_xeb.py

@@ -11,7 +11,7 @@
 # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 # See the License for the specific language governing permissions and
 # limitations under the License.
-from typing import Sequence, TYPE_CHECKING, Optional, Tuple, Dict
+from typing import Sequence, TYPE_CHECKING, Optional, Tuple, Dict, cast
 
 from dataclasses import dataclass
 import itertools
@@ -22,25 +22,24 @@
 import numpy as np
 import pandas as pd
 
-from cirq import ops, devices, value, vis
+from cirq import ops, value, vis
 from cirq.experiments.xeb_sampling import sample_2q_xeb_circuits
 from cirq.experiments.xeb_fitting import benchmark_2q_xeb_fidelities
 from cirq.experiments.xeb_fitting import fit_exponential_decays, exponential_decay
 from cirq.experiments import random_quantum_circuit_generation as rqcg
+from cirq.experiments.qubit_characterizations import ParallelRandomizedBenchmarkingResult
+from cirq.qis import noise_utils
+from cirq._compat import cached_method
 
 if TYPE_CHECKING:
     import cirq
 
 
-def _grid_qubits_for_sampler(sampler: 'cirq.Sampler'):
+def _grid_qubits_for_sampler(sampler: 'cirq.Sampler') -> Optional[Sequence['cirq.GridQubit']]:
     if hasattr(sampler, 'processor'):
         device = sampler.processor.get_device()
         return sorted(device.metadata.qubit_set)
-    else:
-        qubits = devices.GridQubit.rect(3, 2, 4, 3)
-        # Delete one qubit from the rectangular arangement to
-        # 1) make it irregular 2) simplify simulation.
-        return qubits[:-1]
+    return None
 
 
 def _manhattan_distance(qubit1: 'cirq.GridQubit', qubit2: 'cirq.GridQubit') -> int:
@@ -65,7 +64,7 @@ def all_qubit_pairs(self) -> Tuple[Tuple['cirq.GridQubit', 'cirq.GridQubit'], ..
         return tuple(sorted(self._qubit_pair_map.keys()))
 
     def plot_heatmap(self, ax: Optional[plt.Axes] = None, **plot_kwargs) -> plt.Axes:
-        """plot the heatmap for xeb error.
+        """plot the heatmap of XEB errors.
 
         Args:
             ax: the plt.Axes to plot on. If not given, a new figure is created,
@@ -75,7 +74,6 @@ def plot_heatmap(self, ax: Optional[plt.Axes] = None, **plot_kwargs) -> plt.Axes
         show_plot = not ax
         if not isinstance(ax, plt.Axes):
             fig, ax = plt.subplots(1, 1, figsize=(8, 8))
-
         heatmap_data: Dict[Tuple['cirq.GridQubit', ...], float] = {
             pair: self.xeb_error(*pair) for pair in self.all_qubit_pairs
         }
@@ -131,10 +129,13 @@ def _record(self, q0, q1) -> pd.Series:
             q0, q1 = q1, q0
         return self.fidelities.iloc[self._qubit_pair_map[(q0, q1)]]
 
+    def xeb_fidelity(self, q0: 'cirq.GridQubit', q1: 'cirq.GridQubit') -> float:
+        """Return the XEB fidelity of a qubit pair."""
+        return self._record(q0, q1).layer_fid
+
     def xeb_error(self, q0: 'cirq.GridQubit', q1: 'cirq.GridQubit') -> float:
         """Return the XEB error of a qubit pair."""
-        p = self._record(q0, q1).layer_fid
-        return 1 - p
+        return 1 - self.xeb_fidelity(q0, q1)
 
     def all_errors(self) -> Dict[Tuple['cirq.GridQubit', 'cirq.GridQubit'], float]:
         """Return the XEB error of all qubit pairs."""
@@ -156,9 +157,109 @@ def plot_histogram(self, ax: Optional[plt.Axes] = None, **plot_kwargs) -> plt.Ax
             fig.show(**plot_kwargs)
         return ax
 
+    @cached_method
+    def pauli_error(self) -> Dict[Tuple['cirq.GridQubit', 'cirq.GridQubit'], float]:
+        """Return the Pauli error of a qubit pair."""
+        return {
+            pair: noise_utils.decay_constant_to_pauli_error(
+                noise_utils.xeb_fidelity_to_decay_constant(self.xeb_fidelity(*pair), num_qubits=2),
+                num_qubits=2,
+            )
+            for pair in self.all_qubit_pairs
+        }
+
+
+@dataclass(frozen=True)
+class CombinedXEBRBResult:
+    rb_result: ParallelRandomizedBenchmarkingResult
+    xeb_result: TwoQubitXEBResult
+
+    @property
+    def all_qubit_pairs(self) -> Sequence[Tuple['cirq.GridQubit', 'cirq.GridQubit']]:
+        return self.xeb_result.all_qubit_pairs
+
+    def _pauli_error(self, q0: 'cirq.GridQubit', q1: 'cirq.GridQubit') -> float:
+        """Return the total pauli error."""
+        q0, q1 = sorted([q0, q1])
+        single_q_paulis = self.rb_result.pauli_error()
+        return self.xeb_result.pauli_error()[(q0, q1)] + single_q_paulis[q0] + single_q_paulis[q1]
+
+    @cached_method
+    def pauli_error(self) -> Dict[Tuple['cirq.GridQubit', 'cirq.GridQubit'], float]:
+        """Return the equivalent pauli error for all pairs."""
+        return {pair: self._pauli_error(*pair) for pair in self.all_qubit_pairs}
+
+    @cached_method
+    def decay_constant(self) -> Dict[Tuple['cirq.GridQubit', 'cirq.GridQubit'], float]:
+        """Return the equivalent decay constant for all pairs."""
+        return {
+            pair: noise_utils.pauli_error_to_decay_constant(pauli, 2)
+            for pair, pauli in self.pauli_error().items()
+        }
+
+    @cached_method
+    def xeb_error(self) -> Dict[Tuple['cirq.GridQubit', 'cirq.GridQubit'], float]:
+        """Return the equivalent XEB error for all pairs."""
+        return {
+            pair: 1 - noise_utils.decay_constant_to_xeb_fidelity(decay, 2)
+            for pair, decay in self.decay_constant().items()
+        }
+
+    def _target_errors(self, target_error: str) -> Dict[Tuple['cirq.GridQubit', ...], float]:
+        error_funcs = {
+            'pauli': self.pauli_error,
+            'decay_constant': self.decay_constant,
+            'xeb': self.xeb_error,
+        }
+        return cast(Dict[Tuple['cirq.GridQubit', ...], float], error_funcs[target_error]())
+
+    def plot_heatmap(
+        self, target_error: str = 'pauli', ax: Optional[plt.Axes] = None, **plot_kwargs
+    ) -> plt.Axes:
+        """plot the heatmap of the target errors.
+
+        Args:
+            target_error: The error to draw. Must be one of 'xeb', 'pauli', or 'decay_constant'
+            ax: the plt.Axes to plot on. If not given, a new figure is created,
+                plotted on, and shown.
+            **plot_kwargs: Arguments to be passed to 'plt.Axes.plot'.
+        """
+        show_plot = not ax
+        if not isinstance(ax, plt.Axes):
+            fig, ax = plt.subplots(1, 1, figsize=(8, 8))
+        heatmap_data = self._target_errors(target_error)
+
+        name = f'{target_error} error' if target_error != 'decay_constant' else 'decay constant'
+        ax.set_title(f'device {name} heatmap')
+
+        vis.TwoQubitInteractionHeatmap(heatmap_data).plot(ax=ax, **plot_kwargs)
+        if show_plot:
+            fig.show()
+        return ax
+
+    def plot_histogram(
+        self, target_error: str = 'pauli', ax: Optional[plt.Axes] = None, **plot_kwargs
+    ) -> plt.Axes:
+        """plot a histogram of target error.
+
+        Args:
+            target_error: The error to draw. Must be one of 'xeb', 'pauli', or 'decay_constant'
+            ax: the plt.Axes to plot on. If not given, a new figure is created,
+                plotted on, and shown.
+            **plot_kwargs: Arguments to be passed to 'plt.Axes.plot'.
+        """
+        fig = None
+        if ax is None:
+            fig, ax = plt.subplots(1, 1, figsize=(8, 8))
+        vis.integrated_histogram(data=self._target_errors(target_error), ax=ax, **plot_kwargs)
+        if fig is not None:
+            fig.show(**plot_kwargs)
+        return ax
+
 
 def parallel_two_qubit_xeb(
     sampler: 'cirq.Sampler',
+    qubits: Optional[Sequence['cirq.GridQubit']] = None,
     entangling_gate: 'cirq.Gate' = ops.CZ,
     n_repetitions: int = 10**4,
     n_combinations: int = 10,
@@ -172,6 +273,7 @@ def parallel_two_qubit_xeb(
 
     Args:
         sampler: The quantum engine or simulator to run the circuits.
+        qubits: Qubits under test. If none, uses all qubits on the sampler's device.
         entangling_gate: The entangling gate to use.
         n_repetitions: The number of repetitions to use.
         n_combinations: The number of combinations to generate.
@@ -184,10 +286,17 @@ def parallel_two_qubit_xeb(
 
     Returns:
         A TwoQubitXEBResult object representing the results of the experiment.
+
+    Raises:
+        ValueError: If qubits are not specified and the sampler has no device.
     """
     rs = value.parse_random_state(random_state)
 
-    qubits = _grid_qubits_for_sampler(sampler)
+    if qubits is None:
+        qubits = _grid_qubits_for_sampler(sampler)
+        if qubits is None:
+            raise ValueError("Couldn't determine qubits from sampler. Please specify them.")
+
     graph = nx.Graph(
         pair for pair in itertools.combinations(qubits, 2) if _manhattan_distance(*pair) == 1
     )