From 55d9917e60caa268a06b75df909d2d0d93663a15 Mon Sep 17 00:00:00 2001 From: Toshinari Itoko Date: Fri, 17 Nov 2023 23:46:39 +0900 Subject: [PATCH] WIP: Add layer fidelity experiment --- .../randomized_benchmarking/clifford_utils.py | 11 + .../randomized_benchmarking/layer_fidelity.py | 347 ++++++++++++++++++ .../layer_fidelity_analysis.py | 202 ++++++++++ 3 files changed, 560 insertions(+) create mode 100644 qiskit_experiments/library/randomized_benchmarking/layer_fidelity.py create mode 100644 qiskit_experiments/library/randomized_benchmarking/layer_fidelity_analysis.py diff --git a/qiskit_experiments/library/randomized_benchmarking/clifford_utils.py b/qiskit_experiments/library/randomized_benchmarking/clifford_utils.py index f6ab757b5a..570d3b2f60 100644 --- a/qiskit_experiments/library/randomized_benchmarking/clifford_utils.py +++ b/qiskit_experiments/library/randomized_benchmarking/clifford_utils.py @@ -600,3 +600,14 @@ def _layer_indices_from_num(num: Integral) -> Tuple[Integral, Integral, Integral idx1 = num % _NUM_LAYER_1 idx0 = num // _NUM_LAYER_1 return idx0, idx1, idx2 + + +@lru_cache(maxsize=24 * 24) +def _product_1q_nums(first: Integral, second: Integral) -> Integral: + """Return the 2-qubit Clifford integer that represents the product of 1-qubit Cliffords.""" + qc0 = CliffordUtils.clifford_1_qubit_circuit(first) + qc1 = CliffordUtils.clifford_1_qubit_circuit(second) + qc = QuantumCircuit(2) + qc.compose(qc0, qubits=(0,), inplace=True) + qc.compose(qc1, qubits=(1,), inplace=True) + return num_from_2q_circuit(qc) diff --git a/qiskit_experiments/library/randomized_benchmarking/layer_fidelity.py b/qiskit_experiments/library/randomized_benchmarking/layer_fidelity.py new file mode 100644 index 0000000000..f1ad7813a5 --- /dev/null +++ b/qiskit_experiments/library/randomized_benchmarking/layer_fidelity.py @@ -0,0 +1,347 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2023. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. +""" +Layer Fidelity RB Experiment class. +""" +import logging +from collections import defaultdict +from typing import Union, Iterable, Optional, List, Sequence, Tuple + +import numpy as np +from numpy.random import Generator, default_rng +from numpy.random.bit_generator import BitGenerator, SeedSequence + +from qiskit.circuit import QuantumCircuit, CircuitInstruction, Barrier +from qiskit.circuit.library import get_standard_gate_name_mapping +from qiskit.exceptions import QiskitError +from qiskit.providers import BackendV2Converter +from qiskit.providers.backend import Backend, BackendV1, BackendV2 +from qiskit.pulse.instruction_schedule_map import CalibrationPublisher + +from qiskit_experiments.framework import BaseExperiment, Options +from qiskit_experiments.framework.restless_mixin import RestlessMixin + +from .clifford_utils import ( + CliffordUtils, + compose_1q, + compose_2q, + inverse_1q, + inverse_2q, + _product_1q_nums, + _num_from_2q_gate, + _clifford_1q_int_to_instruction, + _clifford_2q_int_to_instruction, + _decompose_clifford_ops, +) +from .layer_fidelity_analysis import LayerFidelityAnalysis + +LOG = logging.getLogger(__name__) + + +GATE_NAME_MAP = get_standard_gate_name_mapping() +NUM_1Q_CLIFFORD = CliffordUtils.NUM_CLIFFORD_1_QUBIT + + +class LayerFidelity(BaseExperiment, RestlessMixin): + """TODO + + # section: overview + + TODO + + # section: analysis_ref + :class:`LayerFidelityAnalysis` + + # section: reference + .. ref_arxiv:: 1 2311.05933 + """ + + def __init__( + self, + physical_qubits: Sequence[int], + two_qubit_layers: Sequence[Sequence[Tuple[int, int]]], + lengths: Iterable[int], + backend: Optional[Backend] = None, + num_samples: int = 3, + seed: Optional[Union[int, SeedSequence, BitGenerator, Generator]] = None, + # full_sampling: Optional[bool] = False, TODO: can we always do full_sampling and remove the option? + two_qubit_gate: Optional[str] = None, + one_qubit_basis_gates: Optional[Sequence[str]] = None, + ): + """Initialize a standard randomized benchmarking experiment. + + Args: + physical_qubits: List of physical qubits for the experiment. + two_qubit_layers: List of pairs of qubits to run on, will use the direction given here. + lengths: A list of RB sequences lengths. + backend: The backend to run the experiment on. + num_samples: Number of samples to generate for each sequence length. + seed: Optional, seed used to initialize ``numpy.random.default_rng``. + when generating circuits. The ``default_rng`` will be initialized + with this seed value every time :meth:`circuits` is called. + two_qubit_gate: Two-qubit gate name (e.g. "cx", "cz", "ecr") of which the two qubit layers consist. + one_qubit_basis_gates: One-qubit gates to use for implementing 1q Clifford operations. + + Raises: + QiskitError: If any invalid argument is supplied. + """ + # Compute full layers + full_layers = [] + for two_q_layer in two_qubit_layers: + qubits_in_layer = {q for qpair in two_q_layer for q in qpair} + layer = two_q_layer + [q for q in physical_qubits if q not in qubits_in_layer] + full_layers.append(layer) + + # Initialize base experiment + super().__init__(physical_qubits, analysis=LayerFidelityAnalysis(full_layers), backend=backend) + + # Verify parameters + # TODO more checks + if len(set(lengths)) != len(lengths): + raise QiskitError( + f"The lengths list {lengths} should not contain " "duplicate elements." + ) + if num_samples <= 0: + raise QiskitError(f"The number of samples {num_samples} should " "be positive.") + if two_qubit_gate not in GATE_NAME_MAP: + pass # TODO: too restrictive to forbidden custom two qubit gate name? + + # Get parameters from backend + if two_qubit_gate is None: + # TODO: implement and raise an error if backend is None + raise NotImplemented() + if one_qubit_basis_gates is None: + # TODO: implement and raise an error if backend is None + raise NotImplemented() + + # Set configurable options + self.set_experiment_options( + lengths=sorted(lengths), + num_samples=num_samples, + seed=seed, + two_qubit_layers=two_qubit_layers, + two_qubit_gate=two_qubit_gate, + one_qubit_basis_gates=tuple(one_qubit_basis_gates), + ) + # self.analysis.set_options(outcome="0" * self.num_qubits) + + @classmethod + def _default_experiment_options(cls) -> Options: + """Default experiment options. + + Experiment Options: + two_qubit_layers (List[List[Tuple[int, int]]]): List of pairs of qubits to run on. + lengths (List[int]): A list of RB sequences lengths. + num_samples (int): Number of samples to generate for each sequence length. + seed (None or int or SeedSequence or BitGenerator or Generator): A seed + used to initialize ``numpy.random.default_rng`` when generating circuits. + The ``default_rng`` will be initialized with this seed value every time + :meth:`circuits` is called. + two_qubit_gate (str): Two-qubit gate name (e.g. "cx", "cz", "ecr") of which the two qubit layers consist. + one_qubit_basis_gates (Tuple[str]): One-qubit gates to use for implementing 1q Clifford operations. + """ + options = super()._default_experiment_options() + options.update_options( + lengths=None, + num_samples=None, + seed=None, + two_qubit_layers=None, + two_qubit_gate=None, + one_qubit_basis_gates=tuple(), + ) + return options + + def set_experiment_options(self, **fields): + """Set the experiment options. + + Args: + fields: The fields to update the options + + Raises: + AttributeError: If the field passed in is not a supported options + """ + for field in {"two_qubit_layers"}: + if hasattr(self._experiment_options, field) and self._experiment_options[field] is not None: + raise AttributeError( + f"Options field {field} is not allowed to update." + ) + super().set_experiment_options(**fields) + + @classmethod + def _default_transpile_options(cls) -> Options: + """Default transpiler options for transpiling RB circuits.""" + return Options(optimization_level=1) + + def set_transpile_options(self, **fields): + """Transpile options is not supported for LayerFidelity experiments. + + Raises: + QiskitError: If `set_transpile_options` is called. + """ + raise QiskitError( + "Custom transpile options is not supported for LayerFidelity experiments." + ) + + def _set_backend(self, backend: Backend): + """Set the backend V2 for RB experiments since RB experiments only support BackendV2 + except for simulators. If BackendV1 is provided, it is converted to V2 and stored. + """ + if isinstance(backend, BackendV1) and "simulator" not in backend.name(): + super()._set_backend(BackendV2Converter(backend, add_delay=True)) + else: + super()._set_backend(backend) + + def __residual_qubits(self, two_qubit_layer): + qubits_in_layer = {q for qpair in two_qubit_layer for q in qpair} + return [q for q in self.physical_qubits if q not in qubits_in_layer] + + def circuits(self) -> List[QuantumCircuit]: + """Return a list of physical circuits to measure layer fidelity. + + Returns: + A list of :class:`QuantumCircuit`. + """ + opts = self.experiment_options + rng = default_rng(seed=opts.seed) + basis_gates = (opts.two_qubit_gate,) + opts.one_qubit_basis_gates + GATE2Q = GATE_NAME_MAP[opts.two_qubit_gate] + GATE2Q_CLIFF = _num_from_2q_gate(GATE2Q) + residal_qubits_by_layer = [self.__residual_qubits(layer) for layer in opts.two_qubit_layers] + # Circuit generation + circuits = [] + num_qubits = max(self.physical_qubits) + 1 + for i_sample in range(opts.num_samples): + for i_set, (two_qubit_layer, one_qubits) in enumerate(zip(opts.two_qubit_layers, residal_qubits_by_layer)): + num_2q_gates = len(two_qubit_layer) + num_1q_gates = len(one_qubits) + composite_qubits = two_qubit_layer + [(q,) for q in one_qubits] + composite_clbits = [(2*c, 2*c+1) for c in range(num_2q_gates)] + composite_clbits.extend([(c,) for c in range(2*num_2q_gates, 2*num_2q_gates+num_1q_gates)]) + for length in opts.lengths: + # initialize cliffords and a ciruit (0: identity clifford) + cliffs_2q = [0] * num_2q_gates + cliffs_1q = [0] * num_1q_gates + circ = QuantumCircuit(num_qubits) + for _ in range(length): + # sample random 1q-Clifford layer + for j, qpair in enumerate(two_qubit_layer): + # sample product of two 1q-Cliffords as 2q interger Clifford + samples = rng.integers(NUM_1Q_CLIFFORD, size=2) + cliffs_2q[j] = compose_2q(cliffs_2q[j], _product_1q_nums(*samples)) + for sample, q in zip(samples, qpair): + circ._append( + _clifford_1q_int_to_instruction( + sample, opts.one_qubit_basis_gates + ), + (circ.qubits[q],), + tuple(), + ) + for k, q in enumerate(one_qubits): + sample = rng.integers(NUM_1Q_CLIFFORD) + cliffs_1q[k] = compose_1q(cliffs_1q[k], sample) + circ._append( + _clifford_1q_int_to_instruction(sample, opts.one_qubit_basis_gates), + (circ.qubits[q],), + tuple(), + ) + circ.barrier(self.physical_qubits) + # add two qubit gates + for j, qpair in enumerate(two_qubit_layer): + circ._append(GATE2Q, tuple(circ.qubits[q] for q in qpair), tuple()) + cliffs_2q[j] = compose_2q(cliffs_2q[j], GATE2Q_CLIFF) + # TODO: add dd if necessary + for k, q in enumerate(one_qubits): + # TODO: add dd if necessary + pass + circ.barrier(self.physical_qubits) + # add the last inverse + for j, qpair in enumerate(two_qubit_layer): + inv = inverse_2q(cliffs_2q[j]) + circ._append( + _clifford_2q_int_to_instruction(inv, basis_gates), + tuple(circ.qubits[q] for q in qpair), + tuple(), + ) + for k, q in enumerate(one_qubits): + inv = inverse_1q(cliffs_1q[k]) + circ._append( + _clifford_1q_int_to_instruction(inv, opts.one_qubit_basis_gates), + (circ.qubits[q],), + tuple(), + ) + + circ.measure_active() # includes insertion of the barrier before measurement + # store composite structure in metadata + circ.metadata = { + 'experiment_type': 'BatchExperiment', 'composite_metadata': [ + { + 'experiment_type': 'ParallelExperiment', + 'composite_index': list(range(num_2q_gates + num_1q_gates)), + 'composite_metadata': [ + {'experiment_type': 'SubLayerFidelity', 'physical_qubits': qpair, 'sample': i_sample, 'xval': length} + for qpair in two_qubit_layer + ] + [ + {'experiment_type': 'SubLayerFidelity', 'physical_qubits': (q,), 'sample': i_sample, 'xval': length} + for q in one_qubits + ], + 'composite_qubits': composite_qubits, + 'composite_clbits': composite_clbits + } + ], + 'composite_index': [i_set] + } + circuits.append(circ) + + return circuits + + def _transpiled_circuits(self) -> List[QuantumCircuit]: + """Return a list of experiment circuits, transpiled.""" + transpiled = [_decompose_clifford_ops(circ) for circ in self.circuits()] + # Set custom calibrations provided in backend + if isinstance(self.backend, BackendV2): + instructions = [] # (op_name, qargs) for each element where qargs means qubit tuple + for two_qubit_layer in self.experiment_options.two_qubit_layers: + for qpair in two_qubit_layer: + instructions.append((self.experiment_options.two_qubit_gate, tuple(qpair))) + for q in self.__residual_qubits(two_qubit_layer): + for gate_1q in self.experiment_options.one_qubit_basis_gates: + instructions.append((gate_1q, (q,))) + + common_calibrations = defaultdict(dict) + for op_name, qargs in instructions: + inst_prop = self.backend.target[op_name].get(qargs, None) + if inst_prop is None: + continue + schedule = inst_prop.calibration + if schedule is None: + continue + publisher = schedule.metadata.get("publisher", CalibrationPublisher.QISKIT) + if publisher != CalibrationPublisher.BACKEND_PROVIDER: + common_calibrations[op_name][(qargs, tuple())] = schedule + + for circ in transpiled: + # This logic is inefficient in terms of payload size and backend compilation + # because this binds every custom pulse to a circuit regardless of + # its existence. It works but redundant calibration must be removed -- NK. + circ.calibrations = common_calibrations + + return transpiled + + def _metadata(self): + metadata = super()._metadata() + # Store measurement level and meas return if they have been + # set for the experiment + for run_opt in ["meas_level", "meas_return"]: + if hasattr(self.run_options, run_opt): + metadata[run_opt] = getattr(self.run_options, run_opt) + + return metadata diff --git a/qiskit_experiments/library/randomized_benchmarking/layer_fidelity_analysis.py b/qiskit_experiments/library/randomized_benchmarking/layer_fidelity_analysis.py new file mode 100644 index 0000000000..fa2010582d --- /dev/null +++ b/qiskit_experiments/library/randomized_benchmarking/layer_fidelity_analysis.py @@ -0,0 +1,202 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2023. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. +""" +Layer Fidelity RB analysis class. +""" +from typing import List, Tuple, Union + +import lmfit + +import qiskit_experiments.curve_analysis as curve +from qiskit_experiments.exceptions import AnalysisError +from qiskit_experiments.framework import CompositeAnalysis, AnalysisResultData, ExperimentData + + +class _SubLayerFidelityAnalysis(curve.CurveAnalysis): + r"""A class to analyze a sub-experiment for estimating layer fidelity, + i.e. one of 1Q/2Q simultaneous direct RBs. + + # section: overview + This analysis takes only single series. + This series is fit by the exponential decay function. + From the fit :math:`\alpha` value this analysis estimates the error per gate (EPG). + + # section: fit_model + .. math:: + + F(x) = a \alpha^x + b + + # section: fit_parameters + defpar a: + desc: Height of decay curve. + init_guess: Determined by :math:`1 - b`. + bounds: [0, 1] + defpar b: + desc: Base line. + init_guess: Determined by :math:`(1/2)^n` where :math:`n` is number of qubit. + bounds: [0, 1] + defpar \alpha: + desc: Depolarizing parameter. + init_guess: Determined by :func:`~.guess.rb_decay`. + bounds: [0, 1] + + # section: reference + .. ref_arxiv:: 1 2311.05933 + """ + + def __init__(self, physical_qubits): + super().__init__( + models=[ + lmfit.models.ExpressionModel( + expr="a * alpha ** x + b", + name="rb_decay", + ) + ] + ) + self._physical_qubits = physical_qubits + + @classmethod + def _default_options(cls): + """Default analysis options. + """ + default_options = super()._default_options() + default_options.plotter.set_figure_options( + xlabel="Layer Length", + ylabel="P(0)", + ) + default_options.plot_raw_data = True + default_options.result_parameters = ["alpha"] + default_options.average_method = "sample" + + return default_options + + def _generate_fit_guesses( + self, + user_opt: curve.FitOptions, + curve_data: curve.ScatterTable, + ) -> Union[curve.FitOptions, List[curve.FitOptions]]: + """Create algorithmic initial fit guess from analysis options and curve data. + + Args: + user_opt: Fit options filled with user provided guess and bounds. + curve_data: Formatted data collection to fit. + + Returns: + List of fit options that are passed to the fitter function. + """ + user_opt.bounds.set_if_empty( + a=(0, 1), + alpha=(0, 1), + b=(0, 1), + ) + + b_guess = 1 / 2 ** len(self._physical_qubits) + alpha_guess = curve.guess.rb_decay(curve_data.x, curve_data.y, b=b_guess) + a_guess = (curve_data.y[0] - b_guess) / (alpha_guess ** curve_data.x[0]) + + user_opt.p0.set_if_empty( + b=b_guess, + a=a_guess, + alpha=alpha_guess, + ) + + return user_opt + + def _create_analysis_results( + self, + fit_data: curve.CurveFitResult, + quality: str, + **metadata, + ) -> List[AnalysisResultData]: + """Create analysis results for important fit parameters. + + Args: + fit_data: Fit outcome. + quality: Quality of fit outcome. + + Returns: + List of analysis result data. + """ + outcomes = super()._create_analysis_results(fit_data, quality, **metadata) + num_qubits = len(self._physical_qubits) + + # Calculate EPC + alpha = fit_data.ufloat_params["alpha"] + scale = (2**num_qubits - 1) / (2**num_qubits) + epg = scale * (1 - alpha) + + outcomes.append( + AnalysisResultData( + name="EPG", + value=epg, + chisq=fit_data.reduced_chisq, + quality=quality, + extra=metadata, + ) + ) + return outcomes + + + +class LayerFidelityAnalysis(CompositeAnalysis): + r"""A class to analyze layer fidelity experiments. + + # section: see_also + * :py:class:`qiskit_experiments.library.characterization.analysis.SubLayerFidelityAnalysis` + + # section: reference + .. ref_arxiv:: 1 2311.05933 + """ + + def __init__(self, layers, analyses=None): + if analyses: + # TODO: Validation + pass + else: + analyses = [] + for a_layer in layers: + a_layer_analyses = [_SubLayerFidelityAnalysis(qubits) for qubits in a_layer] + analyses.append(CompositeAnalysis(a_layer_analyses, flatten_results=True)) + + super().__init__(analyses, flatten_results=True) + + def _run_analysis( + self, experiment_data: ExperimentData + ) -> Tuple[List[AnalysisResultData], List["matplotlib.figure.Figure"]]: + r"""Run analysis for Layer Fidelity experiment. + It invokes CompositeAnalysis._run_analysis that will invoke + _run_analysis for the sub-experiments (1Q/2Q simultaneous direct RBs). + Based on the results, it computes the result for Layer Fidelity. + """ + + # Run composite analysis and extract sub-experiments results + analysis_results, figures = super()._run_analysis(experiment_data) + + # Calculate Layer Fidelity from EPGs + lf = None # TODO + quality_lf = ( + "good" if all(sub.quality == "good" for sub in analysis_results) else "bad" + ) + lf_result = AnalysisResultData( + name="LF", + value=lf, + chisq=None, + quality=quality_lf, + extra={}, + ) + + # TODO: Plot LF by chain length for a full 2q-gate chain + + # Return combined results + analysis_results = [lf_result] + analysis_results + # figures = [lf_plot] + figures + return analysis_results, figures