Simultaneous readout

cirq-core/cirq/__init__.py

@@ -95,6 +95,7 @@
 )
 
 from cirq.experiments import (
+    estimate_correlated_single_qubit_readout_errors,
     estimate_single_qubit_readout_errors,
     hog_score_xeb_fidelity_from_probabilities,
     least_squares_xeb_fidelity_from_expectations,

cirq-core/cirq/experiments/__init__.py

@@ -72,6 +72,7 @@
 )
 
 from cirq.experiments.single_qubit_readout_calibration import (
+    estimate_correlated_single_qubit_readout_errors,
     estimate_single_qubit_readout_errors,
     SingleQubitReadoutCalibrationResult,
 )

cirq-core/cirq/experiments/single_qubit_readout_calibration.py

@@ -11,21 +11,23 @@
 # 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 Any, Dict, Iterable, TYPE_CHECKING
-
+"""Module for supporting single qubit readout experiments using
+either correlated or uncorrelated readout statistics.
+"""
 import dataclasses
+import random
 import time
+from typing import Any, Dict, Iterable, List, Optional, TYPE_CHECKING
 
+import sympy
 import numpy as np
-
-from cirq import circuits, ops
+from cirq import circuits, ops, study
 
 if TYPE_CHECKING:
     import cirq
 
 
-@dataclasses.dataclass(frozen=True)
+@dataclasses.dataclass
 class SingleQubitReadoutCalibrationResult:
     """Result of estimating single qubit readout error.
 
@@ -96,19 +98,131 @@ def estimate_single_qubit_readout_errors(
         the probabilities. Also stores a timestamp indicating the time when
         data was finished being collected from the sampler.
     """
+    return estimate_correlated_single_qubit_readout_errors(
+        sampler=sampler,
+        qubits=qubits,
+        repetitions=repetitions,
+        trials=2,
+        bit_strings=np.array([[0, 1] for q in qubits]),
+    )
+
+
+def estimate_correlated_single_qubit_readout_errors(
+    sampler: 'cirq.Sampler',
+    *,
+    qubits: Iterable['cirq.Qid'],
+    trials: int = 20,
+    repetitions: int = 1000,
+    trials_per_batch: Optional[int] = None,
+    bit_strings: np.ndarray = None,
+) -> SingleQubitReadoutCalibrationResult:
+    """Estimate single qubit readout error using parallel operations.
+
+    For each trial, prepare and then measure a random computational basis
+    bitstring on qubits using gates in parallel.
+    Returns a SingleQubitReadoutCalibrationResult which can be used to
+    compute readout errors for each qubit.
+
+    Args:
+        sampler: The `cirq.Sampler` used to run the circuits.
+        qubits: The qubits being tested.
+        repetitions: The number of measurement repetitions to perform for
+            each trial.
+        trials: The number of bitstrings to prepare.
+        trials_per_batch:  If provided, split the experiment into batches
+            with this number of trials in each batch.
+        bit_strings: Optional numpy array of shape (qubits, trials) where the
+            first dimension is the qubit (ordered by the qubit order from
+            the qubits parameter) and the second dimension is the number of
+            trials.  Each value should be a 0 or 1 which specifies which
+            state the qubit should be prepared into during that trial.
+            If not provided, the function will generate random bit strings
+            for you.
+
+    Returns:
+        A SingleQubitReadoutCalibrationResult storing the readout error
+        probabilities as well as the number of repetitions used to estimate
+        the probabilities. Also stores a timestamp indicating the time when
+        data was finished being collected from the sampler.  Note that,
+        if there did not exist a trial where a given qubit was set to |0〉,
+        the zero-state error will be set to `nan` (not a number).  Likewise
+        for qubits with no |1〉trial and one-state error.
+    """
     qubits = list(qubits)
 
-    zeros_circuit = circuits.Circuit(ops.measure_each(*qubits, key_func=repr))
-    ones_circuit = circuits.Circuit(
-        ops.X.on_each(*qubits), ops.measure_each(*qubits, key_func=repr)
+    if trials <= 0:
+        raise ValueError("Must provide non-zero trials for readout calibration.")
+    if repetitions <= 0:
+        raise ValueError("Must provide non-zero repetition for readout calibration.")
+    if bit_strings is None:
+        bit_strings = np.array([[random.randint(0, 1) for _ in range(trials)] for _ in qubits])
+    if not hasattr(bit_strings, 'shape') or bit_strings.shape != (len(qubits), trials):
+        raise ValueError(
+            'bit_strings must be numpy array '
+            f'of shape (qubits, trials) ({len(qubits)}, {trials}) '
+            f"but was {bit_strings.shape if hasattr(bit_strings, 'shape') else None}"
+        )
+    if trials_per_batch is None:
+        trials_per_batch = trials
+    if trials_per_batch <= 0:
+        raise ValueError("Must provide non-zero trials_per_batch for readout calibration.")
+
+    all_circuits = []
+    all_sweeps: List[study.Sweepable] = []
+    num_batches = (trials + trials_per_batch - 1) // trials_per_batch
+
+    # Initialize circuits
+    flip_symbols = sympy.symbols(f'flip_0:{len(qubits)}')
+    flip_circuit = circuits.Circuit(
+        [ops.X(q) ** s for q, s in zip(qubits, flip_symbols)],
+        [ops.measure_each(*qubits, key_func=repr)],
     )
-
-    zeros_result = sampler.run(zeros_circuit, repetitions=repetitions)
-    ones_result = sampler.run(ones_circuit, repetitions=repetitions)
+    all_circuits = [flip_circuit] * num_batches
+
+    # Initialize sweeps
+    for batch in range(num_batches):
+        single_sweeps = []
+        for qubit_idx, q in enumerate(qubits):
+            trial_range = range(
+                batch * trials_per_batch, min((batch + 1) * trials_per_batch, trials)
+            )
+            single_sweeps.append(
+                study.Points(
+                    key=f'flip_{qubit_idx}',
+                    points=[bit_strings[qubit_idx][bit] for bit in trial_range],
+                )
+            )
+        total_sweeps = study.Zip(*single_sweeps)
+        all_sweeps.append(total_sweeps)
+
+    # Execute circuits
+    results = sampler.run_batch(all_circuits, all_sweeps, repetitions=repetitions)
     timestamp = time.time()
 
-    zero_state_errors = {q: np.mean(zeros_result.measurements[repr(q)]) for q in qubits}
-    one_state_errors = {q: 1 - np.mean(ones_result.measurements[repr(q)]) for q in qubits}
+    # Analyze results
+    zero_state_trials: Dict[cirq.Qid, int] = {q: 0 for q in qubits}
+    one_state_trials: Dict[cirq.Qid, int] = {q: 0 for q in qubits}
+    zero_state_totals: Dict[cirq.Qid, int] = {q: 0 for q in qubits}
+    one_state_totals: Dict[cirq.Qid, int] = {q: 0 for q in qubits}
+    for batch_result in results:
+        for trial_idx, trial_result in enumerate(batch_result):
+            for qubit_idx, q in enumerate(qubits):
+                had_x_gate = bit_strings[qubit_idx][trial_idx]
+                if had_x_gate:
+                    one_state_trials[q] += repetitions - np.sum(trial_result.measurements[repr(q)])
+                    one_state_totals[q] += repetitions
+                else:
+                    zero_state_trials[q] += np.sum(trial_result.measurements[repr(q)])
+                    zero_state_totals[q] += repetitions
+
+    zero_state_errors = {
+        q: zero_state_trials[q] / zero_state_totals[q] if zero_state_totals[q] > 0 else np.nan
+        for q in qubits
+    }
+    one_state_errors = {
+        q: one_state_trials[q] / one_state_totals[q] if one_state_totals[q] > 0 else np.nan
+        for q in qubits
+    }
 
     return SingleQubitReadoutCalibrationResult(
         zero_state_errors=zero_state_errors,

cirq-core/cirq/experiments/single_qubit_readout_calibration_test.py

@@ -11,14 +11,24 @@
 # 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 List
+import pytest
 
 import numpy as np
 
 import cirq
 
 
+def test_single_qubit_readout_result_repr():
+    result = cirq.experiments.SingleQubitReadoutCalibrationResult(
+        zero_state_errors={cirq.LineQubit(0): 0.1},
+        one_state_errors={cirq.LineQubit(0): 0.2},
+        repetitions=1000,
+        timestamp=0.3,
+    )
+    cirq.testing.assert_equivalent_repr(result)
+
+
 class NoisySingleQubitReadoutSampler(cirq.Sampler):
     def __init__(self, p0: float, p1: float, seed: 'cirq.RANDOM_STATE_OR_SEED_LIKE' = None):
         """Sampler that flips some bits upon readout.
@@ -67,7 +77,7 @@ def test_estimate_single_qubit_readout_errors_no_noise():
 def test_estimate_single_qubit_readout_errors_with_noise():
     qubits = cirq.LineQubit.range(5)
     sampler = NoisySingleQubitReadoutSampler(p0=0.1, p1=0.2, seed=1234)
-    repetitions = 1000
+    repetitions = 2000
     result = cirq.estimate_single_qubit_readout_errors(
         sampler, qubits=qubits, repetitions=repetitions
     )
@@ -79,11 +89,121 @@ def test_estimate_single_qubit_readout_errors_with_noise():
     assert isinstance(result.timestamp, float)
 
 
-def test_single_qubit_readout_calibration_result_repr():
-    result = cirq.experiments.SingleQubitReadoutCalibrationResult(
-        zero_state_errors={cirq.LineQubit(0): 0.1},
-        one_state_errors={cirq.LineQubit(0): 0.2},
-        repetitions=1000,
-        timestamp=0.3,
+def test_estimate_correlated_readout_errors_no_noise():
+    qubits = cirq.LineQubit.range(10)
+    sampler = cirq.Simulator()
+    repetitions = 1000
+    result = cirq.estimate_correlated_single_qubit_readout_errors(
+        sampler, qubits=qubits, repetitions=repetitions
     )
-    cirq.testing.assert_equivalent_repr(result)
+    assert result.zero_state_errors == {q: 0 for q in qubits}
+    assert result.one_state_errors == {q: 0 for q in qubits}
+    assert result.repetitions == repetitions
+    assert isinstance(result.timestamp, float)
+
+
+def test_estimate_correlated_readout_errors_all_zeros():
+    qubits = cirq.LineQubit.range(10)
+    sampler = cirq.ZerosSampler()
+    repetitions = 1000
+    result = cirq.estimate_correlated_single_qubit_readout_errors(
+        sampler, qubits=qubits, repetitions=repetitions
+    )
+    assert result.zero_state_errors == {q: 0 for q in qubits}
+    assert result.one_state_errors == {q: 1 for q in qubits}
+    assert result.repetitions == repetitions
+    assert isinstance(result.timestamp, float)
+
+
+def test_estimate_correlated_readout_errors_bad_bit_string():
+    qubits = cirq.LineQubit.range(10)
+    with pytest.raises(ValueError, match='but was None'):
+        _ = cirq.estimate_correlated_single_qubit_readout_errors(
+            cirq.ZerosSampler(),
+            qubits=qubits,
+            repetitions=1000,
+            trials=35,
+            trials_per_batch=10,
+            bit_strings=[1, 1, 1, 1],
+        )
+
+
+def test_estimate_correlated_readout_errors_zero_reps():
+    qubits = cirq.LineQubit.range(10)
+    with pytest.raises(ValueError, match='non-zero repetition'):
+        _ = cirq.estimate_correlated_single_qubit_readout_errors(
+            cirq.ZerosSampler(),
+            qubits=qubits,
+            repetitions=0,
+            trials=35,
+            trials_per_batch=10,
+        )
+
+
+def test_estimate_correlated_readout_errors_zero_trials():
+    qubits = cirq.LineQubit.range(10)
+    with pytest.raises(ValueError, match='non-zero trials'):
+        _ = cirq.estimate_correlated_single_qubit_readout_errors(
+            cirq.ZerosSampler(),
+            qubits=qubits,
+            repetitions=1000,
+            trials=0,
+            trials_per_batch=10,
+        )
+
+
+def test_estimate_correlated_readout_errors_zero_batch():
+    qubits = cirq.LineQubit.range(10)
+    with pytest.raises(ValueError, match='non-zero trials_per_batch'):
+        _ = cirq.estimate_correlated_single_qubit_readout_errors(
+            cirq.ZerosSampler(),
+            qubits=qubits,
+            repetitions=1000,
+            trials=10,
+            trials_per_batch=0,
+        )
+
+
+def test_estimate_correlated_readout_errors_batching():
+    qubits = cirq.LineQubit.range(10)
+    sampler = cirq.ZerosSampler()
+    repetitions = 1000
+    result = cirq.estimate_correlated_single_qubit_readout_errors(
+        sampler, qubits=qubits, repetitions=repetitions, trials=35, trials_per_batch=10
+    )
+    assert result.zero_state_errors == {q: 0 for q in qubits}
+    assert result.one_state_errors == {q: 1 for q in qubits}
+    assert result.repetitions == repetitions
+    assert isinstance(result.timestamp, float)
+
+
+def test_estimate_correlated_readout_errors_with_noise():
+    qubits = cirq.LineQubit.range(5)
+    sampler = NoisySingleQubitReadoutSampler(p0=0.1, p1=0.2, seed=1234)
+    repetitions = 1000
+    result = cirq.estimate_correlated_single_qubit_readout_errors(
+        sampler, qubits=qubits, repetitions=repetitions, trials=40
+    )
+    for error in result.one_state_errors.values():
+        assert 0.17 < error < 0.23
+    for error in result.zero_state_errors.values():
+        assert 0.07 < error < 0.13
+    assert result.repetitions == repetitions
+    assert isinstance(result.timestamp, float)
+
+
+def test_estimate_correlated_readout_errors_missing_qubits():
+    qubits = cirq.LineQubit.range(4)
+
+    result = cirq.estimate_correlated_single_qubit_readout_errors(
+        cirq.ZerosSampler(),
+        qubits=qubits,
+        repetitions=2000,
+        trials=1,
+        bit_strings=np.array([[0], [0], [0], [0]]),
+    )
+    assert result.zero_state_errors == {q: 0 for q in qubits}
+    # Trial did not include a one-state
+    assert all(np.isnan(result.one_state_errors[q]) for q in qubits)
+    assert result.repetitions == 2000
+    assert isinstance(result.timestamp, float)

cirq-core/cirq/protocols/json_test_data/SingleQubitReadoutCalibrationResult.repr

@@ -1 +1 @@
-cirq.experiments.SingleQubitReadoutCalibrationResult(zero_state_errors={cirq.LineQubit(0): 0.1}, one_state_errors={cirq.LineQubit(0): 0.2}, repetitions=1000, timestamp=0.3)
+cirq.experiments.SingleQubitReadoutCalibrationResult(zero_state_errors={cirq.LineQubit(0): 0.1}, one_state_errors={cirq.LineQubit(0): 0.2}, repetitions=1000, timestamp=0.3)