From af2bd56696ffcf54d8b922e314aa1be93101a736 Mon Sep 17 00:00:00 2001 From: Naoki Kanazawa Date: Wed, 17 Jan 2024 22:39:26 +0900 Subject: [PATCH] Rework - reorganize; move all experiments and analyses to own library.driven_freq_tuning - add util; define StarkCoefficients dataclasses and util functions - separation; delegated the role of coefficient manipulation from the analysis class to util functions --- docs/apidocs/index.rst | 1 + docs/apidocs/mod_driven_freq_tuning.rst | 6 + qiskit_experiments/__init__.py | 1 + qiskit_experiments/library/__init__.py | 13 +- .../library/characterization/__init__.py | 11 +- .../characterization/analysis/__init__.py | 4 +- .../analysis/ramsey_xy_analysis.py | 398 ----------- .../characterization/analysis/t1_analysis.py | 269 +------- .../library/characterization/ramsey_xy.py | 628 +----------------- .../library/characterization/t1.py | 349 +--------- .../library/driven_freq_tuning/__init__.py | 66 ++ .../library/driven_freq_tuning/analyses.py | 584 ++++++++++++++++ .../driven_freq_tuning/coefficient_utils.py | 232 +++++++ .../library/driven_freq_tuning/p1_spect.py | 270 ++++++++ .../library/driven_freq_tuning/ramsey.py | 359 ++++++++++ .../driven_freq_tuning/ramsey_amp_scan.py | 311 +++++++++ test/library/driven_freq_tuning/__init__.py | 12 + .../test_stark_p1_spect.py | 263 +++----- .../test_stark_ramsey_xy.py | 66 +- test/library/driven_freq_tuning/test_utils.py | 161 +++++ 20 files changed, 2134 insertions(+), 1870 deletions(-) create mode 100644 docs/apidocs/mod_driven_freq_tuning.rst create mode 100644 qiskit_experiments/library/driven_freq_tuning/__init__.py create mode 100644 qiskit_experiments/library/driven_freq_tuning/analyses.py create mode 100644 qiskit_experiments/library/driven_freq_tuning/coefficient_utils.py create mode 100644 qiskit_experiments/library/driven_freq_tuning/p1_spect.py create mode 100644 qiskit_experiments/library/driven_freq_tuning/ramsey.py create mode 100644 qiskit_experiments/library/driven_freq_tuning/ramsey_amp_scan.py create mode 100644 test/library/driven_freq_tuning/__init__.py rename test/library/{characterization => driven_freq_tuning}/test_stark_p1_spect.py (53%) rename test/library/{characterization => driven_freq_tuning}/test_stark_ramsey_xy.py (84%) create mode 100644 test/library/driven_freq_tuning/test_utils.py diff --git a/docs/apidocs/index.rst b/docs/apidocs/index.rst index 01efc9c174..fac3299b4d 100644 --- a/docs/apidocs/index.rst +++ b/docs/apidocs/index.rst @@ -30,6 +30,7 @@ Experiment Modules mod_calibration mod_characterization + mod_driven_freq_tuning mod_randomized_benchmarking mod_tomography mod_quantum_volume diff --git a/docs/apidocs/mod_driven_freq_tuning.rst b/docs/apidocs/mod_driven_freq_tuning.rst new file mode 100644 index 0000000000..bdc7fa1462 --- /dev/null +++ b/docs/apidocs/mod_driven_freq_tuning.rst @@ -0,0 +1,6 @@ +.. _qiskit-experiments-driven-freq-tuning: + +.. automodule:: qiskit_experiments.library.driven_freq_tuning + :no-members: + :no-inherited-members: + :no-special-members: diff --git a/qiskit_experiments/__init__.py b/qiskit_experiments/__init__.py index e9c613b259..32532928b8 100644 --- a/qiskit_experiments/__init__.py +++ b/qiskit_experiments/__init__.py @@ -49,6 +49,7 @@ - :mod:`qiskit_experiments.library.calibration` - :mod:`qiskit_experiments.library.characterization` +- :mod:`qiskit_experiments.library.driven_freq_tuning` - :mod:`qiskit_experiments.library.randomized_benchmarking` - :mod:`qiskit_experiments.library.tomography` """ diff --git a/qiskit_experiments/library/__init__.py b/qiskit_experiments/library/__init__.py index 6737402863..29770ed7b3 100644 --- a/qiskit_experiments/library/__init__.py +++ b/qiskit_experiments/library/__init__.py @@ -76,8 +76,9 @@ ~characterization.FineXDrag ~characterization.FineSXDrag ~characterization.MultiStateDiscrimination - ~characterization.StarkRamseyXY - ~characterization.StarkRamseyXYAmpScan + ~driven_freq_tuning.StarkRamseyXY + ~driven_freq_tuning.StarkRamseyXYAmpScan + ~driven_freq_tuning.StarkP1Spectroscopy .. _characterization two qubits: @@ -160,7 +161,6 @@ class instance to manage parameters and pulse schedules. ) from .characterization import ( T1, - StarkP1Spectroscopy, T2Hahn, T2Ramsey, Tphi, @@ -187,8 +187,6 @@ class instance to manage parameters and pulse schedules. CorrelatedReadoutError, ZZRamsey, MultiStateDiscrimination, - StarkRamseyXY, - StarkRamseyXYAmpScan, ) from .randomized_benchmarking import StandardRB, InterleavedRB from .tomography import ( @@ -199,6 +197,11 @@ class instance to manage parameters and pulse schedules. MitigatedProcessTomography, ) from .quantum_volume import QuantumVolume +from .driven_freq_tuning import ( + StarkRamseyXY, + StarkRamseyXYAmpScan, + StarkP1Spectroscopy, +) # Experiment Sub-modules from . import calibration diff --git a/qiskit_experiments/library/characterization/__init__.py b/qiskit_experiments/library/characterization/__init__.py index 8aaaceee4c..daa29bb4a4 100644 --- a/qiskit_experiments/library/characterization/__init__.py +++ b/qiskit_experiments/library/characterization/__init__.py @@ -24,7 +24,6 @@ :template: autosummary/experiment.rst T1 - StarkP1Spectroscopy T2Ramsey T2Hahn Tphi @@ -50,8 +49,6 @@ ResonatorSpectroscopy MultiStateDiscrimination ZZRamsey - StarkRamseyXY - StarkRamseyXYAmpScan Analysis @@ -63,7 +60,6 @@ T1Analysis T1KerneledAnalysis - StarkP1SpectAnalysis T2RamseyAnalysis T2HahnAnalysis TphiAnalysis @@ -71,7 +67,6 @@ DragCalAnalysis FineAmplitudeAnalysis RamseyXYAnalysis - StarkRamseyXYAmpScanAnalysis ReadoutAngleAnalysis ResonatorSpectroscopyAnalysis LocalReadoutErrorAnalysis @@ -85,8 +80,6 @@ DragCalAnalysis, FineAmplitudeAnalysis, RamseyXYAnalysis, - StarkRamseyXYAmpScanAnalysis, - StarkP1SpectAnalysis, T2RamseyAnalysis, T1Analysis, T1KerneledAnalysis, @@ -101,7 +94,7 @@ MultiStateDiscriminationAnalysis, ) -from .t1 import T1, StarkP1Spectroscopy +from .t1 import T1 from .qubit_spectroscopy import QubitSpectroscopy from .ef_spectroscopy import EFSpectroscopy from .t2ramsey import T2Ramsey @@ -111,7 +104,7 @@ from .rabi import Rabi, EFRabi from .half_angle import HalfAngle from .fine_amplitude import FineAmplitude, FineXAmplitude, FineSXAmplitude, FineZXAmplitude -from .ramsey_xy import RamseyXY, StarkRamseyXY, StarkRamseyXYAmpScan +from .ramsey_xy import RamseyXY from .fine_frequency import FineFrequency from .drag import RoughDrag from .readout_angle import ReadoutAngle diff --git a/qiskit_experiments/library/characterization/analysis/__init__.py b/qiskit_experiments/library/characterization/analysis/__init__.py index 8520060772..f249dcd9be 100644 --- a/qiskit_experiments/library/characterization/analysis/__init__.py +++ b/qiskit_experiments/library/characterization/analysis/__init__.py @@ -14,10 +14,10 @@ from .drag_analysis import DragCalAnalysis from .fine_amplitude_analysis import FineAmplitudeAnalysis -from .ramsey_xy_analysis import RamseyXYAnalysis, StarkRamseyXYAmpScanAnalysis +from .ramsey_xy_analysis import RamseyXYAnalysis from .t2ramsey_analysis import T2RamseyAnalysis from .t2hahn_analysis import T2HahnAnalysis -from .t1_analysis import T1Analysis, T1KerneledAnalysis, StarkP1SpectAnalysis +from .t1_analysis import T1Analysis, T1KerneledAnalysis from .tphi_analysis import TphiAnalysis from .cr_hamiltonian_analysis import CrossResonanceHamiltonianAnalysis from .readout_angle_analysis import ReadoutAngleAnalysis diff --git a/qiskit_experiments/library/characterization/analysis/ramsey_xy_analysis.py b/qiskit_experiments/library/characterization/analysis/ramsey_xy_analysis.py index 1a0ac92837..27a588a550 100644 --- a/qiskit_experiments/library/characterization/analysis/ramsey_xy_analysis.py +++ b/qiskit_experiments/library/characterization/analysis/ramsey_xy_analysis.py @@ -16,11 +16,8 @@ import lmfit import numpy as np -from uncertainties import unumpy as unp import qiskit_experiments.curve_analysis as curve -import qiskit_experiments.visualization as vis -from qiskit_experiments.framework import ExperimentData class RamseyXYAnalysis(curve.CurveAnalysis): @@ -209,398 +206,3 @@ def _evaluate_quality(self, fit_data: curve.CurveFitResult) -> Union[str, None]: return "good" return "bad" - - -class StarkRamseyXYAmpScanAnalysis(curve.CurveAnalysis): - r"""Ramsey XY analysis for the Stark shifted phase sweep. - - # section: overview - - This analysis is a variant of :class:`RamseyXYAnalysis`. In both cases, the X and Y - data are treated as the real and imaginary parts of a complex oscillating signal. - In :class:`RamseyXYAnalysis`, the data are fit assuming a phase varying linearly with - the x-data corresponding to a constant frequency and assuming an exponentially - decaying amplitude. By contrast, in this model, the phase is assumed to be - a third order polynomial :math:`\theta(x)` of the x-data. - Additionally, the amplitude is not assumed to follow a specific form. - Techniques to compute a good initial guess for the polynomial coefficients inside - a trigonometric function like this are not trivial. Instead, this analysis extracts the - raw phase and runs fits on the extracted data to a polynomial :math:`\theta(x)` directly. - - The measured P1 values for a Ramsey X and Y experiment can be written in the form of - a trignometric function taking the phase polynomial :math:`\theta(x)`: - - .. math:: - - P_X = \text{amp}(x) \cdot \cos \theta(x) + \text{offset},\\ - P_Y = \text{amp}(x) \cdot \sin \theta(x) + \text{offset}. - - Hence the phase polynomial can be extracted as follows - - .. math:: - - \theta(x) = \tan^{-1} \frac{P_Y}{P_X}. - - Because the arctangent is implemented by the ``atan2`` function - defined in :math:`[-\pi, \pi]`, the computed :math:`\theta(x)` is unwrapped to - ensure continuous phase evolution. - - We call attention to the fact that :math:`\text{amp}(x)` is also Stark tone amplitude - dependent because of the qubit frequency dependence of the dephasing rate. - In general :math:`\text{amp}(x)` is unpredictable due to dephasing from - two-level systems distributed randomly in frequency - or potentially due to qubit heating. This prevents us from precisely fitting - the raw :math:`P_X`, :math:`P_Y` data. Fitting only the phase data makes the - analysis robust to amplitude dependent dephasing. - - In this analysis, the phase polynomial is defined as - - .. math:: - - \theta(x) = 2 \pi t_S f_S(x) - - where - - .. math:: - - f_S(x) = c_1 x + c_2 x^2 + c_3 x^3 + f_{\rm err}, - - denotes the Stark shift. For the lowest order perturbative expansion of a single driven qubit, - the Stark shift is a quadratic function of :math:`x`, but linear and cubic terms - and a constant offset are also considered to account for - other effects, e.g. strong drive, collisions, TLS, and so forth, - and frequency mis-calibration, respectively. - - # section: fit_model - - .. math:: - - \theta^\nu(x) = c_1^\nu x + c_2^\nu x^2 + c_3^\nu x^3 + f_{\rm err}, - - where :math:`\nu \in \{+, -\}`. - The Stark shift is asymmetric with respect to :math:`x=0`, because of the - anti-crossings of higher energy levels. In a typical transmon qubit, - these levels appear only in :math:`f_S < 0` because of the negative anharmonicity. - To precisely fit the results, this analysis uses different model parameters - for positive (:math:`x > 0`) and negative (:math:`x < 0`) shift domains. - - # section: fit_parameters - - defpar c_1^+: - desc: The linear term coefficient of the positive Stark shift - (fit parameter: ``stark_pos_coef_o1``). - init_guess: 0. - bounds: None - - defpar c_2^+: - desc: The quadratic term coefficient of the positive Stark shift. - This parameter must be positive because Stark amplitude is chosen to - induce blue shift when its sign is positive. - Note that the quadratic term is the primary term - (fit parameter: ``stark_pos_coef_o2``). - init_guess: 1e6. - bounds: [0, inf] - - defpar c_3^+: - desc: The cubic term coefficient of the positive Stark shift - (fit parameter: ``stark_pos_coef_o3``). - init_guess: 0. - bounds: None - - defpar c_1^-: - desc: The linear term coefficient of the negative Stark shift. - (fit parameter: ``stark_neg_coef_o1``). - init_guess: 0. - bounds: None - - defpar c_2^-: - desc: The quadratic term coefficient of the negative Stark shift. - This parameter must be negative because Stark amplitude is chosen to - induce red shift when its sign is negative. - Note that the quadratic term is the primary term - (fit parameter: ``stark_neg_coef_o2``). - init_guess: -1e6. - bounds: [-inf, 0] - - defpar c_3^-: - desc: The cubic term coefficient of the negative Stark shift - (fit parameter: ``stark_neg_coef_o3``). - init_guess: 0. - bounds: None - - defpar f_{\rm err}: - desc: Constant phase accumulation which is independent of the Stark tone amplitude. - (fit parameter: ``stark_ferr``). - init_guess: 0 - bounds: None - - # section: see_also - - :class:`qiskit_experiments.library.characterization.analysis.ramsey_xy_analysis.RamseyXYAnalysis` - - """ - - def __init__(self): - - models = [ - lmfit.models.ExpressionModel( - expr="c1_pos * x + c2_pos * x**2 + c3_pos * x**3 + f_err", - name="FREQpos", - ), - lmfit.models.ExpressionModel( - expr="c1_neg * x + c2_neg * x**2 + c3_neg * x**3 + f_err", - name="FREQneg", - ), - ] - super().__init__(models=models) - - @classmethod - def _default_options(cls): - """Default analysis options. - - Analysis Options: - pulse_len (float): Duration of effective Stark pulse in units of sec. - """ - ramsey_plotter = vis.CurvePlotter(vis.MplDrawer()) - ramsey_plotter.set_figure_options( - xlabel="Stark tone amplitude", - ylabel=["Stark shift", "P1"], - yval_unit=["Hz", None], - series_params={ - "Fpos": { - "color": "#123FE8", - "symbol": "^", - "label": "", - "canvas": 0, - }, - "Fneg": { - "color": "#123FE8", - "symbol": "v", - "label": "", - "canvas": 0, - }, - "Xpos": { - "color": "#123FE8", - "symbol": "o", - "label": "Ramsey X", - "canvas": 1, - }, - "Ypos": { - "color": "#6312E8", - "symbol": "^", - "label": "Ramsey Y", - "canvas": 1, - }, - "Xneg": { - "color": "#E83812", - "symbol": "o", - "label": "Ramsey X", - "canvas": 1, - }, - "Yneg": { - "color": "#E89012", - "symbol": "^", - "label": "Ramsey Y", - "canvas": 1, - }, - }, - sharey=False, - ) - ramsey_plotter.set_options(subplots=(2, 1), style=vis.PlotStyle({"figsize": (10, 8)})) - - options = super()._default_options() - options.update_options( - data_subfit_map={ - "Xpos": {"series": "X", "direction": "pos"}, - "Ypos": {"series": "Y", "direction": "pos"}, - "Xneg": {"series": "X", "direction": "neg"}, - "Yneg": {"series": "Y", "direction": "neg"}, - }, - result_parameters=[ - curve.ParameterRepr("c1_pos", "stark_pos_coef_o1", "Hz"), - curve.ParameterRepr("c2_pos", "stark_pos_coef_o2", "Hz"), - curve.ParameterRepr("c3_pos", "stark_pos_coef_o3", "Hz"), - curve.ParameterRepr("c1_neg", "stark_neg_coef_o1", "Hz"), - curve.ParameterRepr("c2_neg", "stark_neg_coef_o2", "Hz"), - curve.ParameterRepr("c3_neg", "stark_neg_coef_o3", "Hz"), - curve.ParameterRepr("f_err", "stark_ferr", "Hz"), - ], - plotter=ramsey_plotter, - fit_category="freq", - pulse_len=None, - ) - - return options - - def _freq_phase_coef(self) -> float: - """Return a coefficient to convert frequency into phase value.""" - try: - return 2 * np.pi * self.options.pulse_len - except TypeError as ex: - raise TypeError( - "A float-value duration in units of sec of the Stark pulse must be provided. " - f"The pulse_len option value {self.options.pulse_len} is not valid." - ) from ex - - def _format_data( - self, - curve_data: curve.ScatterTable, - category: str = "freq", - ) -> curve.ScatterTable: - - curve_data = super()._format_data(curve_data, category="ramsey_xy") - ramsey_xy = curve_data[curve_data.category == "ramsey_xy"] - - # Create phase data by arctan(Y/X) - columns = list(curve_data.columns) - phase_data = np.empty((0, len(columns))) - y_mean = ramsey_xy.yval.mean() - - grouped = ramsey_xy.groupby("name") - for m_id, direction in enumerate(("pos", "neg")): - x_quadrature = grouped.get_group(f"X{direction}") - y_quadrature = grouped.get_group(f"Y{direction}") - if not np.array_equal(x_quadrature.xval, y_quadrature.xval): - raise ValueError( - "Amplitude values of X and Y quadrature are different. " - "Same values must be used." - ) - x_uarray = unp.uarray(x_quadrature.yval, x_quadrature.yerr) - y_uarray = unp.uarray(y_quadrature.yval, y_quadrature.yerr) - - amplitudes = x_quadrature.xval.to_numpy() - - # pylint: disable=no-member - phase = unp.arctan2(y_uarray - y_mean, x_uarray - y_mean) - phase_n = unp.nominal_values(phase) - phase_s = unp.std_devs(phase) - - # Unwrap phase - # We assume a smooth slope and correct 2pi phase jump to minimize the change of the slope. - unwrapped_phase = np.unwrap(phase_n) - if amplitudes[0] < 0: - # Preserve phase value closest to 0 amplitude - unwrapped_phase = unwrapped_phase + (phase_n[-1] - unwrapped_phase[-1]) - - # Store new data - tmp = np.empty((len(amplitudes), len(columns)), dtype=object) - tmp[:, columns.index("xval")] = amplitudes - tmp[:, columns.index("yval")] = unwrapped_phase / self._freq_phase_coef() - tmp[:, columns.index("yerr")] = phase_s / self._freq_phase_coef() - tmp[:, columns.index("name")] = f"FREQ{direction}" - tmp[:, columns.index("class_id")] = m_id - tmp[:, columns.index("shots")] = x_quadrature.shots + y_quadrature.shots - tmp[:, columns.index("category")] = category - phase_data = np.r_[phase_data, tmp] - - return curve_data.append_list_values(other=phase_data) - - 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(c2_pos=(0, np.inf), c2_neg=(-np.inf, 0)) - user_opt.p0.set_if_empty( - c1_pos=0, c2_pos=1e6, c3_pos=0, c1_neg=0, c2_neg=-1e6, c3_neg=0, f_err=0 - ) - return user_opt - - def _create_figures( - self, - curve_data: curve.ScatterTable, - ) -> List["matplotlib.figure.Figure"]: - - # plot unwrapped phase on first axis - for d in ("pos", "neg"): - sub_data = curve_data[(curve_data.name == f"FREQ{d}") & (curve_data.category == "freq")] - self.plotter.set_series_data( - series_name=f"F{d}", - x_formatted=sub_data.xval.to_numpy(), - y_formatted=sub_data.yval.to_numpy(), - y_formatted_err=sub_data.yerr.to_numpy(), - ) - - # plot raw RamseyXY plot on second axis - for name in ("Xpos", "Ypos", "Xneg", "Yneg"): - sub_data = curve_data[(curve_data.name == name) & (curve_data.category == "ramsey_xy")] - self.plotter.set_series_data( - series_name=name, - x_formatted=sub_data.xval.to_numpy(), - y_formatted=sub_data.yval.to_numpy(), - y_formatted_err=sub_data.yerr.to_numpy(), - ) - - # find base and amplitude guess - ramsey_xy = curve_data[curve_data.category == "ramsey_xy"] - offset_guess = 0.5 * (ramsey_xy.yval.min() + ramsey_xy.yval.max()) - amp_guess = 0.5 * np.ptp(ramsey_xy.yval) - - # plot frequency and Ramsey fit lines - line_data = curve_data[curve_data.category == "fitted"] - for direction in ("pos", "neg"): - sub_data = line_data[line_data.name == f"FREQ{direction}"] - if len(sub_data) == 0: - continue - xval = sub_data.xval.to_numpy() - yn = sub_data.yval.to_numpy() - ys = sub_data.yerr.to_numpy() - yval = unp.uarray(yn, ys) * self._freq_phase_coef() - - # Ramsey fit lines are predicted from the phase fit line. - # Note that this line doesn't need to match with the expeirment data - # because Ramsey P1 data may fluctuate due to phase damping. - - # pylint: disable=no-member - ramsey_cos = amp_guess * unp.cos(yval) + offset_guess - ramsey_sin = amp_guess * unp.sin(yval) + offset_guess - - self.plotter.set_series_data( - series_name=f"F{direction}", - x_interp=xval, - y_interp=yn, - ) - self.plotter.set_series_data( - series_name=f"X{direction}", - x_interp=xval, - y_interp=unp.nominal_values(ramsey_cos), - ) - self.plotter.set_series_data( - series_name=f"Y{direction}", - x_interp=xval, - y_interp=unp.nominal_values(ramsey_sin), - ) - - if np.isfinite(ys).all(): - self.plotter.set_series_data( - series_name=f"F{direction}", - y_interp_err=ys, - ) - self.plotter.set_series_data( - series_name=f"X{direction}", - y_interp_err=unp.std_devs(ramsey_cos), - ) - self.plotter.set_series_data( - series_name=f"Y{direction}", - y_interp_err=unp.std_devs(ramsey_sin), - ) - return [self.plotter.figure()] - - def _initialize( - self, - experiment_data: ExperimentData, - ): - super()._initialize(experiment_data) - - # Set scaling factor to convert phase to frequency - if "stark_length" in experiment_data.metadata: - self.set_options(pulse_len=experiment_data.metadata["stark_length"]) diff --git a/qiskit_experiments/library/characterization/analysis/t1_analysis.py b/qiskit_experiments/library/characterization/analysis/t1_analysis.py index e5887704a5..9ef0ed3bc3 100644 --- a/qiskit_experiments/library/characterization/analysis/t1_analysis.py +++ b/qiskit_experiments/library/characterization/analysis/t1_analysis.py @@ -12,19 +12,12 @@ """ T1 Analysis class. """ -from typing import Union, Tuple, List, Dict +from typing import Union import numpy as np -from qiskit_ibm_experiment import IBMExperimentService -from qiskit_ibm_experiment.exceptions import IBMApiError -from uncertainties import unumpy as unp import qiskit_experiments.curve_analysis as curve -import qiskit_experiments.data_processing as dp -import qiskit_experiments.visualization as vis -from qiskit_experiments.data_processing.exceptions import DataProcessorError -from qiskit_experiments.database_service.device_component import Qubit -from qiskit_experiments.framework import BaseAnalysis, ExperimentData, AnalysisResultData, Options +from qiskit_experiments.framework import Options class T1Analysis(curve.DecayAnalysis): @@ -139,261 +132,3 @@ def _format_data( if avg_slope > 0: curve_data.yval = 1 - curve_data.yval return super()._format_data(curve_data) - - -class StarkP1SpectAnalysis(BaseAnalysis): - """Analysis class for StarkP1Spectroscopy. - - # section: overview - - The P1 landscape is hardly predictable because of the random appearance of - lossy TLS notches, and hence this analysis doesn't provide any - generic mathematical model to fit the measurement data. - A developer may subclass this to conduct own analysis. - - This analysis just visualizes the measured P1 values against Stark tone amplitudes. - The tone amplitudes can be converted into the amount of Stark shift - when the calibrated coefficients are provided in the analysis option, - or the calibration experiment results are available in the result database. - - # section: see_also - :class:`qiskit_experiments.library.characterization.ramsey_xy.StarkRamseyXYAmpScan` - - """ - - stark_coefficients_names = [ - "stark_pos_coef_o1", - "stark_pos_coef_o2", - "stark_pos_coef_o3", - "stark_neg_coef_o1", - "stark_neg_coef_o2", - "stark_neg_coef_o3", - "stark_ferr", - ] - - @property - def plotter(self) -> vis.CurvePlotter: - """Curve plotter instance.""" - return self.options.plotter - - @classmethod - def _default_options(cls) -> Options: - """Default analysis options. - - Analysis Options: - plotter (Plotter): Plotter to visualize P1 landscape. - data_processor (DataProcessor): Data processor to compute P1 value. - stark_coefficients (Union[Dict, str]): Dictionary of Stark shift coefficients to - convert tone amplitudes into amount of Stark shift. This dictionary must include - all keys defined in :attr:`.StarkP1SpectAnalysis.stark_coefficients_names`, - which are calibrated with :class:`.StarkRamseyXYAmpScan`. - Alternatively, it searches for these coefficients in the result database - when "latest" is set. This requires having the experiment service set in - the experiment data to analyze. - x_key (str): Key of the circuit metadata to represent x value. - """ - options = super()._default_options() - - p1spect_plotter = vis.CurvePlotter(vis.MplDrawer()) - p1spect_plotter.set_figure_options( - xlabel="Stark amplitude", - ylabel="P(1)", - xscale="quadratic", - ) - - options.update_options( - plotter=p1spect_plotter, - data_processor=dp.DataProcessor("counts", [dp.Probability("1")]), - stark_coefficients="latest", - x_key="xval", - ) - return options - - # pylint: disable=unused-argument - def _run_spect_analysis( - self, - xdata: np.ndarray, - ydata: np.ndarray, - ydata_err: np.ndarray, - ) -> List[AnalysisResultData]: - """Run further analysis on the spectroscopy data. - - .. note:: - A subclass can overwrite this method to conduct analysis. - - Args: - xdata: X values. This is either amplitudes or frequencies. - ydata: Y values. This is P1 values measured at different Stark tones. - ydata_err: Sampling error of the Y values. - - Returns: - A list of analysis results. - """ - return [] - - @classmethod - def retrieve_coefficients_from_service( - cls, - service: IBMExperimentService, - qubit: int, - backend: str, - ) -> Dict: - """Retrieve stark coefficient dictionary from the experiment service. - - Args: - service: A valid experiment service instance. - qubit: Qubit index. - backend: Name of the backend. - - Returns: - A dictionary of Stark coefficients to convert amplitude to frequency. - None value is returned when the dictionary is incomplete. - """ - out = {} - try: - for name in cls.stark_coefficients_names: - results = service.analysis_results( - device_components=[str(Qubit(qubit))], - result_type=name, - backend_name=backend, - sort_by=["creation_datetime:desc"], - ) - if len(results) == 0: - return None - result_data = getattr(results[0], "result_data") - out[name] = result_data["value"] - except (IBMApiError, ValueError, KeyError, AttributeError): - return None - return out - - @classmethod - def estimate_minmax_frequencies( - cls, - coefficients: Dict[str, float], - max_amplitudes: Tuple[float, float] = (-0.9, 0.9), - ) -> Tuple[float, float]: - """Inquire maximum and minimum Stark shift available within specified amplitude range. - - Args: - coefficients: A dictionary of Stark coefficients. - max_amplitudes: Minimum and maximum amplitude. - - Returns: - Tuple of minimum and maximum frequency. - - Raises: - KeyError: When coefficients are incomplete. - """ - missing = set(cls.stark_coefficients_names) - coefficients.keys() - if any(missing): - raise KeyError( - "Following coefficient data is missing in the " - f"stark 'coefficients' dictionary: {missing}." - ) - - names = cls.stark_coefficients_names # alias - pos_idxs = [2, 1, 0] - neg_idxs = [5, 4, 3] - - freqs = [] - for idxs, max_amp in zip((neg_idxs, pos_idxs), max_amplitudes): - # Solve for inflection points by computing the point where derivertive becomes zero. - solutions = np.roots( - [deriv * coefficients[names[idx]] for deriv, idx in zip([3, 2, 1], idxs)] - ) - inflection_points = solutions[ - (solutions.imag == 0) & (np.sign(solutions) == np.sign(max_amp)) - ] - if len(inflection_points) == 0: - amp = max_amp - else: - # When multiple inflection points are found, use the most outer one. - # There could be a small inflection point around amp=0, - # when the first order term is significant. - amp = min(max_amp, max(inflection_points, key=abs), key=abs) - polyfun = np.poly1d([coefficients[names[idx]] for idx in [*idxs, 6]]) - freqs.append(polyfun(amp)) - return tuple(freqs) - - def _convert_axis( - self, - xdata: np.ndarray, - coefficients: Dict[str, float], - ) -> np.ndarray: - """A helper method to convert x-axis. - - Args: - xdata: An array of Stark tone amplitude. - coefficients: Stark coefficients to convert amplitudes into frequencies. - - Returns: - An array of amount of Stark shift. - """ - names = self.stark_coefficients_names # alias - positive = np.poly1d([coefficients[names[idx]] for idx in [2, 1, 0, 6]]) - negative = np.poly1d([coefficients[names[idx]] for idx in [5, 4, 3, 6]]) - - new_xdata = np.where(xdata > 0, positive(xdata), negative(xdata)) - self.plotter.set_figure_options( - xlabel="Stark shift", - xval_unit="Hz", - xscale="linear", - ) - return new_xdata - - def _run_analysis( - self, - experiment_data: ExperimentData, - ) -> Tuple[List[AnalysisResultData], List["matplotlib.figure.Figure"]]: - - x_key = self.options.x_key - - # Get calibrated Stark tone coefficients - if self.options.stark_coefficients == "latest" and experiment_data.service is not None: - # Get value from service - stark_coeffs = self.retrieve_coefficients_from_service( - service=experiment_data.service, - qubit=experiment_data.metadata["physical_qubits"][0], - backend=experiment_data.backend_name, - ) - elif isinstance(self.options.stark_coefficients, dict): - # Get value from experiment options - missing = set(self.stark_coefficients_names) - self.options.stark_coefficients.keys() - if any(missing): - raise KeyError( - "Following coefficient data is missing in the " - f"'stark_coefficients' dictionary: {missing}." - ) - stark_coeffs = self.options.stark_coefficients - else: - # No calibration is available - stark_coeffs = None - - # Compute P1 value and sampling error - data = experiment_data.data() - try: - xdata = np.asarray([datum["metadata"][x_key] for datum in data], dtype=float) - except KeyError as ex: - raise DataProcessorError( - f"X value key {x_key} is not defined in circuit metadata." - ) from ex - ydata_ufloat = self.options.data_processor(data) - ydata = unp.nominal_values(ydata_ufloat) - ydata_err = unp.std_devs(ydata_ufloat) - - # Convert x-axis of amplitudes into Stark shift by consuming calibrated parameters. - if stark_coeffs: - xdata = self._convert_axis(xdata, stark_coeffs) - - # Draw figures and create analysis results. - self.plotter.set_series_data( - series_name="stark_p1", - x_formatted=xdata, - y_formatted=ydata, - y_formatted_err=ydata_err, - x_interp=xdata, - y_interp=ydata, - ) - analysis_results = self._run_spect_analysis(xdata, ydata, ydata_err) - - return analysis_results, [self.plotter.figure()] diff --git a/qiskit_experiments/library/characterization/ramsey_xy.py b/qiskit_experiments/library/characterization/ramsey_xy.py index 1d0238b652..ccb1987481 100644 --- a/qiskit_experiments/library/characterization/ramsey_xy.py +++ b/qiskit_experiments/library/characterization/ramsey_xy.py @@ -1,6 +1,6 @@ # This code is part of Qiskit. # -# (C) Copyright IBM 2021, 2023. +# (C) Copyright IBM 2021. # # 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 @@ -11,27 +11,16 @@ # that they have been altered from the originals. """Ramsey XY frequency characterization experiment.""" -import warnings -from typing import List, Tuple, Dict, Optional, Sequence +from typing import List, Optional, Sequence import numpy as np -from qiskit import pulse -from qiskit.circuit import QuantumCircuit, Gate, ParameterExpression, Parameter +from qiskit.circuit import QuantumCircuit, Parameter from qiskit.providers.backend import Backend from qiskit.qobj.utils import MeasLevel -from qiskit.utils import optionals as _optional from qiskit_experiments.framework import BaseExperiment, Options, BackendTiming from qiskit_experiments.framework.restless_mixin import RestlessMixin -from qiskit_experiments.library.characterization.analysis import ( - RamseyXYAnalysis, - StarkRamseyXYAmpScanAnalysis, -) - -if _optional.HAS_SYMENGINE: - import symengine as sym -else: - import sympy as sym +from qiskit_experiments.library.characterization.analysis import RamseyXYAnalysis class RamseyXY(BaseExperiment, RestlessMixin): @@ -208,612 +197,3 @@ def _metadata(self): if hasattr(self.run_options, run_opt): metadata[run_opt] = getattr(self.run_options, run_opt) return metadata - - -class StarkRamseyXY(BaseExperiment): - """A sign-sensitive experiment to measure the frequency of a qubit under a pulsed Stark tone. - - # section: overview - - This experiment is a variant of :class:`.RamseyXY` with a pulsed Stark tone - and consists of the following two circuits: - - .. parsed-literal:: - - (Ramsey X) The pulse before measurement rotates by pi-half around the X axis - - ┌────┐┌────────┐┌───┐┌───────────────┐┌────────┐┌────┐┌─┐ - q: ┤ √X ├┤ StarkV ├┤ X ├┤ StarkU(delay) ├┤ Rz(-π) ├┤ √X ├┤M├ - └────┘└────────┘└───┘└───────────────┘└────────┘└────┘└╥┘ - c: 1/═══════════════════════════════════════════════════════╩═ - 0 - - (Ramsey Y) The pulse before measurement rotates by pi-half around the Y axis - - ┌────┐┌────────┐┌───┐┌───────────────┐┌───────────┐┌────┐┌─┐ - q: ┤ √X ├┤ StarkV ├┤ X ├┤ StarkU(delay) ├┤ Rz(-3π/2) ├┤ √X ├┤M├ - └────┘└────────┘└───┘└───────────────┘└───────────┘└────┘└╥┘ - c: 1/══════════════════════════════════════════════════════════╩═ - 0 - - In principle, the sequence is a variant of :class:`.RamseyXY` circuit. - However, the delay in between √X gates is replaced with an off-resonant drive. - This off-resonant drive shifts the qubit frequency due to the - Stark effect and causes it to accumulate phase during the - Ramsey sequence. This frequency shift is a function of the - offset of the Stark tone frequency from the qubit frequency - and of the magnitude of the tone. - - Note that the Stark tone pulse (StarkU) takes the form of a flat-topped Gaussian envelope. - The magnitude of the pulse varies in time during its rising and falling edges. - It is difficult to characterize the net phase accumulation of the qubit during the - edges of the pulse when the frequency shift is varying with the pulse amplitude. - In order to simplify the analysis, an additional pulse (StarkV) - involving only the edges of StarkU is added to the sequence. - The sign of the phase accumulation is inverted for StarkV relative to that of StarkU - by inserting an X gate in between the two pulses. - - This technique allows the experiment to accumulate only the net phase - during the flat-top part of the StarkU pulse with constant magnitude. - - # section: analysis_ref - :py:class:`RamseyXYAnalysis` - - # section: see_also - :class:`qiskit_experiments.library.characterization.ramsey_xy.RamseyXY` - - # section: manual - :doc:`/manuals/characterization/stark_experiment` - - """ - - def __init__( - self, - physical_qubits: Sequence[int], - backend: Optional[Backend] = None, - **experiment_options, - ): - """Create new experiment. - - Args: - physical_qubits: Index of physical qubit. - backend: Optional, the backend to run the experiment on. - experiment_options: Experiment options. See the class documentation or - ``self._default_experiment_options`` for descriptions. - """ - self._timing = None - - super().__init__( - physical_qubits=physical_qubits, - analysis=RamseyXYAnalysis(), - backend=backend, - ) - self.set_experiment_options(**experiment_options) - - @classmethod - def _default_experiment_options(cls) -> Options: - """Default experiment options. - - Experiment Options: - stark_amp (float): A single float parameter to represent the magnitude of Stark tone - and the sign of expected Stark shift. - See :ref:`stark_tone_implementation` for details. - stark_channel (PulseChannel): Pulse channel to apply Stark tones. - If not provided, the same channel with the qubit drive is used. - See :ref:`stark_channel_consideration` for details. - stark_freq_offset (float): Offset of Stark tone frequency from the qubit frequency. - This offset should be sufficiently large so that the Stark pulse - does not Rabi drive the qubit. - See :ref:`stark_frequency_consideration` for details. - stark_sigma (float): Gaussian sigma of the rising and falling edges - of the Stark tone, in seconds. - stark_risefall (float): Ratio of sigma to the duration of - the rising and falling edges of the Stark tone. - min_freq (float): Minimum frequency that this experiment is guaranteed to resolve. - Note that fitter algorithm :class:`.RamseyXYAnalysis` of this experiment - is still capable of fitting experiment data with lower frequency. - max_freq (float): Maximum frequency that this experiment can resolve. - delays (list[float]): The list of delays if set that will be scanned in the - experiment. If not set, then evenly spaced delays with interval - computed from ``min_freq`` and ``max_freq`` are used. - See :meth:`StarkRamseyXY.delays` for details. - """ - options = super()._default_experiment_options() - options.update_options( - stark_amp=0.0, - stark_channel=None, - stark_freq_offset=80e6, - stark_sigma=15e-9, - stark_risefall=2, - min_freq=5e6, - max_freq=100e6, - delays=None, - ) - options.set_validator("stark_freq_offset", (0, np.inf)) - options.set_validator("stark_channel", pulse.channels.PulseChannel) - return options - - def _set_backend(self, backend: Backend): - super()._set_backend(backend) - self._timing = BackendTiming(backend) - - def set_experiment_options(self, **fields): - _warning_circuit_length = 300 - - # Do validation for circuit number - min_freq = fields.get("min_freq", None) - max_freq = fields.get("max_freq", None) - delays = fields.get("delays", None) - if min_freq is not None and max_freq is not None: - if delays: - warnings.warn( - "Experiment option 'min_freq' and 'max_freq' are ignored " - "when 'delays' are explicitly specified.", - UserWarning, - ) - else: - n_expr_circs = 2 * int(2 * max_freq / min_freq) # delays * (x, y) - max_circs_per_job = None - if self._backend_data: - max_circs_per_job = self._backend_data.max_circuits() - if n_expr_circs > (max_circs_per_job or _warning_circuit_length): - warnings.warn( - f"Provided configuration generates {n_expr_circs} circuits. " - "You can set smaller 'max_freq' or larger 'min_freq' to reduce this number. " - "This experiment is still executable but your execution may consume " - "unnecessary long quantum device time, and result may suffer " - "device parameter drift in consequence of the long execution time.", - UserWarning, - ) - # Do validation for spectrum overlap to avoid real excitation - stark_freq_offset = fields.get("stark_freq_offset", None) - stark_sigma = fields.get("stark_sigma", None) - if stark_freq_offset is not None and stark_sigma is not None: - if stark_freq_offset < 1 / stark_sigma: - warnings.warn( - "Provided configuration may induce coherent state exchange between qubit levels " - "because of the potential spectrum overlap. You can avoid this by " - "increasing the 'stark_sigma' or 'stark_freq_offset'. " - "Note that this experiment is still executable.", - UserWarning, - ) - pass - - super().set_experiment_options(**fields) - - def parameters(self) -> np.ndarray: - """Delay values to use in circuits. - - .. note:: - - The delays are computed with the ``min_freq`` and ``max_freq`` experiment options. - The maximum point is computed from the ``min_freq`` to guarantee the result - contains at least one Ramsey oscillation cycle at this frequency. - The interval is computed from the ``max_freq`` to sample with resolution - such that the Nyquist frequency is twice ``max_freq``. - - Returns: - The list of delays to use for the different circuits based on the - experiment options. - - Raises: - ValueError: When ``min_freq`` is larger than ``max_freq``. - """ - opt = self.experiment_options # alias - - if opt.delays is None: - if opt.min_freq > opt.max_freq: - raise ValueError("Experiment option 'min_freq' must be smaller than 'max_freq'.") - # Delay is longer enough to capture 1 cycle of the minimum frequency. - # Fitter can still accurately fit samples shorter than 1 cycle. - max_period = 1 / opt.min_freq - # Inverse of interval should be greater than Nyquist frequency. - sampling_freq = 2 * opt.max_freq - interval = 1 / sampling_freq - return np.arange(0, max_period, interval) - return opt.delays - - def parameterized_circuits(self) -> Tuple[QuantumCircuit, ...]: - """Create circuits with parameters for Ramsey XY experiment with Stark tone. - - Returns: - Quantum template circuits for Ramsey X and Ramsey Y experiment. - """ - opt = self.experiment_options # alias - param = Parameter("delay") - - # Pulse gates - stark_v = Gate("StarkV", 1, []) - stark_u = Gate("StarkU", 1, [param]) - - # Note that Stark tone yields negative (positive) frequency shift - # when the Stark tone frequency is higher (lower) than qubit f01 frequency. - # This choice gives positive frequency shift with positive Stark amplitude. - qubit_f01 = self._backend_data.drive_freqs[self.physical_qubits[0]] - stark_freq = qubit_f01 - np.sign(opt.stark_amp) * opt.stark_freq_offset - stark_amp = np.abs(opt.stark_amp) - stark_channel = opt.stark_channel or pulse.DriveChannel(self.physical_qubits[0]) - ramps_dt = self._timing.round_pulse(time=2 * opt.stark_risefall * opt.stark_sigma) - sigma_dt = ramps_dt / 2 / opt.stark_risefall - - with pulse.build() as stark_v_schedule: - pulse.set_frequency(stark_freq, stark_channel) - pulse.play( - pulse.Gaussian( - duration=ramps_dt, - amp=stark_amp, - sigma=sigma_dt, - name="StarkV", - ), - stark_channel, - ) - - with pulse.build() as stark_u_schedule: - pulse.set_frequency(stark_freq, stark_channel) - pulse.play( - pulse.GaussianSquare( - duration=ramps_dt + param, - amp=stark_amp, - sigma=sigma_dt, - risefall_sigma_ratio=opt.stark_risefall, - name="StarkU", - ), - stark_channel, - ) - - ram_x = QuantumCircuit(1, 1) - ram_x.sx(0) - ram_x.append(stark_v, [0]) - ram_x.x(0) - ram_x.append(stark_u, [0]) - ram_x.rz(-np.pi, 0) - ram_x.sx(0) - ram_x.measure(0, 0) - ram_x.metadata = {"series": "X"} - ram_x.add_calibration( - gate=stark_v, - qubits=self.physical_qubits, - schedule=stark_v_schedule, - ) - ram_x.add_calibration( - gate=stark_u, - qubits=self.physical_qubits, - schedule=stark_u_schedule, - ) - - ram_y = QuantumCircuit(1, 1) - ram_y.sx(0) - ram_y.append(stark_v, [0]) - ram_y.x(0) - ram_y.append(stark_u, [0]) - ram_y.rz(-np.pi * 3 / 2, 0) - ram_y.sx(0) - ram_y.measure(0, 0) - ram_y.metadata = {"series": "Y"} - ram_y.add_calibration( - gate=stark_v, - qubits=self.physical_qubits, - schedule=stark_v_schedule, - ) - ram_y.add_calibration( - gate=stark_u, - qubits=self.physical_qubits, - schedule=stark_u_schedule, - ) - - return ram_x, ram_y - - def circuits(self) -> List[QuantumCircuit]: - """Create circuits. - - Returns: - A list of circuits with a variable delay. - """ - timing = BackendTiming(self.backend, min_length=0) - - ramx_circ, ramy_circ = self.parameterized_circuits() - param = next(iter(ramx_circ.parameters)) - - circs = [] - for delay in self.parameters(): - valid_delay_dt = timing.round_pulse(time=delay) - net_delay_sec = timing.pulse_time(time=delay) - - ramx_circ_assigned = ramx_circ.assign_parameters({param: valid_delay_dt}, inplace=False) - ramx_circ_assigned.metadata["xval"] = net_delay_sec - - ramy_circ_assigned = ramy_circ.assign_parameters({param: valid_delay_dt}, inplace=False) - ramy_circ_assigned.metadata["xval"] = net_delay_sec - - circs.extend([ramx_circ_assigned, ramy_circ_assigned]) - - return circs - - def _metadata(self) -> Dict[str, any]: - """Return experiment metadata for ExperimentData.""" - metadata = super()._metadata() - metadata["stark_amp"] = self.experiment_options.stark_amp - metadata["stark_freq_offset"] = self.experiment_options.stark_freq_offset - - return metadata - - -class StarkRamseyXYAmpScan(BaseExperiment): - r"""A fast characterization of Stark frequency shift by varying Stark tone amplitude. - - # section: overview - - This experiment scans Stark tone amplitude at a fixed tone duration. - The experimental circuits are identical to the :class:`.StarkRamseyXY` experiment - except that the Stark pulse amplitude is the scanned parameter rather than the pulse width. - - .. parsed-literal:: - - (Ramsey X) The pulse before measurement rotates by pi-half around the X axis - - ┌────┐┌───────────────────┐┌───┐┌───────────────────┐┌────────┐┌────┐┌─┐ - q: ┤ √X ├┤ StarkV(stark_amp) ├┤ X ├┤ StarkU(stark_amp) ├┤ Rz(-π) ├┤ √X ├┤M├ - └────┘└───────────────────┘└───┘└───────────────────┘└────────┘└────┘└╥┘ - c: 1/══════════════════════════════════════════════════════════════════════╩═ - 0 - - (Ramsey Y) The pulse before measurement rotates by pi-half around the Y axis - - ┌────┐┌───────────────────┐┌───┐┌───────────────────┐┌───────────┐┌────┐┌─┐ - q: ┤ √X ├┤ StarkV(stark_amp) ├┤ X ├┤ StarkU(stark_amp) ├┤ Rz(-3π/2) ├┤ √X ├┤M├ - └────┘└───────────────────┘└───┘└───────────────────┘└───────────┘└────┘└╥┘ - c: 1/═════════════════════════════════════════════════════════════════════════╩═ - 0 - - The AC Stark effect can be used to shift the frequency of a qubit with a microwave drive. - To calibrate a specific frequency shift, the :class:`.StarkRamseyXY` experiment can be run - to scan the Stark pulse duration at every amplitude, but such a two dimensional scan of - the tone duration and amplitude may require many circuit executions. - To avoid this overhead, the :class:`.StarkRamseyXYAmpScan` experiment fixes the - tone duration and scans only amplitude. - - Recall that an observed Ramsey oscillation in each quadrature may be represented by - - .. math:: - - {\cal E}_X(\Omega, t_S) = A e^{-t_S/\tau} \cos \left( 2\pi f_S(\Omega) t_S \right), \\ - {\cal E}_Y(\Omega, t_S) = A e^{-t_S/\tau} \sin \left( 2\pi f_S(\Omega) t_S \right), - - where :math:`f_S(\Omega)` denotes the amount of Stark shift - at a constant tone amplitude :math:`\Omega`, and :math:`t_S` is the duration of the - applied tone. For a fixed tone duration, - one can still observe the Ramsey oscillation by scanning the tone amplitude. - However, since :math:`f_S` is usually a higher order polynomial of :math:`\Omega`, - one must manage to fit the y-data for trigonometric functions with - phase which non-linearly changes with the x-data. - The :class:`.StarkRamseyXYAmpScan` experiment thus drastically reduces the number of - circuits to run in return for greater complexity in the fitting model. - - # section: analysis_ref - :py:class:`StarkRamseyXYAmpScanAnalysis` - - # section: see_also - :class:`qiskit_experiments.library.characterization.ramsey_xy.StarkRamseyXY` - :class:`qiskit_experiments.library.characterization.ramsey_xy.RamseyXY` - - # section: manual - :doc:`/manuals/characterization/stark_experiment` - - """ - - def __init__( - self, - physical_qubits: Sequence[int], - backend: Optional[Backend] = None, - **experiment_options, - ): - """Create new experiment. - - Args: - physical_qubits: Sequence with the index of the physical qubit. - backend: Optional, the backend to run the experiment on. - experiment_options: Experiment options. See the class documentation or - ``self._default_experiment_options`` for descriptions. - """ - self._timing = None - - super().__init__( - physical_qubits=physical_qubits, - analysis=StarkRamseyXYAmpScanAnalysis(), - backend=backend, - ) - self.set_experiment_options(**experiment_options) - - @classmethod - def _default_experiment_options(cls) -> Options: - """Default experiment options. - - Experiment Options: - stark_channel (PulseChannel): Pulse channel on which to apply Stark tones. - If not provided, the same channel with the qubit drive is used. - See :ref:`stark_channel_consideration` for details. - stark_freq_offset (float): Offset of Stark tone frequency from the qubit frequency. - This offset should be sufficiently large so that the Stark pulse - does not Rabi drive the qubit. - See :ref:`stark_frequency_consideration` for details. - stark_sigma (float): Gaussian sigma of the rising and falling edges - of the Stark tone, in seconds. - stark_risefall (float): Ratio of sigma to the duration of - the rising and falling edges of the Stark tone. - stark_length (float): Time to accumulate Stark shifted phase in seconds. - min_stark_amp (float): Minimum Stark tone amplitude. - max_stark_amp (float): Maximum Stark tone amplitude. - num_stark_amps (int): Number of Stark tone amplitudes to scan. - stark_amps (list[float]): The list of amplitude that will be scanned in the experiment. - If not set, then ``num_stark_amps`` evenly spaced amplitudes - between ``min_stark_amp`` and ``max_stark_amp`` are used. If ``stark_amps`` - is set, these parameters are ignored. - """ - options = super()._default_experiment_options() - options.update_options( - stark_channel=None, - stark_freq_offset=80e6, - stark_sigma=15e-9, - stark_risefall=2, - stark_length=50e-9, - min_stark_amp=-1.0, - max_stark_amp=1.0, - num_stark_amps=101, - stark_amps=None, - ) - options.set_validator("stark_freq_offset", (0, np.inf)) - options.set_validator("stark_channel", pulse.channels.PulseChannel) - return options - - def _set_backend(self, backend: Backend): - super()._set_backend(backend) - self._timing = BackendTiming(backend) - - def parameters(self) -> np.ndarray: - """Stark tone amplitudes to use in circuits. - - Returns: - The list of amplitudes to use for the different circuits based on the - experiment options. - """ - opt = self.experiment_options # alias - - if opt.stark_amps is None: - params = np.linspace(opt.min_stark_amp, opt.max_stark_amp, opt.num_stark_amps) - else: - params = np.asarray(opt.stark_amps, dtype=float) - - return params - - def parameterized_circuits(self) -> Tuple[QuantumCircuit, ...]: - """Create circuits with parameters for Ramsey XY experiment with Stark tone. - - Returns: - Quantum template circuits for Ramsey X and Ramsey Y experiment. - """ - opt = self.experiment_options # alias - param = Parameter("stark_amp") - sym_param = param._symbol_expr - - # Pulse gates - stark_v = Gate("StarkV", 1, [param]) - stark_u = Gate("StarkU", 1, [param]) - - # Note that Stark tone yields negative (positive) frequency shift - # when the Stark tone frequency is higher (lower) than qubit f01 frequency. - # This choice gives positive frequency shift with positive Stark amplitude. - qubit_f01 = self._backend_data.drive_freqs[self.physical_qubits[0]] - neg_sign_of_amp = ParameterExpression( - symbol_map={param: sym_param}, - expr=-sym.sign(sym_param), - ) - abs_of_amp = ParameterExpression( - symbol_map={param: sym_param}, - expr=sym.Abs(sym_param), - ) - stark_freq = qubit_f01 + neg_sign_of_amp * opt.stark_freq_offset - stark_channel = opt.stark_channel or pulse.DriveChannel(self.physical_qubits[0]) - ramps_dt = self._timing.round_pulse(time=2 * opt.stark_risefall * opt.stark_sigma) - sigma_dt = ramps_dt / 2 / opt.stark_risefall - width_dt = self._timing.round_pulse(time=opt.stark_length) - - with pulse.build() as stark_v_schedule: - pulse.set_frequency(stark_freq, stark_channel) - pulse.play( - pulse.Gaussian( - duration=ramps_dt, - amp=abs_of_amp, - sigma=sigma_dt, - ), - stark_channel, - ) - - with pulse.build() as stark_u_schedule: - pulse.set_frequency(stark_freq, stark_channel) - pulse.play( - pulse.GaussianSquare( - duration=ramps_dt + width_dt, - amp=abs_of_amp, - sigma=sigma_dt, - risefall_sigma_ratio=opt.stark_risefall, - ), - stark_channel, - ) - - ram_x = QuantumCircuit(1, 1) - ram_x.sx(0) - ram_x.append(stark_v, [0]) - ram_x.x(0) - ram_x.append(stark_u, [0]) - ram_x.rz(-np.pi, 0) - ram_x.sx(0) - ram_x.measure(0, 0) - ram_x.metadata = {"series": "X"} - ram_x.add_calibration( - gate=stark_v, - qubits=self.physical_qubits, - schedule=stark_v_schedule, - ) - ram_x.add_calibration( - gate=stark_u, - qubits=self.physical_qubits, - schedule=stark_u_schedule, - ) - - ram_y = QuantumCircuit(1, 1) - ram_y.sx(0) - ram_y.append(stark_v, [0]) - ram_y.x(0) - ram_y.append(stark_u, [0]) - ram_y.rz(-np.pi * 3 / 2, 0) - ram_y.sx(0) - ram_y.measure(0, 0) - ram_y.metadata = {"series": "Y"} - ram_y.add_calibration( - gate=stark_v, - qubits=self.physical_qubits, - schedule=stark_v_schedule, - ) - ram_y.add_calibration( - gate=stark_u, - qubits=self.physical_qubits, - schedule=stark_u_schedule, - ) - - return ram_x, ram_y - - def circuits(self) -> List[QuantumCircuit]: - """Create circuits. - - Returns: - A list of circuits with a variable Stark tone amplitudes. - """ - ramx_circ, ramy_circ = self.parameterized_circuits() - param = next(iter(ramx_circ.parameters)) - - circs = [] - for amp in self.parameters(): - # Add metadata "direction" to ease the filtering of the data - # by curve analysis. Indeed, the fit parameters are amplitude sign dependent. - - ramx_circ_assigned = ramx_circ.assign_parameters({param: amp}, inplace=False) - ramx_circ_assigned.metadata["xval"] = amp - ramx_circ_assigned.metadata["direction"] = "pos" if amp > 0 else "neg" - - ramy_circ_assigned = ramy_circ.assign_parameters({param: amp}, inplace=False) - ramy_circ_assigned.metadata["xval"] = amp - ramy_circ_assigned.metadata["direction"] = "pos" if amp > 0 else "neg" - - circs.extend([ramx_circ_assigned, ramy_circ_assigned]) - - return circs - - def _metadata(self) -> Dict[str, any]: - """Return experiment metadata for ExperimentData.""" - metadata = super()._metadata() - metadata["stark_length"] = self._timing.pulse_time( - time=self.experiment_options.stark_length - ) - metadata["stark_freq_offset"] = self.experiment_options.stark_freq_offset - - return metadata diff --git a/qiskit_experiments/library/characterization/t1.py b/qiskit_experiments/library/characterization/t1.py index 3873d03f60..0554599a25 100644 --- a/qiskit_experiments/library/characterization/t1.py +++ b/qiskit_experiments/library/characterization/t1.py @@ -13,24 +13,14 @@ T1 Experiment class. """ -from typing import List, Tuple, Dict, Optional, Union, Sequence +from typing import List, Optional, Union, Sequence import numpy as np -from qiskit import pulse -from qiskit.circuit import QuantumCircuit, Gate, Parameter, ParameterExpression +from qiskit.circuit import QuantumCircuit from qiskit.providers.backend import Backend -from qiskit.utils import optionals as _optional -from qiskit_experiments.framework import BackendTiming, BaseExperiment, ExperimentData, Options -from qiskit_experiments.library.characterization.analysis.t1_analysis import ( - T1Analysis, - StarkP1SpectAnalysis, -) - -if _optional.HAS_SYMENGINE: - import symengine as sym -else: - import sympy as sym +from qiskit_experiments.framework import BackendTiming, BaseExperiment, Options +from qiskit_experiments.library.characterization.analysis.t1_analysis import T1Analysis class T1(BaseExperiment): @@ -119,334 +109,3 @@ def _metadata(self): if hasattr(self.run_options, run_opt): metadata[run_opt] = getattr(self.run_options, run_opt) return metadata - - -class StarkP1Spectroscopy(BaseExperiment): - """P1 spectroscopy experiment with Stark tone. - - # section: overview - - This experiment measures a probability of the excitation state of the qubit - with a certain delay after excitation. - A Stark tone is applied during this delay to move the - qubit frequency to conduct a spectroscopy of qubit relaxation quantity. - - .. parsed-literal:: - - ┌───┐┌──────────────────┐┌─┐ - q: ┤ X ├┤ Stark(stark_amp) ├┤M├ - └───┘└──────────────────┘└╥┘ - c: 1/══════════════════════════╩═ - 0 - - Since the qubit relaxation rate may depend on the qubit frequency due to the - coupling to nearby energy levels, this experiment is useful to find out - lossy operation frequency that might be harmful to the gate fidelity [1]. - - # section: analysis_ref - :py:class:`.StarkP1SpectAnalysis` - - # section: reference - .. ref_arxiv:: 1 2105.15201 - - # section: see_also - :class:`qiskit_experiments.library.characterization.ramsey_xy.StarkRamseyXY` - :class:`qiskit_experiments.library.characterization.ramsey_xy.StarkRamseyXYAmpScan` - - # section: manual - :doc:`/manuals/characterization/stark_experiment` - - """ - - def __init__( - self, - physical_qubits: Sequence[int], - backend: Optional[Backend] = None, - **experiment_options, - ): - """ - Initialize the T1 experiment class. - - Args: - physical_qubits: Sequence with the index of the physical qubit. - backend: Optional, the backend to run the experiment on. - experiment_options: Experiment options. See the class documentation or - ``self._default_experiment_options`` for descriptions. - """ - self._timing = None - - super().__init__( - physical_qubits=physical_qubits, - analysis=StarkP1SpectAnalysis(), - backend=backend, - ) - self.set_experiment_options(**experiment_options) - - @classmethod - def _default_experiment_options(cls) -> Options: - """Default experiment options. - - Experiment Options: - t1_delay (float): The T1 delay time after excitation pulse. The delay must be - sufficiently greater than the edge duration determined by the stark_sigma. - stark_channel (PulseChannel): Pulse channel to apply Stark tones. - If not provided, the same channel with the qubit drive is used. - stark_freq_offset (float): Offset of Stark tone frequency from the qubit frequency. - This must be greater than zero not to apply Rabi drive. - stark_sigma (float): Gaussian sigma of the rising and falling edges - of the Stark tone, in seconds. - stark_risefall (float): Ratio of sigma to the duration of - the rising and falling edges of the Stark tone. - min_xval (float): Minimum x value. - max_xval (float): Maximum x value. - num_xvals (int): Number of x-values to scan. - xval_type (str): Type of x-value. Either ``amplitude`` or ``frequency``. - Setting to frequency requires pre-calibration of Stark shift coefficients. - spacing (str): A policy for the spacing to create an amplitude list from - ``min_stark_amp`` to ``max_stark_amp``. Either ``linear`` or ``quadratic`` - must be specified. - xvals (list[float]): The list of x-values that will be scanned in the experiment. - If not set, then ``num_xvals`` parameters spaced according to - the ``spacing`` policy between ``min_xval`` and ``max_xval`` are used. - If ``xvals`` is set, these parameters are ignored. - service (IBMExperimentService): A valid experiment service instance that can - provide the Stark coefficients for the qubit to run experiment. - This is required only when ``stark_coefficients`` is ``latest`` and - ``xval_type`` is ``frequency``. This value is automatically set when - a backend is attached to this experiment instance. - stark_coefficients (Union[Dict, str]): Dictionary of Stark shift coefficients to - convert tone amplitudes into amount of Stark shift. This dictionary must include - all keys defined in :attr:`.StarkP1SpectAnalysis.stark_coefficients_names`, - which are calibrated with :class:`.StarkRamseyXYAmpScan`. - Alternatively, a search for these coefficients in the result database is run - when "latest" is set. This requires having the experiment service available - in the ``backend`` set for the experiment. - """ - options = super()._default_experiment_options() - options.update_options( - t1_delay=20e-6, - stark_channel=None, - stark_freq_offset=80e6, - stark_sigma=15e-9, - stark_risefall=2, - min_xval=-1.0, - max_xval=1.0, - num_xvals=201, - xval_type="amplitude", - spacing="quadratic", - xvals=None, - service=None, - stark_coefficients="latest", - ) - options.set_validator("spacing", ["linear", "quadratic"]) - options.set_validator("xval_type", ["amplitude", "frequency"]) - options.set_validator("stark_freq_offset", (0, np.inf)) - options.set_validator("stark_channel", pulse.channels.PulseChannel) - return options - - def _set_backend(self, backend: Backend): - super()._set_backend(backend) - self._timing = BackendTiming(backend) - if self.experiment_options.service is None: - self.set_experiment_options( - service=ExperimentData.get_service_from_backend(backend), - ) - - def parameters(self) -> np.ndarray: - """Stark tone amplitudes to use in circuits. - - Returns: - The list of amplitudes to use for the different circuits based on the - experiment options. - """ - opt = self.experiment_options # alias - - if opt.xvals is None: - if opt.spacing == "linear": - params = np.linspace(opt.min_xval, opt.max_xval, opt.num_xvals) - elif opt.spacing == "quadratic": - min_sqrt = np.sign(opt.min_xval) * np.sqrt(np.abs(opt.min_xval)) - max_sqrt = np.sign(opt.max_xval) * np.sqrt(np.abs(opt.max_xval)) - lin_params = np.linspace(min_sqrt, max_sqrt, opt.num_xvals) - params = np.sign(lin_params) * lin_params**2 - else: - raise ValueError(f"Spacing option {opt.spacing} is not valid.") - else: - params = np.asarray(opt.xvals, dtype=float) - - if opt.xval_type == "frequency": - return self._frequencies_to_amplitudes(params) - return params - - def _frequencies_to_amplitudes(self, params: np.ndarray) -> np.ndarray: - """A helper method to convert frequency values to amplitude. - - Args: - params: Parameters representing a frequency of Stark shift. - - Returns: - Corresponding Stark tone amplitudes. - - Raises: - RuntimeError: When service or analysis results for Stark coefficients are not available. - TypeError: When attached analysis class is not valid. - KeyError: When stark_coefficients dictionary is provided but keys are missing. - ValueError: When specified Stark shift is not available. - """ - opt = self.experiment_options # alias - - if not isinstance(self.analysis, StarkP1SpectAnalysis): - raise TypeError( - f"Analysis class {self.analysis.__class__.__name__} is not a subclass of " - "StarkP1SpectAnalysis. Use proper analysis class to scan frequencies." - ) - coef_names = self.analysis.stark_coefficients_names - - if opt.stark_coefficients == "latest": - if opt.service is None: - raise RuntimeError( - "Experiment service is not available. Provide a dictionary of " - "Stark coefficients in the experiment options." - ) - coefficients = self.analysis.retrieve_coefficients_from_service( - service=opt.service, - qubit=self.physical_qubits[0], - backend=self._backend_data.name, - ) - if coefficients is None: - raise RuntimeError( - "Experiment results for the coefficients of the Stark shift is not found " - f"for the backend {self._backend_data.name} qubit {self.physical_qubits}." - ) - else: - missing = set(coef_names) - opt.stark_coefficients.keys() - if any(missing): - raise KeyError( - f"Following coefficient data is missing in the 'stark_coefficients': {missing}." - ) - coefficients = opt.stark_coefficients - positive = np.asarray([coefficients[coef_names[idx]] for idx in [2, 1, 0]]) - negative = np.asarray([coefficients[coef_names[idx]] for idx in [5, 4, 3]]) - offset = coefficients[coef_names[6]] - - amplitudes = np.zeros_like(params) - for idx, tgt_freq in enumerate(params): - stark_shift = tgt_freq - offset - if np.isclose(stark_shift, 0): - amplitudes[idx] = 0 - continue - if np.sign(stark_shift) > 0: - fit_coeffs = [*positive, -stark_shift] - else: - fit_coeffs = [*negative, -stark_shift] - amp_candidates = np.roots(fit_coeffs) - # Because the fit function is third order, we get three solutions here. - # Only one valid solution must exist because we assume - # a monotonic trend for Stark shift against tone amplitude in domain of definition. - criteria = np.all( - [ - # Frequency shift and tone have the same sign by definition - np.sign(amp_candidates.real) == np.sign(stark_shift), - # Tone amplitude is a real value - np.isclose(amp_candidates.imag, 0.0), - # The absolute value of tone amplitude must be less than 1.0 - np.abs(amp_candidates.real) < 1.0, - ], - axis=0, - ) - valid_amps = amp_candidates[criteria] - if len(valid_amps) == 0: - raise ValueError( - f"Stark shift at frequency value of {tgt_freq} Hz is not available on " - f"the backend {self._backend_data.name} qubit {self.physical_qubits}." - ) - if len(valid_amps) > 1: - # We assume a monotonic trend but sometimes a large third-order term causes - # inflection point and inverts the trend in larger amplitudes. - # In this case we would have more than one solutions, but we can - # take the smallerst amplitude before reaching to the inflection point. - before_inflection = np.argmin(np.abs(valid_amps.real)) - valid_amp = float(valid_amps[before_inflection].real) - else: - valid_amp = float(valid_amps.real) - amplitudes[idx] = valid_amp - - return amplitudes - - def parameterized_circuits(self) -> Tuple[QuantumCircuit, ...]: - """Create circuits with parameters for P1 experiment with Stark shift. - - Returns: - Quantum template circuit for P1 experiment. - """ - opt = self.experiment_options # alias - param = Parameter("stark_amp") - sym_param = param._symbol_expr - - # Pulse gates - stark = Gate("Stark", 1, [param]) - - # Note that Stark tone yields negative (positive) frequency shift - # when the Stark tone frequency is higher (lower) than qubit f01 frequency. - # This choice gives positive frequency shift with positive Stark amplitude. - qubit_f01 = self._backend_data.drive_freqs[self.physical_qubits[0]] - neg_sign_of_amp = ParameterExpression( - symbol_map={param: sym_param}, - expr=-sym.sign(sym_param), - ) - abs_of_amp = ParameterExpression( - symbol_map={param: sym_param}, - expr=sym.Abs(sym_param), - ) - stark_freq = qubit_f01 + neg_sign_of_amp * opt.stark_freq_offset - stark_channel = opt.stark_channel or pulse.DriveChannel(self.physical_qubits[0]) - sigma_dt = opt.stark_sigma / self._backend_data.dt - delay_dt = self._timing.round_pulse(time=opt.t1_delay) - - with pulse.build() as stark_schedule: - pulse.set_frequency(stark_freq, stark_channel) - pulse.play( - pulse.GaussianSquare( - duration=delay_dt, - amp=abs_of_amp, - sigma=sigma_dt, - risefall_sigma_ratio=opt.stark_risefall, - ), - stark_channel, - ) - - temp_t1 = QuantumCircuit(1, 1) - temp_t1.x(0) - temp_t1.append(stark, [0]) - temp_t1.measure(0, 0) - temp_t1.add_calibration( - gate=stark, - qubits=self.physical_qubits, - schedule=stark_schedule, - ) - - return (temp_t1,) - - def circuits(self) -> List[QuantumCircuit]: - """Create circuits. - - Returns: - A list of P1 circuits with a variable Stark tone amplitudes. - """ - (t1_circ,) = self.parameterized_circuits() - param = next(iter(t1_circ.parameters)) - - circs = [] - for amp in self.parameters(): - t1_assigned = t1_circ.assign_parameters({param: amp}, inplace=False) - t1_assigned.metadata = {"xval": amp} - circs.append(t1_assigned) - - return circs - - def _metadata(self) -> Dict[str, any]: - """Return experiment metadata for ExperimentData.""" - metadata = super()._metadata() - metadata["stark_freq_offset"] = self.experiment_options.stark_freq_offset - - return metadata diff --git a/qiskit_experiments/library/driven_freq_tuning/__init__.py b/qiskit_experiments/library/driven_freq_tuning/__init__.py new file mode 100644 index 0000000000..9cdd2faf99 --- /dev/null +++ b/qiskit_experiments/library/driven_freq_tuning/__init__.py @@ -0,0 +1,66 @@ +# 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. + +""" +=============================================================================================== +Driven Frequency Tuning (:mod:`qiskit_experiments.library.driven_freq_tuning`) +=============================================================================================== + +.. currentmodule:: qiskit_experiments.library.driven_freq_tuning + +Experiments +=========== +.. autosummary:: + :toctree: ../stubs/ + :template: autosummary/experiment.rst + + StarkRamseyXY + StarkRamseyXYAmpScan + StarkP1Spectroscopy + + +Analysis +======== + +.. autosummary:: + :toctree: ../stubs/ + :template: autosummary/analysis.rst + + StarkRamseyXYAmpScanAnalysis + StarkP1SpectAnalysis + + +Utilities +========= + +.. autosummary:: + :toctree: ../stubs/ + + StarkCoefficients + convert_amp_to_freq + convert_freq_to_amp + retrieve_coefficients_from_backend + retrieve_coefficients_from_service +""" + +from .analyses import StarkRamseyXYAmpScanAnalysis, StarkP1SpectAnalysis +from .ramsey import StarkRamseyXY +from .ramsey_amp_scan import StarkRamseyXYAmpScan +from .p1_spect import StarkP1Spectroscopy + +from .coefficient_utils import ( + StarkCoefficients, + convert_amp_to_freq, + convert_freq_to_amp, + retrieve_coefficients_from_backend, + retrieve_coefficients_from_service, +) diff --git a/qiskit_experiments/library/driven_freq_tuning/analyses.py b/qiskit_experiments/library/driven_freq_tuning/analyses.py new file mode 100644 index 0000000000..67f510a6f6 --- /dev/null +++ b/qiskit_experiments/library/driven_freq_tuning/analyses.py @@ -0,0 +1,584 @@ +# 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. +"""Stark shift analyses.""" + +from __future__ import annotations + +from typing import List, Union + +import lmfit +import numpy as np +from uncertainties import unumpy as unp + +import qiskit_experiments.curve_analysis as curve +import qiskit_experiments.data_processing as dp +import qiskit_experiments.visualization as vis +from qiskit_experiments.data_processing.exceptions import DataProcessorError +from qiskit_experiments.framework import BaseAnalysis, ExperimentData, AnalysisResultData, Options +from .coefficient_utils import ( + StarkCoefficients, + convert_amp_to_freq, + retrieve_coefficients_from_service, +) + + +class StarkRamseyXYAmpScanAnalysis(curve.CurveAnalysis): + r"""Ramsey XY analysis for the Stark shifted phase sweep. + + # section: overview + + This analysis is a variant of :class:`RamseyXYAnalysis`. In both cases, the X and Y + data are treated as the real and imaginary parts of a complex oscillating signal. + In :class:`RamseyXYAnalysis`, the data are fit assuming a phase varying linearly with + the x-data corresponding to a constant frequency and assuming an exponentially + decaying amplitude. By contrast, in this model, the phase is assumed to be + a third order polynomial :math:`\theta(x)` of the x-data. + Additionally, the amplitude is not assumed to follow a specific form. + Techniques to compute a good initial guess for the polynomial coefficients inside + a trigonometric function like this are not trivial. Instead, this analysis extracts the + raw phase and runs fits on the extracted data to a polynomial :math:`\theta(x)` directly. + + The measured P1 values for a Ramsey X and Y experiment can be written in the form of + a trignometric function taking the phase polynomial :math:`\theta(x)`: + + .. math:: + + P_X = \text{amp}(x) \cdot \cos \theta(x) + \text{offset},\\ + P_Y = \text{amp}(x) \cdot \sin \theta(x) + \text{offset}. + + Hence the phase polynomial can be extracted as follows + + .. math:: + + \theta(x) = \tan^{-1} \frac{P_Y}{P_X}. + + Because the arctangent is implemented by the ``atan2`` function + defined in :math:`[-\pi, \pi]`, the computed :math:`\theta(x)` is unwrapped to + ensure continuous phase evolution. + + We call attention to the fact that :math:`\text{amp}(x)` is also Stark tone amplitude + dependent because of the qubit frequency dependence of the dephasing rate. + In general :math:`\text{amp}(x)` is unpredictable due to dephasing from + two-level systems distributed randomly in frequency + or potentially due to qubit heating. This prevents us from precisely fitting + the raw :math:`P_X`, :math:`P_Y` data. Fitting only the phase data makes the + analysis robust to amplitude dependent dephasing. + + In this analysis, the phase polynomial is defined as + + .. math:: + + \theta(x) = 2 \pi t_S f_S(x) + + where + + .. math:: + + f_S(x) = c_1 x + c_2 x^2 + c_3 x^3 + f_{\rm err}, + + denotes the Stark shift. For the lowest order perturbative expansion of a single driven qubit, + the Stark shift is a quadratic function of :math:`x`, but linear and cubic terms + and a constant offset are also considered to account for + other effects, e.g. strong drive, collisions, TLS, and so forth, + and frequency mis-calibration, respectively. + + # section: fit_model + + .. math:: + + \theta^\nu(x) = c_1^\nu x + c_2^\nu x^2 + c_3^\nu x^3 + f_{\rm err}, + + where :math:`\nu \in \{+, -\}`. + The Stark shift is asymmetric with respect to :math:`x=0`, because of the + anti-crossings of higher energy levels. In a typical transmon qubit, + these levels appear only in :math:`f_S < 0` because of the negative anharmonicity. + To precisely fit the results, this analysis uses different model parameters + for positive (:math:`x > 0`) and negative (:math:`x < 0`) shift domains. + + # section: fit_parameters + + defpar c_1^+: + desc: The linear term coefficient of the positive Stark shift + (fit parameter: ``stark_pos_coef_o1``). + init_guess: 0. + bounds: None + + defpar c_2^+: + desc: The quadratic term coefficient of the positive Stark shift. + This parameter must be positive because Stark amplitude is chosen to + induce blue shift when its sign is positive. + Note that the quadratic term is the primary term + (fit parameter: ``stark_pos_coef_o2``). + init_guess: 1e6. + bounds: [0, inf] + + defpar c_3^+: + desc: The cubic term coefficient of the positive Stark shift + (fit parameter: ``stark_pos_coef_o3``). + init_guess: 0. + bounds: None + + defpar c_1^-: + desc: The linear term coefficient of the negative Stark shift. + (fit parameter: ``stark_neg_coef_o1``). + init_guess: 0. + bounds: None + + defpar c_2^-: + desc: The quadratic term coefficient of the negative Stark shift. + This parameter must be negative because Stark amplitude is chosen to + induce red shift when its sign is negative. + Note that the quadratic term is the primary term + (fit parameter: ``stark_neg_coef_o2``). + init_guess: -1e6. + bounds: [-inf, 0] + + defpar c_3^-: + desc: The cubic term coefficient of the negative Stark shift + (fit parameter: ``stark_neg_coef_o3``). + init_guess: 0. + bounds: None + + defpar f_{\rm err}: + desc: Constant phase accumulation which is independent of the Stark tone amplitude. + (fit parameter: ``stark_ferr``). + init_guess: 0 + bounds: None + + # section: see_also + + :class:`qiskit_experiments.library.characterization.analysis.ramsey_xy_analysis.RamseyXYAnalysis` + + """ + + def __init__(self): + + models = [ + lmfit.models.ExpressionModel( + expr="c1_pos * x + c2_pos * x**2 + c3_pos * x**3 + f_err", + name="FREQpos", + ), + lmfit.models.ExpressionModel( + expr="c1_neg * x + c2_neg * x**2 + c3_neg * x**3 + f_err", + name="FREQneg", + ), + ] + super().__init__(models=models) + + @classmethod + def _default_options(cls): + """Default analysis options. + + Analysis Options: + pulse_len (float): Duration of effective Stark pulse in units of sec. + """ + ramsey_plotter = vis.CurvePlotter(vis.MplDrawer()) + ramsey_plotter.set_figure_options( + xlabel="Stark tone amplitude", + ylabel=["Stark shift", "P1"], + yval_unit=["Hz", None], + series_params={ + "Fpos": { + "color": "#123FE8", + "symbol": "^", + "label": "", + "canvas": 0, + }, + "Fneg": { + "color": "#123FE8", + "symbol": "v", + "label": "", + "canvas": 0, + }, + "Xpos": { + "color": "#123FE8", + "symbol": "o", + "label": "Ramsey X", + "canvas": 1, + }, + "Ypos": { + "color": "#6312E8", + "symbol": "^", + "label": "Ramsey Y", + "canvas": 1, + }, + "Xneg": { + "color": "#E83812", + "symbol": "o", + "label": "Ramsey X", + "canvas": 1, + }, + "Yneg": { + "color": "#E89012", + "symbol": "^", + "label": "Ramsey Y", + "canvas": 1, + }, + }, + sharey=False, + ) + ramsey_plotter.set_options(subplots=(2, 1), style=vis.PlotStyle({"figsize": (10, 8)})) + + options = super()._default_options() + options.update_options( + data_subfit_map={ + "Xpos": {"series": "X", "direction": "pos"}, + "Ypos": {"series": "Y", "direction": "pos"}, + "Xneg": {"series": "X", "direction": "neg"}, + "Yneg": {"series": "Y", "direction": "neg"}, + }, + plotter=ramsey_plotter, + fit_category="freq", + pulse_len=None, + ) + + return options + + def _freq_phase_coef(self) -> float: + """Return a coefficient to convert frequency into phase value.""" + try: + return 2 * np.pi * self.options.pulse_len + except TypeError as ex: + raise TypeError( + "A float-value duration in units of sec of the Stark pulse must be provided. " + f"The pulse_len option value {self.options.pulse_len} is not valid." + ) from ex + + def _format_data( + self, + curve_data: curve.ScatterTable, + category: str = "freq", + ) -> curve.ScatterTable: + + curve_data = super()._format_data(curve_data, category="ramsey_xy") + ramsey_xy = curve_data[curve_data.category == "ramsey_xy"] + + # Create phase data by arctan(Y/X) + columns = list(curve_data.columns) + phase_data = np.empty((0, len(columns))) + y_mean = ramsey_xy.yval.mean() + + grouped = ramsey_xy.groupby("name") + for m_id, direction in enumerate(("pos", "neg")): + x_quadrature = grouped.get_group(f"X{direction}") + y_quadrature = grouped.get_group(f"Y{direction}") + if not np.array_equal(x_quadrature.xval, y_quadrature.xval): + raise ValueError( + "Amplitude values of X and Y quadrature are different. " + "Same values must be used." + ) + x_uarray = unp.uarray(x_quadrature.yval, x_quadrature.yerr) + y_uarray = unp.uarray(y_quadrature.yval, y_quadrature.yerr) + + amplitudes = x_quadrature.xval.to_numpy() + + # pylint: disable=no-member + phase = unp.arctan2(y_uarray - y_mean, x_uarray - y_mean) + phase_n = unp.nominal_values(phase) + phase_s = unp.std_devs(phase) + + # Unwrap phase + # We assume a smooth slope and correct 2pi phase jump to minimize the change of the slope. + unwrapped_phase = np.unwrap(phase_n) + if amplitudes[0] < 0: + # Preserve phase value closest to 0 amplitude + unwrapped_phase = unwrapped_phase + (phase_n[-1] - unwrapped_phase[-1]) + + # Store new data + tmp = np.empty((len(amplitudes), len(columns)), dtype=object) + tmp[:, columns.index("xval")] = amplitudes + tmp[:, columns.index("yval")] = unwrapped_phase / self._freq_phase_coef() + tmp[:, columns.index("yerr")] = phase_s / self._freq_phase_coef() + tmp[:, columns.index("name")] = f"FREQ{direction}" + tmp[:, columns.index("class_id")] = m_id + tmp[:, columns.index("shots")] = x_quadrature.shots + y_quadrature.shots + tmp[:, columns.index("category")] = category + phase_data = np.r_[phase_data, tmp] + + return curve_data.append_list_values(other=phase_data) + + 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(c2_pos=(0, np.inf), c2_neg=(-np.inf, 0)) + user_opt.p0.set_if_empty( + c1_pos=0, c2_pos=1e6, c3_pos=0, c1_neg=0, c2_neg=-1e6, c3_neg=0, f_err=0 + ) + return user_opt + + def _create_analysis_results( + self, + fit_data: curve.CurveFitResult, + quality: str, + **metadata, + ) -> List[AnalysisResultData]: + outcomes = super()._create_analysis_results(fit_data, quality, **metadata) + + # Combine fit coefficients + coeffs = StarkCoefficients( + pos_coef_o1=fit_data.ufloat_params["c1_pos"].nominal_value, + pos_coef_o2=fit_data.ufloat_params["c2_pos"].nominal_value, + pos_coef_o3=fit_data.ufloat_params["c3_pos"].nominal_value, + neg_coef_o1=fit_data.ufloat_params["c1_neg"].nominal_value, + neg_coef_o2=fit_data.ufloat_params["c2_neg"].nominal_value, + neg_coef_o3=fit_data.ufloat_params["c3_neg"].nominal_value, + offset=fit_data.ufloat_params["f_err"].nominal_value, + ) + outcomes.append( + AnalysisResultData( + name="stark_coefficients", + value=coeffs, + chisq=fit_data.reduced_chisq, + quality=quality, + extra=metadata, + ) + ) + return outcomes + + def _create_figures( + self, + curve_data: curve.ScatterTable, + ) -> List["matplotlib.figure.Figure"]: + + # plot unwrapped phase on first axis + for d in ("pos", "neg"): + sub_data = curve_data[(curve_data.name == f"FREQ{d}") & (curve_data.category == "freq")] + self.plotter.set_series_data( + series_name=f"F{d}", + x_formatted=sub_data.xval.to_numpy(), + y_formatted=sub_data.yval.to_numpy(), + y_formatted_err=sub_data.yerr.to_numpy(), + ) + + # plot raw RamseyXY plot on second axis + for name in ("Xpos", "Ypos", "Xneg", "Yneg"): + sub_data = curve_data[(curve_data.name == name) & (curve_data.category == "ramsey_xy")] + self.plotter.set_series_data( + series_name=name, + x_formatted=sub_data.xval.to_numpy(), + y_formatted=sub_data.yval.to_numpy(), + y_formatted_err=sub_data.yerr.to_numpy(), + ) + + # find base and amplitude guess + ramsey_xy = curve_data[curve_data.category == "ramsey_xy"] + offset_guess = 0.5 * (ramsey_xy.yval.min() + ramsey_xy.yval.max()) + amp_guess = 0.5 * np.ptp(ramsey_xy.yval) + + # plot frequency and Ramsey fit lines + line_data = curve_data[curve_data.category == "fitted"] + for direction in ("pos", "neg"): + sub_data = line_data[line_data.name == f"FREQ{direction}"] + if len(sub_data) == 0: + continue + xval = sub_data.xval.to_numpy() + yn = sub_data.yval.to_numpy() + ys = sub_data.yerr.to_numpy() + yval = unp.uarray(yn, ys) * self._freq_phase_coef() + + # Ramsey fit lines are predicted from the phase fit line. + # Note that this line doesn't need to match with the expeirment data + # because Ramsey P1 data may fluctuate due to phase damping. + + # pylint: disable=no-member + ramsey_cos = amp_guess * unp.cos(yval) + offset_guess + ramsey_sin = amp_guess * unp.sin(yval) + offset_guess + + self.plotter.set_series_data( + series_name=f"F{direction}", + x_interp=xval, + y_interp=yn, + ) + self.plotter.set_series_data( + series_name=f"X{direction}", + x_interp=xval, + y_interp=unp.nominal_values(ramsey_cos), + ) + self.plotter.set_series_data( + series_name=f"Y{direction}", + x_interp=xval, + y_interp=unp.nominal_values(ramsey_sin), + ) + + if np.isfinite(ys).all(): + self.plotter.set_series_data( + series_name=f"F{direction}", + y_interp_err=ys, + ) + self.plotter.set_series_data( + series_name=f"X{direction}", + y_interp_err=unp.std_devs(ramsey_cos), + ) + self.plotter.set_series_data( + series_name=f"Y{direction}", + y_interp_err=unp.std_devs(ramsey_sin), + ) + return [self.plotter.figure()] + + def _initialize( + self, + experiment_data: ExperimentData, + ): + super()._initialize(experiment_data) + + # Set scaling factor to convert phase to frequency + if "stark_length" in experiment_data.metadata: + self.set_options(pulse_len=experiment_data.metadata["stark_length"]) + + +class StarkP1SpectAnalysis(BaseAnalysis): + """Analysis class for StarkP1Spectroscopy. + + # section: overview + + The P1 landscape is hardly predictable because of the random appearance of + lossy TLS notches, and hence this analysis doesn't provide any + generic mathematical model to fit the measurement data. + A developer may subclass this to conduct own analysis. + + This analysis just visualizes the measured P1 values against Stark tone amplitudes. + The tone amplitudes can be converted into the amount of Stark shift + when the calibrated coefficients are provided in the analysis option, + or the calibration experiment results are available in the result database. + + # section: see_also + :class:`qiskit_experiments.library.driven_freq_tuning.StarkRamseyXYAmpScan` + + """ + + @property + def plotter(self) -> vis.CurvePlotter: + """Curve plotter instance.""" + return self.options.plotter + + @classmethod + def _default_options(cls) -> Options: + """Default analysis options. + + Analysis Options: + plotter (Plotter): Plotter to visualize P1 landscape. + data_processor (DataProcessor): Data processor to compute P1 value. + stark_coefficients (Union[Dict, str]): Dictionary of Stark shift coefficients to + convert tone amplitudes into amount of Stark shift. This dictionary must include + all keys defined in :attr:`.StarkP1SpectAnalysis.stark_coefficients_names`, + which are calibrated with :class:`.StarkRamseyXYAmpScan`. + Alternatively, it searches for these coefficients in the result database + when "latest" is set. This requires having the experiment service set in + the experiment data to analyze. + x_key (str): Key of the circuit metadata to represent x value. + """ + options = super()._default_options() + + p1spect_plotter = vis.CurvePlotter(vis.MplDrawer()) + p1spect_plotter.set_figure_options( + xlabel="Stark amplitude", + ylabel="P(1)", + xscale="quadratic", + ) + + options.update_options( + plotter=p1spect_plotter, + data_processor=dp.DataProcessor("counts", [dp.Probability("1")]), + stark_coefficients=None, + x_key="xval", + ) + options.set_validator("stark_coefficients", StarkCoefficients) + + return options + + # pylint: disable=unused-argument + def _run_spect_analysis( + self, + xdata: np.ndarray, + ydata: np.ndarray, + ydata_err: np.ndarray, + ) -> list[AnalysisResultData]: + """Run further analysis on the spectroscopy data. + + .. note:: + A subclass can overwrite this method to conduct analysis. + + Args: + xdata: X values. This is either amplitudes or frequencies. + ydata: Y values. This is P1 values measured at different Stark tones. + ydata_err: Sampling error of the Y values. + + Returns: + A list of analysis results. + """ + return [] + + def _run_analysis( + self, + experiment_data: ExperimentData, + ) -> tuple[list[AnalysisResultData], list["matplotlib.figure.Figure"]]: + + x_key = self.options.x_key + + # Get calibrated Stark tone coefficients + if self.options.stark_coefficients is None and experiment_data.service is not None: + # Get value from service + stark_coeffs = retrieve_coefficients_from_service( + service=experiment_data.service, + backend_name=experiment_data.backend_name, + qubit=experiment_data.metadata["physical_qubits"][0], + ) + else: + stark_coeffs = self.options.stark_coefficients + + # Compute P1 value and sampling error + data = experiment_data.data() + try: + xdata = np.asarray([datum["metadata"][x_key] for datum in data], dtype=float) + except KeyError as ex: + raise DataProcessorError( + f"X value key {x_key} is not defined in circuit metadata." + ) from ex + ydata_ufloat = self.options.data_processor(data) + ydata = unp.nominal_values(ydata_ufloat) + ydata_err = unp.std_devs(ydata_ufloat) + + # Convert x-axis of amplitudes into Stark shift by consuming calibrated parameters. + if isinstance(stark_coeffs, StarkCoefficients): + xdata = convert_amp_to_freq( + amps=xdata, + coeffs=stark_coeffs, + ) + self.plotter.set_figure_options( + xlabel="Stark shift", + xval_unit="Hz", + xscale="linear", + ) + + # Draw figures and create analysis results. + self.plotter.set_series_data( + series_name="stark_p1", + x_formatted=xdata, + y_formatted=ydata, + y_formatted_err=ydata_err, + x_interp=xdata, + y_interp=ydata, + ) + analysis_results = self._run_spect_analysis(xdata, ydata, ydata_err) + + return analysis_results, [self.plotter.figure()] diff --git a/qiskit_experiments/library/driven_freq_tuning/coefficient_utils.py b/qiskit_experiments/library/driven_freq_tuning/coefficient_utils.py new file mode 100644 index 0000000000..5d47a88faa --- /dev/null +++ b/qiskit_experiments/library/driven_freq_tuning/coefficient_utils.py @@ -0,0 +1,232 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2024. +# +# 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. +"""Coefficients characterizing Stark shift.""" + +from dataclasses import dataclass + +import numpy as np + +from qiskit.providers.backend import Backend +from qiskit_ibm_experiment.service import IBMExperimentService +from qiskit_ibm_experiment.exceptions import IBMApiError + +from qiskit_experiments.framework.backend_data import BackendData +from qiskit_experiments.framework.experiment_data import ExperimentData + + +@dataclass +class StarkCoefficients: + """A dataclass representing a set of coefficients characterizing Stark shift.""" + + pos_coef_o1: float + """The first order shift coefficient on positive amplitude.""" + + pos_coef_o2: float + """The second order shift coefficient on positive amplitude.""" + + pos_coef_o3: float + """The third order shift coefficient on positive amplitude.""" + + neg_coef_o1: float + """The first order shift coefficient on negative amplitude.""" + + neg_coef_o2: float + """The second order shift coefficient on negative amplitude.""" + + neg_coef_o3: float + """The third order shift coefficient on negative amplitude.""" + + offset: float + """Offset frequency.""" + + def positive_coeffs(self) -> list[float]: + """Positive coefficients.""" + return [self.pos_coef_o3, self.pos_coef_o2, self.pos_coef_o1] + + def negative_coeffs(self) -> list[float]: + """Negative coefficients.""" + return [self.neg_coef_o3, self.neg_coef_o2, self.neg_coef_o1] + + def __str__(self): + # Short representation for dataframe + return "StarkCoefficients(...)" + + +def convert_freq_to_amp( + freqs: np.ndarray, + coeffs: StarkCoefficients, +) -> np.ndarray: + """A helper function to convert Stark frequency to amplitude. + + Args: + freqs: Target frequency shifts to compute required Stark amplitude. + coeffs: Calibrated Stark coefficients. + + Returns: + Estimated Stark amplitudes to induce input frequency shifts. + + Raises: + ValueError: When amplitude value cannot be solved. + """ + amplitudes = np.zeros_like(freqs) + for idx, freq in enumerate(freqs): + shift = freq - coeffs.offset + if np.isclose(shift, 0.0): + amplitudes[idx] = 0.0 + continue + if shift > 0: + fit = [*coeffs.positive_coeffs(), -shift] + else: + fit = [*coeffs.negative_coeffs(), -shift] + amp_candidates = np.roots(fit) + # Because the fit function is third order, we get three solutions here. + # Only one valid solution must exist because we assume + # a monotonic trend for Stark shift against tone amplitude in domain of definition. + criteria = np.all( + [ + # Frequency shift and tone have the same sign by definition + np.sign(amp_candidates.real) == np.sign(shift), + # Tone amplitude is a real value + np.isclose(amp_candidates.imag, 0.0), + # The absolute value of tone amplitude must be less than 1.0 + 10 mp + np.abs(amp_candidates.real) < 1.0 + 10 * np.finfo(float).eps, + ], + axis=0, + ) + valid_amps = amp_candidates[criteria] + if len(valid_amps) == 0: + raise ValueError(f"Stark shift at frequency value of {freq} Hz is not available.") + if len(valid_amps) > 1: + # We assume a monotonic trend but sometimes a large third-order term causes + # inflection point and inverts the trend in larger amplitudes. + # In this case we would have more than one solutions, but we can + # take the smallerst amplitude before reaching to the inflection point. + before_inflection = np.argmin(np.abs(valid_amps.real)) + valid_amp = float(valid_amps[before_inflection].real) + else: + valid_amp = float(valid_amps[0].real) + amplitudes[idx] = min(valid_amp, 1.0) + return amplitudes + + +def convert_amp_to_freq( + amps: np.ndarray, + coeffs: StarkCoefficients, +) -> np.ndarray: + """A helper function to convert Stark amplitude to frequency shift. + + Args: + amps: Amplitude values to convert into frequency shift. + coeffs: Calibrated Stark coefficients. + + Returns: + Calculated frequency shift at given Stark amplitude. + """ + pos_fit = np.poly1d([*coeffs.positive_coeffs(), coeffs.offset]) + neg_fit = np.poly1d([*coeffs.negative_coeffs(), coeffs.offset]) + + return np.where(amps > 0, pos_fit(amps), neg_fit(amps)) + + +def find_min_max_frequency( + min_amp: float, + max_amp: float, + coeffs: StarkCoefficients, +) -> tuple[float, float]: + """A helper function to estimate maximum frequency shift within given amplitude budget. + + Args: + min_amp: Minimum Stark amplitude. + max_amp: Maximum Stark amplitude. + coeffs: Calibrated Stark coefficients. + + Returns: + Minimum and maximum frequency shift available within the amplitude range. + """ + trim_amps = [] + for amp in [min_amp, max_amp]: + if amp > 0: + fit = coeffs.positive_coeffs() + else: + fit = coeffs.negative_coeffs() + # Solve for inflection points by computing the point where derivative becomes zero. + solutions = np.roots([deriv * coeff for deriv, coeff in zip((3, 2, 1), fit)]) + inflection_points = solutions[(solutions.imag == 0) & (np.sign(solutions) == np.sign(amp))] + if len(inflection_points) > 0: + # When multiple inflection points are found, use the most outer one. + # There could be a small inflection point around amp=0, + # when the first order term is significant. + amp = min([amp, max(inflection_points, key=abs)], key=abs) + trim_amps.append(amp) + return tuple(convert_amp_to_freq(np.asarray(trim_amps), coeffs)) + + +def retrieve_coefficients_from_service( + service: IBMExperimentService, + backend_name: str, + qubit: int, +) -> StarkCoefficients: + """Retrieve StarkCoefficients object from experiment service. + + Args: + service: IBM Experiment service instance interfacing with result database. + backend_name: Name of target backend. + qubit: Index of qubit. + + Returns: + StarkCoefficients object. + + Raises: + RuntimeError: When stark_coefficients entry doesn't exist in the service. + """ + try: + retrieved = service.analysis_results( + device_components=[f"Q{qubit}"], + result_type="stark_coefficients", + backend_name=backend_name, + sort_by=["creation_datetime:desc"], + ) + except (IBMApiError, ValueError) as ex: + raise RuntimeError( + f"Failed to retrieve the result of stark_coefficients: {ex.message}" + ) from ex + if len(retrieved) == 0: + raise RuntimeError( + "Analysis result of stark_coefficients is not found in the " + "experiment service. Run and save the result of StarkRamseyXYAmpScan." + ) + return retrieved[0].result_data["value"] + + +def retrieve_coefficients_from_backend( + backend: Backend, + qubit: int, +) -> StarkCoefficients: + """Retrieve StarkCoefficients object from the Qiskit backend. + + Args: + backend: Qiskit backend object. + qubit: Index of qubit. + + Returns: + StarkCoefficients object. + + Raises: + RuntimeError: When experiment service cannot be loaded from backend. + """ + name = BackendData(backend).name + service = ExperimentData.get_service_from_backend(backend) + + if service is None: + raise RuntimeError(f"Valid experiment service is not found for the backend {name}.") + + return retrieve_coefficients_from_service(service, name, qubit) diff --git a/qiskit_experiments/library/driven_freq_tuning/p1_spect.py b/qiskit_experiments/library/driven_freq_tuning/p1_spect.py new file mode 100644 index 0000000000..a104087853 --- /dev/null +++ b/qiskit_experiments/library/driven_freq_tuning/p1_spect.py @@ -0,0 +1,270 @@ +# 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. +"""P1 experiment at various qubit frequencies.""" + +from __future__ import annotations + +from collections.abc import Sequence + +import numpy as np +from qiskit import pulse +from qiskit.circuit import QuantumCircuit, Gate, Parameter, ParameterExpression +from qiskit.providers.backend import Backend +from qiskit.utils import optionals as _optional + +from qiskit_experiments.framework import BackendTiming, BaseExperiment, Options +from .analyses import StarkP1SpectAnalysis + +from .coefficient_utils import ( + StarkCoefficients, + convert_freq_to_amp, + retrieve_coefficients_from_backend, +) + +if _optional.HAS_SYMENGINE: + import symengine as sym +else: + import sympy as sym + + +class StarkP1Spectroscopy(BaseExperiment): + """P1 spectroscopy experiment with Stark tone. + + # section: overview + + This experiment measures a probability of the excitation state of the qubit + with a certain delay after excitation. + A Stark tone is applied during this delay to move the + qubit frequency to conduct a spectroscopy of qubit relaxation quantity. + + .. parsed-literal:: + + ┌───┐┌──────────────────┐┌─┐ + q: ┤ X ├┤ Stark(stark_amp) ├┤M├ + └───┘└──────────────────┘└╥┘ + c: 1/══════════════════════════╩═ + 0 + + Since the qubit relaxation rate may depend on the qubit frequency due to the + coupling to nearby energy levels, this experiment is useful to find out + lossy operation frequency that might be harmful to the gate fidelity [1]. + + # section: analysis_ref + :class:`qiskit_experiments.library.driven_freq_tuning.StarkP1SpectAnalysis` + + # section: reference + .. ref_arxiv:: 1 2105.15201 + + # section: see_also + :class:`qiskit_experiments.library.driven_freq_tuning.ramsey.StarkRamseyXY` + :class:`qiskit_experiments.library.driven_freq_tuning.ramsey_amp_scan.StarkRamseyXYAmpScan` + + # section: manual + :doc:`/manuals/characterization/stark_experiment` + + """ + + def __init__( + self, + physical_qubits: Sequence[int], + backend: Backend | None = None, + **experiment_options, + ): + """ + Initialize new experiment class. + + Args: + physical_qubits: Sequence with the index of the physical qubit. + backend: Optional, the backend to run the experiment on. + experiment_options: Experiment options. See the class documentation or + ``self._default_experiment_options`` for descriptions. + """ + self._timing = None + + super().__init__( + physical_qubits=physical_qubits, + analysis=StarkP1SpectAnalysis(), + backend=backend, + ) + self.set_experiment_options(**experiment_options) + + @classmethod + def _default_experiment_options(cls) -> Options: + """Default experiment options. + + Experiment Options: + t1_delay (float): The T1 delay time after excitation pulse. The delay must be + sufficiently greater than the edge duration determined by the stark_sigma. + stark_channel (PulseChannel): Pulse channel to apply Stark tones. + If not provided, the same channel with the qubit drive is used. + stark_freq_offset (float): Offset of Stark tone frequency from the qubit frequency. + This must be greater than zero not to apply Rabi drive. + stark_sigma (float): Gaussian sigma of the rising and falling edges + of the Stark tone, in seconds. + stark_risefall (float): Ratio of sigma to the duration of + the rising and falling edges of the Stark tone. + min_xval (float): Minimum x value. + max_xval (float): Maximum x value. + num_xvals (int): Number of x-values to scan. + xval_type (str): Type of x-value. Either ``amplitude`` or ``frequency``. + Setting to frequency requires pre-calibration of Stark shift coefficients. + spacing (str): A policy for the spacing to create an amplitude list from + ``min_stark_amp`` to ``max_stark_amp``. Either ``linear`` or ``quadratic`` + must be specified. + xvals (list[float]): The list of x-values that will be scanned in the experiment. + If not set, then ``num_xvals`` parameters spaced according to + the ``spacing`` policy between ``min_xval`` and ``max_xval`` are used. + If ``xvals`` is set, these parameters are ignored. + stark_coefficients (StarkCoefficients): Calibrated Stark shift coefficients. + This value is necessary when xval_type is "frequency". + When this value is None, a search for the "stark_coefficients" in the + result database is run. This requires having the experiment service + available in the backend set for the experiment. + """ + options = super()._default_experiment_options() + options.update_options( + t1_delay=20e-6, + stark_channel=None, + stark_freq_offset=80e6, + stark_sigma=15e-9, + stark_risefall=2, + min_xval=-1.0, + max_xval=1.0, + num_xvals=201, + xval_type="amplitude", + spacing="quadratic", + xvals=None, + stark_coefficients=None, + ) + options.set_validator("spacing", ["linear", "quadratic"]) + options.set_validator("xval_type", ["amplitude", "frequency"]) + options.set_validator("stark_freq_offset", (0, np.inf)) + options.set_validator("stark_channel", pulse.channels.PulseChannel) + options.set_validator("stark_coefficients", StarkCoefficients) + return options + + def _set_backend(self, backend: Backend): + super()._set_backend(backend) + self._timing = BackendTiming(backend) + + def parameters(self) -> np.ndarray: + """Stark tone amplitudes to use in circuits. + + Returns: + The list of amplitudes to use for the different circuits based on the + experiment options. + + Raises: + ValueError: When invalid xval spacing is specified. + """ + opt = self.experiment_options # alias + + if opt.xvals is None: + if opt.spacing == "linear": + params = np.linspace(opt.min_xval, opt.max_xval, opt.num_xvals) + elif opt.spacing == "quadratic": + min_sqrt = np.sign(opt.min_xval) * np.sqrt(np.abs(opt.min_xval)) + max_sqrt = np.sign(opt.max_xval) * np.sqrt(np.abs(opt.max_xval)) + lin_params = np.linspace(min_sqrt, max_sqrt, opt.num_xvals) + params = np.sign(lin_params) * lin_params**2 + else: + raise ValueError(f"Spacing option {opt.spacing} is not valid.") + else: + params = np.asarray(opt.xvals, dtype=float) + + if opt.xval_type == "frequency": + coeffs = opt.stark_coefficients + if coeffs is None: + coeffs = retrieve_coefficients_from_backend( + backend=self.backend, + qubit=self.physical_qubits[0], + ) + return convert_freq_to_amp(freqs=params, coeffs=coeffs) + return params + + def parameterized_circuits(self) -> tuple[QuantumCircuit, ...]: + """Create circuits with parameters for P1 experiment with Stark shift. + + Returns: + Quantum template circuit for P1 experiment. + """ + opt = self.experiment_options # alias + param = Parameter("stark_amp") + sym_param = param._symbol_expr + + # Pulse gates + stark = Gate("Stark", 1, [param]) + + # Note that Stark tone yields negative (positive) frequency shift + # when the Stark tone frequency is higher (lower) than qubit f01 frequency. + # This choice gives positive frequency shift with positive Stark amplitude. + qubit_f01 = self._backend_data.drive_freqs[self.physical_qubits[0]] + neg_sign_of_amp = ParameterExpression( + symbol_map={param: sym_param}, + expr=-sym.sign(sym_param), + ) + abs_of_amp = ParameterExpression( + symbol_map={param: sym_param}, + expr=sym.Abs(sym_param), + ) + stark_freq = qubit_f01 + neg_sign_of_amp * opt.stark_freq_offset + stark_channel = opt.stark_channel or pulse.DriveChannel(self.physical_qubits[0]) + sigma_dt = opt.stark_sigma / self._backend_data.dt + delay_dt = self._timing.round_pulse(time=opt.t1_delay) + + with pulse.build() as stark_schedule: + pulse.set_frequency(stark_freq, stark_channel) + pulse.play( + pulse.GaussianSquare( + duration=delay_dt, + amp=abs_of_amp, + sigma=sigma_dt, + risefall_sigma_ratio=opt.stark_risefall, + ), + stark_channel, + ) + + temp_t1 = QuantumCircuit(1, 1) + temp_t1.x(0) + temp_t1.append(stark, [0]) + temp_t1.measure(0, 0) + temp_t1.add_calibration( + gate=stark, + qubits=self.physical_qubits, + schedule=stark_schedule, + ) + + return (temp_t1,) + + def circuits(self) -> list[QuantumCircuit]: + """Create circuits. + + Returns: + A list of P1 circuits with a variable Stark tone amplitudes. + """ + (t1_circ,) = self.parameterized_circuits() + param = next(iter(t1_circ.parameters)) + + circs = [] + for amp in self.parameters(): + t1_assigned = t1_circ.assign_parameters({param: amp}, inplace=False) + t1_assigned.metadata = {"xval": amp} + circs.append(t1_assigned) + + return circs + + def _metadata(self) -> dict[str, any]: + """Return experiment metadata for ExperimentData.""" + metadata = super()._metadata() + metadata["stark_freq_offset"] = self.experiment_options.stark_freq_offset + + return metadata diff --git a/qiskit_experiments/library/driven_freq_tuning/ramsey.py b/qiskit_experiments/library/driven_freq_tuning/ramsey.py new file mode 100644 index 0000000000..e214a5c4b4 --- /dev/null +++ b/qiskit_experiments/library/driven_freq_tuning/ramsey.py @@ -0,0 +1,359 @@ +# 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. +"""Stark Ramsey experiment.""" + +from __future__ import annotations + +import warnings +from collections.abc import Sequence + +import numpy as np +from qiskit import pulse +from qiskit.circuit import QuantumCircuit, Gate, Parameter +from qiskit.providers.backend import Backend +from qiskit.utils import optionals as _optional + +from qiskit_experiments.framework import BaseExperiment, Options, BackendTiming +from qiskit_experiments.library.characterization.analysis import RamseyXYAnalysis + +if _optional.HAS_SYMENGINE: + pass +else: + pass + + +class StarkRamseyXY(BaseExperiment): + """A sign-sensitive experiment to measure the frequency of a qubit under a pulsed Stark tone. + + # section: overview + + This experiment is a variant of :class:`.RamseyXY` with a pulsed Stark tone + and consists of the following two circuits: + + .. parsed-literal:: + + (Ramsey X) The pulse before measurement rotates by pi-half around the X axis + + ┌────┐┌────────┐┌───┐┌───────────────┐┌────────┐┌────┐┌─┐ + q: ┤ √X ├┤ StarkV ├┤ X ├┤ StarkU(delay) ├┤ Rz(-π) ├┤ √X ├┤M├ + └────┘└────────┘└───┘└───────────────┘└────────┘└────┘└╥┘ + c: 1/═══════════════════════════════════════════════════════╩═ + 0 + + (Ramsey Y) The pulse before measurement rotates by pi-half around the Y axis + + ┌────┐┌────────┐┌───┐┌───────────────┐┌───────────┐┌────┐┌─┐ + q: ┤ √X ├┤ StarkV ├┤ X ├┤ StarkU(delay) ├┤ Rz(-3π/2) ├┤ √X ├┤M├ + └────┘└────────┘└───┘└───────────────┘└───────────┘└────┘└╥┘ + c: 1/══════════════════════════════════════════════════════════╩═ + 0 + + In principle, the sequence is a variant of :class:`.RamseyXY` circuit. + However, the delay in between √X gates is replaced with an off-resonant drive. + This off-resonant drive shifts the qubit frequency due to the + Stark effect and causes it to accumulate phase during the + Ramsey sequence. This frequency shift is a function of the + offset of the Stark tone frequency from the qubit frequency + and of the magnitude of the tone. + + Note that the Stark tone pulse (StarkU) takes the form of a flat-topped Gaussian envelope. + The magnitude of the pulse varies in time during its rising and falling edges. + It is difficult to characterize the net phase accumulation of the qubit during the + edges of the pulse when the frequency shift is varying with the pulse amplitude. + In order to simplify the analysis, an additional pulse (StarkV) + involving only the edges of StarkU is added to the sequence. + The sign of the phase accumulation is inverted for StarkV relative to that of StarkU + by inserting an X gate in between the two pulses. + + This technique allows the experiment to accumulate only the net phase + during the flat-top part of the StarkU pulse with constant magnitude. + + # section: analysis_ref + :class:`qiskit_experiments.library.characterization.RamseyXYAnalysis` + + # section: see_also + :class:`qiskit_experiments.library.characterization.ramsey_xy.RamseyXY` + + # section: manual + :doc:`/manuals/characterization/stark_experiment` + + """ + + def __init__( + self, + physical_qubits: Sequence[int], + backend: Backend | None = None, + **experiment_options, + ): + """Create new experiment. + + Args: + physical_qubits: Index of physical qubit. + backend: Optional, the backend to run the experiment on. + experiment_options: Experiment options. See the class documentation or + ``self._default_experiment_options`` for descriptions. + """ + self._timing = None + + super().__init__( + physical_qubits=physical_qubits, + analysis=RamseyXYAnalysis(), + backend=backend, + ) + self.set_experiment_options(**experiment_options) + + @classmethod + def _default_experiment_options(cls) -> Options: + """Default experiment options. + + Experiment Options: + stark_amp (float): A single float parameter to represent the magnitude of Stark tone + and the sign of expected Stark shift. + See :ref:`stark_tone_implementation` for details. + stark_channel (PulseChannel): Pulse channel to apply Stark tones. + If not provided, the same channel with the qubit drive is used. + See :ref:`stark_channel_consideration` for details. + stark_freq_offset (float): Offset of Stark tone frequency from the qubit frequency. + This offset should be sufficiently large so that the Stark pulse + does not Rabi drive the qubit. + See :ref:`stark_frequency_consideration` for details. + stark_sigma (float): Gaussian sigma of the rising and falling edges + of the Stark tone, in seconds. + stark_risefall (float): Ratio of sigma to the duration of + the rising and falling edges of the Stark tone. + min_freq (float): Minimum frequency that this experiment is guaranteed to resolve. + Note that fitter algorithm :class:`.RamseyXYAnalysis` of this experiment + is still capable of fitting experiment data with lower frequency. + max_freq (float): Maximum frequency that this experiment can resolve. + delays (list[float]): The list of delays if set that will be scanned in the + experiment. If not set, then evenly spaced delays with interval + computed from ``min_freq`` and ``max_freq`` are used. + See :meth:`StarkRamseyXY.delays` for details. + """ + options = super()._default_experiment_options() + options.update_options( + stark_amp=0.0, + stark_channel=None, + stark_freq_offset=80e6, + stark_sigma=15e-9, + stark_risefall=2, + min_freq=5e6, + max_freq=100e6, + delays=None, + ) + options.set_validator("stark_freq_offset", (0, np.inf)) + options.set_validator("stark_channel", pulse.channels.PulseChannel) + return options + + def _set_backend(self, backend: Backend): + super()._set_backend(backend) + self._timing = BackendTiming(backend) + + def set_experiment_options(self, **fields): + _warning_circuit_length = 300 + + # Do validation for circuit number + min_freq = fields.get("min_freq", None) + max_freq = fields.get("max_freq", None) + delays = fields.get("delays", None) + if min_freq is not None and max_freq is not None: + if delays: + warnings.warn( + "Experiment option 'min_freq' and 'max_freq' are ignored " + "when 'delays' are explicitly specified.", + UserWarning, + ) + else: + n_expr_circs = 2 * int(2 * max_freq / min_freq) # delays * (x, y) + max_circs_per_job = None + if self._backend_data: + max_circs_per_job = self._backend_data.max_circuits() + if n_expr_circs > (max_circs_per_job or _warning_circuit_length): + warnings.warn( + f"Provided configuration generates {n_expr_circs} circuits. " + "You can set smaller 'max_freq' or larger 'min_freq' to reduce this number. " + "This experiment is still executable but your execution may consume " + "unnecessary long quantum device time, and result may suffer " + "device parameter drift in consequence of the long execution time.", + UserWarning, + ) + # Do validation for spectrum overlap to avoid real excitation + stark_freq_offset = fields.get("stark_freq_offset", None) + stark_sigma = fields.get("stark_sigma", None) + if stark_freq_offset is not None and stark_sigma is not None: + if stark_freq_offset < 1 / stark_sigma: + warnings.warn( + "Provided configuration may induce coherent state exchange between qubit levels " + "because of the potential spectrum overlap. You can avoid this by " + "increasing the 'stark_sigma' or 'stark_freq_offset'. " + "Note that this experiment is still executable.", + UserWarning, + ) + pass + + super().set_experiment_options(**fields) + + def parameters(self) -> np.ndarray: + """Delay values to use in circuits. + + .. note:: + + The delays are computed with the ``min_freq`` and ``max_freq`` experiment options. + The maximum point is computed from the ``min_freq`` to guarantee the result + contains at least one Ramsey oscillation cycle at this frequency. + The interval is computed from the ``max_freq`` to sample with resolution + such that the Nyquist frequency is twice ``max_freq``. + + Returns: + The list of delays to use for the different circuits based on the + experiment options. + + Raises: + ValueError: When ``min_freq`` is larger than ``max_freq``. + """ + opt = self.experiment_options # alias + + if opt.delays is None: + if opt.min_freq > opt.max_freq: + raise ValueError("Experiment option 'min_freq' must be smaller than 'max_freq'.") + # Delay is longer enough to capture 1 cycle of the minimum frequency. + # Fitter can still accurately fit samples shorter than 1 cycle. + max_period = 1 / opt.min_freq + # Inverse of interval should be greater than Nyquist frequency. + sampling_freq = 2 * opt.max_freq + interval = 1 / sampling_freq + return np.arange(0, max_period, interval) + return opt.delays + + def parameterized_circuits(self) -> tuple[QuantumCircuit, ...]: + """Create circuits with parameters for Ramsey XY experiment with Stark tone. + + Returns: + Quantum template circuits for Ramsey X and Ramsey Y experiment. + """ + opt = self.experiment_options # alias + param = Parameter("delay") + + # Pulse gates + stark_v = Gate("StarkV", 1, []) + stark_u = Gate("StarkU", 1, [param]) + + # Note that Stark tone yields negative (positive) frequency shift + # when the Stark tone frequency is higher (lower) than qubit f01 frequency. + # This choice gives positive frequency shift with positive Stark amplitude. + qubit_f01 = self._backend_data.drive_freqs[self.physical_qubits[0]] + stark_freq = qubit_f01 - np.sign(opt.stark_amp) * opt.stark_freq_offset + stark_amp = np.abs(opt.stark_amp) + stark_channel = opt.stark_channel or pulse.DriveChannel(self.physical_qubits[0]) + ramps_dt = self._timing.round_pulse(time=2 * opt.stark_risefall * opt.stark_sigma) + sigma_dt = ramps_dt / 2 / opt.stark_risefall + + with pulse.build() as stark_v_schedule: + pulse.set_frequency(stark_freq, stark_channel) + pulse.play( + pulse.Gaussian( + duration=ramps_dt, + amp=stark_amp, + sigma=sigma_dt, + name="StarkV", + ), + stark_channel, + ) + + with pulse.build() as stark_u_schedule: + pulse.set_frequency(stark_freq, stark_channel) + pulse.play( + pulse.GaussianSquare( + duration=ramps_dt + param, + amp=stark_amp, + sigma=sigma_dt, + risefall_sigma_ratio=opt.stark_risefall, + name="StarkU", + ), + stark_channel, + ) + + ram_x = QuantumCircuit(1, 1) + ram_x.sx(0) + ram_x.append(stark_v, [0]) + ram_x.x(0) + ram_x.append(stark_u, [0]) + ram_x.rz(-np.pi, 0) + ram_x.sx(0) + ram_x.measure(0, 0) + ram_x.metadata = {"series": "X"} + ram_x.add_calibration( + gate=stark_v, + qubits=self.physical_qubits, + schedule=stark_v_schedule, + ) + ram_x.add_calibration( + gate=stark_u, + qubits=self.physical_qubits, + schedule=stark_u_schedule, + ) + + ram_y = QuantumCircuit(1, 1) + ram_y.sx(0) + ram_y.append(stark_v, [0]) + ram_y.x(0) + ram_y.append(stark_u, [0]) + ram_y.rz(-np.pi * 3 / 2, 0) + ram_y.sx(0) + ram_y.measure(0, 0) + ram_y.metadata = {"series": "Y"} + ram_y.add_calibration( + gate=stark_v, + qubits=self.physical_qubits, + schedule=stark_v_schedule, + ) + ram_y.add_calibration( + gate=stark_u, + qubits=self.physical_qubits, + schedule=stark_u_schedule, + ) + + return ram_x, ram_y + + def circuits(self) -> list[QuantumCircuit]: + """Create circuits. + + Returns: + A list of circuits with a variable delay. + """ + timing = BackendTiming(self.backend, min_length=0) + + ramx_circ, ramy_circ = self.parameterized_circuits() + param = next(iter(ramx_circ.parameters)) + + circs = [] + for delay in self.parameters(): + valid_delay_dt = timing.round_pulse(time=delay) + net_delay_sec = timing.pulse_time(time=delay) + + ramx_circ_assigned = ramx_circ.assign_parameters({param: valid_delay_dt}, inplace=False) + ramx_circ_assigned.metadata["xval"] = net_delay_sec + + ramy_circ_assigned = ramy_circ.assign_parameters({param: valid_delay_dt}, inplace=False) + ramy_circ_assigned.metadata["xval"] = net_delay_sec + + circs.extend([ramx_circ_assigned, ramy_circ_assigned]) + + return circs + + def _metadata(self) -> dict[str, any]: + """Return experiment metadata for ExperimentData.""" + metadata = super()._metadata() + metadata["stark_amp"] = self.experiment_options.stark_amp + metadata["stark_freq_offset"] = self.experiment_options.stark_freq_offset + + return metadata diff --git a/qiskit_experiments/library/driven_freq_tuning/ramsey_amp_scan.py b/qiskit_experiments/library/driven_freq_tuning/ramsey_amp_scan.py new file mode 100644 index 0000000000..d04eb607d1 --- /dev/null +++ b/qiskit_experiments/library/driven_freq_tuning/ramsey_amp_scan.py @@ -0,0 +1,311 @@ +# 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. +"""Stark Ramsey experiment directly scanning Stark amplitude.""" + +from __future__ import annotations + +from collections.abc import Sequence + +import numpy as np +from qiskit import pulse +from qiskit.circuit import QuantumCircuit, Gate, ParameterExpression, Parameter +from qiskit.providers.backend import Backend +from qiskit.utils import optionals as _optional + +from qiskit_experiments.framework import BaseExperiment, Options, BackendTiming +from .analyses import StarkRamseyXYAmpScanAnalysis + +if _optional.HAS_SYMENGINE: + import symengine as sym +else: + import sympy as sym + + +class StarkRamseyXYAmpScan(BaseExperiment): + r"""A fast characterization of Stark frequency shift by varying Stark tone amplitude. + + # section: overview + + This experiment scans Stark tone amplitude at a fixed tone duration. + The experimental circuits are identical to the :class:`.StarkRamseyXY` experiment + except that the Stark pulse amplitude is the scanned parameter rather than the pulse width. + + .. parsed-literal:: + + (Ramsey X) The pulse before measurement rotates by pi-half around the X axis + + ┌────┐┌───────────────────┐┌───┐┌───────────────────┐┌────────┐┌────┐┌─┐ + q: ┤ √X ├┤ StarkV(stark_amp) ├┤ X ├┤ StarkU(stark_amp) ├┤ Rz(-π) ├┤ √X ├┤M├ + └────┘└───────────────────┘└───┘└───────────────────┘└────────┘└────┘└╥┘ + c: 1/══════════════════════════════════════════════════════════════════════╩═ + 0 + + (Ramsey Y) The pulse before measurement rotates by pi-half around the Y axis + + ┌────┐┌───────────────────┐┌───┐┌───────────────────┐┌───────────┐┌────┐┌─┐ + q: ┤ √X ├┤ StarkV(stark_amp) ├┤ X ├┤ StarkU(stark_amp) ├┤ Rz(-3π/2) ├┤ √X ├┤M├ + └────┘└───────────────────┘└───┘└───────────────────┘└───────────┘└────┘└╥┘ + c: 1/═════════════════════════════════════════════════════════════════════════╩═ + 0 + + The AC Stark effect can be used to shift the frequency of a qubit with a microwave drive. + To calibrate a specific frequency shift, the :class:`.StarkRamseyXY` experiment can be run + to scan the Stark pulse duration at every amplitude, but such a two dimensional scan of + the tone duration and amplitude may require many circuit executions. + To avoid this overhead, the :class:`.StarkRamseyXYAmpScan` experiment fixes the + tone duration and scans only amplitude. + + Recall that an observed Ramsey oscillation in each quadrature may be represented by + + .. math:: + + {\cal E}_X(\Omega, t_S) = A e^{-t_S/\tau} \cos \left( 2\pi f_S(\Omega) t_S \right), \\ + {\cal E}_Y(\Omega, t_S) = A e^{-t_S/\tau} \sin \left( 2\pi f_S(\Omega) t_S \right), + + where :math:`f_S(\Omega)` denotes the amount of Stark shift + at a constant tone amplitude :math:`\Omega`, and :math:`t_S` is the duration of the + applied tone. For a fixed tone duration, + one can still observe the Ramsey oscillation by scanning the tone amplitude. + However, since :math:`f_S` is usually a higher order polynomial of :math:`\Omega`, + one must manage to fit the y-data for trigonometric functions with + phase which non-linearly changes with the x-data. + The :class:`.StarkRamseyXYAmpScan` experiment thus drastically reduces the number of + circuits to run in return for greater complexity in the fitting model. + + # section: analysis_ref + :class:`qiskit_experiments.library.driven_freq_tuning.StarkRamseyXYAmpScanAnalysis` + + # section: see_also + :class:`qiskit_experiments.library.driven_freq_tuning.ramsey.StarkRamseyXY` + :class:`qiskit_experiments.library.characterization.ramsey_xy.RamseyXY` + + # section: manual + :doc:`/manuals/characterization/stark_experiment` + + """ + + def __init__( + self, + physical_qubits: Sequence[int], + backend: Backend | None = None, + **experiment_options, + ): + """Create new experiment. + + Args: + physical_qubits: Sequence with the index of the physical qubit. + backend: Optional, the backend to run the experiment on. + experiment_options: Experiment options. See the class documentation or + ``self._default_experiment_options`` for descriptions. + """ + self._timing = None + + super().__init__( + physical_qubits=physical_qubits, + analysis=StarkRamseyXYAmpScanAnalysis(), + backend=backend, + ) + self.set_experiment_options(**experiment_options) + + @classmethod + def _default_experiment_options(cls) -> Options: + """Default experiment options. + + Experiment Options: + stark_channel (PulseChannel): Pulse channel on which to apply Stark tones. + If not provided, the same channel with the qubit drive is used. + See :ref:`stark_channel_consideration` for details. + stark_freq_offset (float): Offset of Stark tone frequency from the qubit frequency. + This offset should be sufficiently large so that the Stark pulse + does not Rabi drive the qubit. + See :ref:`stark_frequency_consideration` for details. + stark_sigma (float): Gaussian sigma of the rising and falling edges + of the Stark tone, in seconds. + stark_risefall (float): Ratio of sigma to the duration of + the rising and falling edges of the Stark tone. + stark_length (float): Time to accumulate Stark shifted phase in seconds. + min_stark_amp (float): Minimum Stark tone amplitude. + max_stark_amp (float): Maximum Stark tone amplitude. + num_stark_amps (int): Number of Stark tone amplitudes to scan. + stark_amps (list[float]): The list of amplitude that will be scanned in the experiment. + If not set, then ``num_stark_amps`` evenly spaced amplitudes + between ``min_stark_amp`` and ``max_stark_amp`` are used. If ``stark_amps`` + is set, these parameters are ignored. + """ + options = super()._default_experiment_options() + options.update_options( + stark_channel=None, + stark_freq_offset=80e6, + stark_sigma=15e-9, + stark_risefall=2, + stark_length=50e-9, + min_stark_amp=-1.0, + max_stark_amp=1.0, + num_stark_amps=101, + stark_amps=None, + ) + options.set_validator("stark_freq_offset", (0, np.inf)) + options.set_validator("stark_channel", pulse.channels.PulseChannel) + return options + + def _set_backend(self, backend: Backend): + super()._set_backend(backend) + self._timing = BackendTiming(backend) + + def parameters(self) -> np.ndarray: + """Stark tone amplitudes to use in circuits. + + Returns: + The list of amplitudes to use for the different circuits based on the + experiment options. + """ + opt = self.experiment_options # alias + + if opt.stark_amps is None: + params = np.linspace(opt.min_stark_amp, opt.max_stark_amp, opt.num_stark_amps) + else: + params = np.asarray(opt.stark_amps, dtype=float) + + return params + + def parameterized_circuits(self) -> tuple[QuantumCircuit, ...]: + """Create circuits with parameters for Ramsey XY experiment with Stark tone. + + Returns: + Quantum template circuits for Ramsey X and Ramsey Y experiment. + """ + opt = self.experiment_options # alias + param = Parameter("stark_amp") + sym_param = param._symbol_expr + + # Pulse gates + stark_v = Gate("StarkV", 1, [param]) + stark_u = Gate("StarkU", 1, [param]) + + # Note that Stark tone yields negative (positive) frequency shift + # when the Stark tone frequency is higher (lower) than qubit f01 frequency. + # This choice gives positive frequency shift with positive Stark amplitude. + qubit_f01 = self._backend_data.drive_freqs[self.physical_qubits[0]] + neg_sign_of_amp = ParameterExpression( + symbol_map={param: sym_param}, + expr=-sym.sign(sym_param), + ) + abs_of_amp = ParameterExpression( + symbol_map={param: sym_param}, + expr=sym.Abs(sym_param), + ) + stark_freq = qubit_f01 + neg_sign_of_amp * opt.stark_freq_offset + stark_channel = opt.stark_channel or pulse.DriveChannel(self.physical_qubits[0]) + ramps_dt = self._timing.round_pulse(time=2 * opt.stark_risefall * opt.stark_sigma) + sigma_dt = ramps_dt / 2 / opt.stark_risefall + width_dt = self._timing.round_pulse(time=opt.stark_length) + + with pulse.build() as stark_v_schedule: + pulse.set_frequency(stark_freq, stark_channel) + pulse.play( + pulse.Gaussian( + duration=ramps_dt, + amp=abs_of_amp, + sigma=sigma_dt, + ), + stark_channel, + ) + + with pulse.build() as stark_u_schedule: + pulse.set_frequency(stark_freq, stark_channel) + pulse.play( + pulse.GaussianSquare( + duration=ramps_dt + width_dt, + amp=abs_of_amp, + sigma=sigma_dt, + risefall_sigma_ratio=opt.stark_risefall, + ), + stark_channel, + ) + + ram_x = QuantumCircuit(1, 1) + ram_x.sx(0) + ram_x.append(stark_v, [0]) + ram_x.x(0) + ram_x.append(stark_u, [0]) + ram_x.rz(-np.pi, 0) + ram_x.sx(0) + ram_x.measure(0, 0) + ram_x.metadata = {"series": "X"} + ram_x.add_calibration( + gate=stark_v, + qubits=self.physical_qubits, + schedule=stark_v_schedule, + ) + ram_x.add_calibration( + gate=stark_u, + qubits=self.physical_qubits, + schedule=stark_u_schedule, + ) + + ram_y = QuantumCircuit(1, 1) + ram_y.sx(0) + ram_y.append(stark_v, [0]) + ram_y.x(0) + ram_y.append(stark_u, [0]) + ram_y.rz(-np.pi * 3 / 2, 0) + ram_y.sx(0) + ram_y.measure(0, 0) + ram_y.metadata = {"series": "Y"} + ram_y.add_calibration( + gate=stark_v, + qubits=self.physical_qubits, + schedule=stark_v_schedule, + ) + ram_y.add_calibration( + gate=stark_u, + qubits=self.physical_qubits, + schedule=stark_u_schedule, + ) + + return ram_x, ram_y + + def circuits(self) -> list[QuantumCircuit]: + """Create circuits. + + Returns: + A list of circuits with a variable Stark tone amplitudes. + """ + ramx_circ, ramy_circ = self.parameterized_circuits() + param = next(iter(ramx_circ.parameters)) + + circs = [] + for amp in self.parameters(): + # Add metadata "direction" to ease the filtering of the data + # by curve analysis. Indeed, the fit parameters are amplitude sign dependent. + + ramx_circ_assigned = ramx_circ.assign_parameters({param: amp}, inplace=False) + ramx_circ_assigned.metadata["xval"] = amp + ramx_circ_assigned.metadata["direction"] = "pos" if amp > 0 else "neg" + + ramy_circ_assigned = ramy_circ.assign_parameters({param: amp}, inplace=False) + ramy_circ_assigned.metadata["xval"] = amp + ramy_circ_assigned.metadata["direction"] = "pos" if amp > 0 else "neg" + + circs.extend([ramx_circ_assigned, ramy_circ_assigned]) + + return circs + + def _metadata(self) -> dict[str, any]: + """Return experiment metadata for ExperimentData.""" + metadata = super()._metadata() + metadata["stark_length"] = self._timing.pulse_time( + time=self.experiment_options.stark_length + ) + metadata["stark_freq_offset"] = self.experiment_options.stark_freq_offset + + return metadata diff --git a/test/library/driven_freq_tuning/__init__.py b/test/library/driven_freq_tuning/__init__.py new file mode 100644 index 0000000000..4575d01965 --- /dev/null +++ b/test/library/driven_freq_tuning/__init__.py @@ -0,0 +1,12 @@ +# 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. +"""Tests for driven frequency tuning.""" diff --git a/test/library/characterization/test_stark_p1_spect.py b/test/library/driven_freq_tuning/test_stark_p1_spect.py similarity index 53% rename from test/library/characterization/test_stark_p1_spect.py rename to test/library/driven_freq_tuning/test_stark_p1_spect.py index 3ae4c17e76..f4fb1f2e13 100644 --- a/test/library/characterization/test_stark_p1_spect.py +++ b/test/library/driven_freq_tuning/test_stark_p1_spect.py @@ -23,7 +23,8 @@ from qiskit_experiments.framework import ExperimentData, AnalysisResultData from qiskit_experiments.library import StarkP1Spectroscopy -from qiskit_experiments.library.characterization.analysis import StarkP1SpectAnalysis +from qiskit_experiments.library.driven_freq_tuning.analyses import StarkP1SpectAnalysis +from qiskit_experiments.library.driven_freq_tuning import coefficient_utils as util from qiskit_experiments.test import FakeService @@ -48,49 +49,6 @@ def _run_spect_analysis( class TestStarkP1Spectroscopy(QiskitExperimentsTestCase): """Test case for the Stark P1 Spectroscopy experiment.""" - @staticmethod - def create_service_helper( - pos_coef_o1: float, - pos_coef_o2: float, - pos_coef_o3: float, - neg_coef_o1: float, - neg_coef_o2: float, - neg_coef_o3: float, - ferr: float, - qubit: int, - backend_name: str, - ): - """A helper method to create service with analysis results.""" - service = FakeService() - - service.create_experiment( - experiment_type="StarkRamseyXYAmpScan", - backend_name=backend_name, - experiment_id="123456789", - ) - - coeffs = { - "stark_pos_coef_o1": pos_coef_o1, - "stark_pos_coef_o2": pos_coef_o2, - "stark_pos_coef_o3": pos_coef_o3, - "stark_neg_coef_o1": neg_coef_o1, - "stark_neg_coef_o2": neg_coef_o2, - "stark_neg_coef_o3": neg_coef_o3, - "stark_ferr": ferr, - } - for i, (key, value) in enumerate(coeffs.items()): - service.create_analysis_result( - experiment_id="123456789", - result_data={"value": value}, - result_type=key, - device_components=[f"Q{qubit}"], - tags=[], - quality="Good", - verified=False, - result_id=str(i), - ) - return service - def test_linear_spaced_parameters(self): """Test generating parameters with linear spacing.""" exp = StarkP1Spectroscopy((0,)) @@ -126,96 +84,59 @@ def test_invalid_spacing(self): exp.set_experiment_options(spacing="invalid_option") def test_raises_scanning_frequency_without_service(self): - """Test raises error when frequency is set without having service or coefficients set.""" + """Test raises error when frequency is set without no coefficients. + + This covers following situations: + - stark_coefficients options is None + - backend object doesn't provide experiment service + """ exp = StarkP1Spectroscopy((0,), backend=FakeHanoiV2()) exp.set_experiment_options( xvals=[-100e6, -50e6, 0, 50e6, 100e6], xval_type="frequency", - stark_coefficients="latest", ) with self.assertRaises(RuntimeError): exp.parameters() - @named_data( - ["ordinary", 5e6, 200e6, -50e6, 5e6, -180e6, -40e6, 100e3], - ["asymmetric_inflection_1st_ord", 10e6, 200e6, -20e6, -50e6, -180e6, -20e6, -10e6], - ["inflection_3st_ord", 10e6, 200e6, -80e6, 80e6, -180e6, -200e6, 100e3], - ) - @unpack - def test_scanning_frequency(self, po1, po2, po3, no1, no2, no3, ferr): - """Test scanning frequency with experiment service. - - This is a sort of round-trip test. - We generate amplitudes from frequencies through the experiment class. - These amplitudes are converted into frequencies again with the same coefficients. - The two sets of frequencies must be consistent. - """ - service = self.create_service_helper(po1, po2, po3, no1, no2, no3, ferr, 0, "fake_hanoi") - - exp = StarkP1Spectroscopy((0,), backend=FakeHanoiV2()) - - ref_freqs = np.linspace(-70e6, 70e6, 31) - exp.set_experiment_options( - xvals=ref_freqs, - xval_type="frequency", - service=service, - ) - - amplitudes = exp.parameters() - - # Compute frequency from parameter values with the same coefficient - analysis = StarkP1SpectAnalysis() - coeffs = analysis.retrieve_coefficients_from_service(service, 0, "fake_hanoi") - frequencies = analysis._convert_axis(amplitudes, coeffs) - - np.testing.assert_array_almost_equal(frequencies, ref_freqs) - def test_scanning_frequency_with_coeffs(self): - """Test scanning frequency with manually provided Stark coefficients. - - This is just a difference of API from the test_scanning_frequency. - Data driven test is omitted here. - """ - po1, po2, po3, no1, no2, no3, ferr = 5e6, 200e6, -50e6, 5e6, -180e6, -40e6, 100e3 + """Test scanning frequency with manually provided Stark coefficients.""" + coeffs = util.StarkCoefficients( + pos_coef_o1=5e6, + pos_coef_o2=200e6, + pos_coef_o3=-50e6, + neg_coef_o1=5e6, + neg_coef_o2=-180e6, + neg_coef_o3=-40e6, + offset=100e3, + ) exp = StarkP1Spectroscopy((0,), backend=FakeHanoiV2()) ref_amps = np.array([-0.50, -0.25, 0.0, 0.25, 0.50], dtype=float) - test_freqs = np.where( - ref_amps > 0, - po1 * ref_amps + po2 * ref_amps**2 + po3 * ref_amps**3 + ferr, - no1 * ref_amps + no2 * ref_amps**2 + no3 * ref_amps**3 + ferr, - ) + test_freqs = util.convert_amp_to_freq(ref_amps, coeffs) exp.set_experiment_options( xvals=test_freqs, xval_type="frequency", - stark_coefficients={ - "stark_pos_coef_o1": po1, - "stark_pos_coef_o2": po2, - "stark_pos_coef_o3": po3, - "stark_neg_coef_o1": no1, - "stark_neg_coef_o2": no2, - "stark_neg_coef_o3": no3, - "stark_ferr": ferr, - }, + stark_coefficients=coeffs, ) params = exp.parameters() np.testing.assert_array_almost_equal(params, ref_amps) def test_scaning_frequency_around_zero(self): """Test scanning frequency around zero.""" + coeffs = util.StarkCoefficients( + pos_coef_o1=5e6, + pos_coef_o2=100e6, + pos_coef_o3=10e6, + neg_coef_o1=-5e6, + neg_coef_o2=-100e6, + neg_coef_o3=-10e6, + offset=500e3, + ) exp = StarkP1Spectroscopy((0,), backend=FakeHanoiV2()) exp.set_experiment_options( xvals=[0, 500e3], xval_type="frequency", - stark_coefficients={ - "stark_pos_coef_o1": 5e6, - "stark_pos_coef_o2": 100e6, - "stark_pos_coef_o3": 10e6, - "stark_neg_coef_o1": -5e6, - "stark_neg_coef_o2": -100e6, - "stark_neg_coef_o3": -10e6, - "stark_ferr": 500e3, - }, + stark_coefficients=coeffs, ) params = exp.parameters() # Frequency offset is 500 kHz and we need negative shift to tune frequency at zero. @@ -278,50 +199,6 @@ def test_circuits(self): self.assertEqual(circs[0], qc1) self.assertEqual(circs[1], qc2) - def test_retrieve_coefficients(self): - """Test retrieving Stark coefficients from the experiment service.""" - po1, po2, po3, no1, no2, no3, ferr = 5e6, 200e6, -50e6, 5e6, -180e6, -40e6, 100e3 - service = self.create_service_helper(po1, po2, po3, no1, no2, no3, ferr, 0, "fake_hanoi") - - retrieved_coeffs = StarkP1SpectAnalysis.retrieve_coefficients_from_service( - service=service, - qubit=0, - backend="fake_hanoi", - ) - self.assertDictEqual( - retrieved_coeffs, - { - "stark_pos_coef_o1": po1, - "stark_pos_coef_o2": po2, - "stark_pos_coef_o3": po3, - "stark_neg_coef_o1": no1, - "stark_neg_coef_o2": no2, - "stark_neg_coef_o3": no3, - "stark_ferr": ferr, - }, - ) - - @named_data( - ["ordinary", 5e6, 200e6, -50e6, 5e6, -180e6, -40e6, 100e3], - ["asymmetric_inflection_1st_ord", 10e6, 200e6, -20e6, -50e6, -180e6, -20e6, -10e6], - ["inflection_3st_ord", 10e6, 200e6, -80e6, 80e6, -180e6, -200e6, 100e3], - ) - @unpack - def test_estimating_minmax_frequency(self, po1, po2, po3, no1, no2, no3, ferr): - """Test computing the minimum and maximum frequency within the amplitude budget.""" - service = self.create_service_helper(po1, po2, po3, no1, no2, no3, ferr, 0, "fake_hanoi") - analysis = StarkP1SpectAnalysis() - - coeffs = analysis.retrieve_coefficients_from_service(service, 0, "fake_hanoi") - # pylint: disable=unbalanced-tuple-unpacking - minf, maxf = analysis.estimate_minmax_frequencies(coeffs, (-0.9, 0.9)) - - amps = np.linspace(-0.9, 0.9, 101) - freqs = analysis._convert_axis(amps, coeffs) - - self.assertAlmostEqual(minf, min(freqs), delta=1e6) - self.assertAlmostEqual(maxf, max(freqs), delta=1e6) - def test_running_analysis_without_service(self): """Test running analysis without setting service to the experiment data. @@ -349,15 +226,42 @@ def test_running_analysis_with_service(self, po1, po2, po3, no1, no2, no3, ferr) This must convert x-axis into frequencies with the Stark coefficients. """ - service = self.create_service_helper(po1, po2, po3, no1, no2, no3, ferr, 0, "fake_hanoi") + mock_experiment_id = "6453f3d1-04ef-4e3b-82c6-1a92e3e066eb" + mock_result_id = "d067ae34-96db-4e8e-adc8-030305d3d404" + mock_backend = FakeHanoiV2().name + + coeffs = util.StarkCoefficients( + pos_coef_o1=po1, + pos_coef_o2=po2, + pos_coef_o3=po3, + neg_coef_o1=no1, + neg_coef_o2=no2, + neg_coef_o3=no3, + offset=ferr, + ) + + service = FakeService() + service.create_experiment( + experiment_type="StarkRamseyXYAmpScan", + backend_name=mock_backend, + experiment_id=mock_experiment_id, + ) + service.create_analysis_result( + experiment_id=mock_experiment_id, + result_data={"value": coeffs}, + result_type="stark_coefficients", + device_components=["Q0"], + tags=[], + quality="Good", + verified=False, + result_id=mock_result_id, + ) + analysis = StarkP1SpectAnalysisReturnXvals() xvals = np.linspace(-1, 1, 11) - ref_xvals = np.where( - xvals > 0, - po1 * xvals + po2 * xvals**2 + po3 * xvals**3 + ferr, - no1 * xvals + no2 * xvals**2 + no3 * xvals**3 + ferr, - ) + ref_fvals = util.convert_amp_to_freq(xvals, coeffs) + exp_data = ExperimentData( service=service, backend=FakeHanoiV2(), @@ -365,9 +269,9 @@ def test_running_analysis_with_service(self, po1, po2, po3, no1, no2, no3, ferr) exp_data.metadata.update({"physical_qubits": [0]}) for x in xvals: exp_data.add_data({"counts": {"0": 1000, "1": 0}, "metadata": {"xval": x}}) - analysis.run(exp_data, replace_results=True) - test_xvals = exp_data.analysis_results("xvals").value - np.testing.assert_array_almost_equal(test_xvals, ref_xvals) + analysis.run(exp_data, replace_results=True).block_for_results() + test_fvals = exp_data.analysis_results("xvals").value + np.testing.assert_array_almost_equal(test_fvals, ref_fvals) def test_running_analysis_with_user_provided_coeffs(self): """Test running analysis by manually providing Stark coefficients. @@ -376,30 +280,25 @@ def test_running_analysis_with_user_provided_coeffs(self): This is just a difference of API from the test_running_analysis_with_service. Data driven test is omitted here. """ - po1, po2, po3, no1, no2, no3, ferr = 5e6, 200e6, -50e6, 5e6, -180e6, -40e6, 100e3 + coeffs = util.StarkCoefficients( + pos_coef_o1=5e6, + pos_coef_o2=200e6, + pos_coef_o3=-50e6, + neg_coef_o1=5e6, + neg_coef_o2=-180e6, + neg_coef_o3=-40e6, + offset=100e3, + ) analysis = StarkP1SpectAnalysisReturnXvals() - analysis.set_options( - stark_coefficients={ - "stark_pos_coef_o1": po1, - "stark_pos_coef_o2": po2, - "stark_pos_coef_o3": po3, - "stark_neg_coef_o1": no1, - "stark_neg_coef_o2": no2, - "stark_neg_coef_o3": no3, - "stark_ferr": ferr, - } - ) + analysis.set_options(stark_coefficients=coeffs) xvals = np.linspace(-1, 1, 11) - ref_xvals = np.where( - xvals > 0, - po1 * xvals + po2 * xvals**2 + po3 * xvals**3 + ferr, - no1 * xvals + no2 * xvals**2 + no3 * xvals**3 + ferr, - ) + ref_fvals = util.convert_amp_to_freq(xvals, coeffs) + exp_data = ExperimentData() for x in xvals: exp_data.add_data({"counts": {"0": 1000, "1": 0}, "metadata": {"xval": x}}) - analysis.run(exp_data, replace_results=True) - test_xvals = exp_data.analysis_results("xvals").value - np.testing.assert_array_almost_equal(test_xvals, ref_xvals) + analysis.run(exp_data, replace_results=True).block_for_results() + test_fvals = exp_data.analysis_results("xvals").value + np.testing.assert_array_almost_equal(test_fvals, ref_fvals) diff --git a/test/library/characterization/test_stark_ramsey_xy.py b/test/library/driven_freq_tuning/test_stark_ramsey_xy.py similarity index 84% rename from test/library/characterization/test_stark_ramsey_xy.py rename to test/library/driven_freq_tuning/test_stark_ramsey_xy.py index 795fde5297..33188dcf92 100644 --- a/test/library/characterization/test_stark_ramsey_xy.py +++ b/test/library/driven_freq_tuning/test_stark_ramsey_xy.py @@ -21,7 +21,8 @@ from qiskit.providers.fake_provider import FakeHanoiV2 from qiskit_experiments.library import StarkRamseyXY, StarkRamseyXYAmpScan -from qiskit_experiments.library.characterization.analysis import StarkRamseyXYAmpScanAnalysis +from qiskit_experiments.library.driven_freq_tuning.analyses import StarkRamseyXYAmpScanAnalysis +from qiskit_experiments.library.driven_freq_tuning import coefficient_utils as util from qiskit_experiments.framework import ExperimentData @@ -242,24 +243,33 @@ def test_ramsey_fast_analysis(self, c1p, c2p, c3p, c1n, c2n, c3n, ferr): exp_data = ExperimentData() exp_data.metadata.update({"stark_length": 50e-9}) + ref_coeffs = util.StarkCoefficients( + pos_coef_o1=c1p, + pos_coef_o2=c2p, + pos_coef_o3=c3p, + neg_coef_o1=c1n, + neg_coef_o2=c2n, + neg_coef_o3=c3n, + offset=ferr, + ) + yvals = util.convert_amp_to_freq(xvals, ref_coeffs) + # Generate fake data based on fit model. - for x in xvals: + for x, y in zip(xvals, yvals): if x >= 0.0: - fs = c1p * x + c2p * x**2 + c3p * x**3 + ferr direction = "pos" else: - fs = c1n * x + c2n * x**2 + c3n * x**3 + ferr direction = "neg" # Add some sampling error - ramx_count = rng.binomial(shots, amp * np.cos(const * fs) + off) + ramx_count = rng.binomial(shots, amp * np.cos(const * y) + off) exp_data.add_data( { "counts": {"0": shots - ramx_count, "1": ramx_count}, "metadata": {"xval": x, "series": "X", "direction": direction}, } ) - ramy_count = rng.binomial(shots, amp * np.sin(const * fs) + off) + ramy_count = rng.binomial(shots, amp * np.sin(const * y) + off) exp_data.add_data( { "counts": {"0": shots - ramy_count, "1": ramy_count}, @@ -271,38 +281,18 @@ def test_ramsey_fast_analysis(self, c1p, c2p, c3p, c1n, c2n, c3n, ferr): analysis.run(exp_data, replace_results=True) self.assertExperimentDone(exp_data) - # Check the fitted parameter can approximate the same polynominal - x_pos = np.linspace(0, 1, 51) - x_neg = np.linspace(-1, 0, 51) - ref_yvals_pos = c1p * x_pos + c2p * x_pos**2 + c3p * x_pos**3 + ferr - ref_yvals_neg = c1n * x_neg + c2n * x_neg**2 + c3n * x_neg**3 + ferr - - # Note that these parameter values are not necessary the same with input values - # as long as they can approximate the original phase polynominal. - c1p_est = exp_data.analysis_results("stark_pos_coef_o1").value.n - c2p_est = exp_data.analysis_results("stark_pos_coef_o2").value.n - c3p_est = exp_data.analysis_results("stark_pos_coef_o3").value.n - c1n_est = exp_data.analysis_results("stark_neg_coef_o1").value.n - c2n_est = exp_data.analysis_results("stark_neg_coef_o2").value.n - c3n_est = exp_data.analysis_results("stark_neg_coef_o3").value.n - ferr_est = exp_data.analysis_results("stark_ferr").value.n - - test_yvals_pos = c1p_est * x_pos + c2p_est * x_pos**2 + c3p_est * x_pos**3 + ferr_est - test_yvals_neg = c1n_est * x_neg + c2n_est * x_neg**2 + c3n_est * x_neg**3 + ferr_est - - # Check similality of reconstructed polynominals - # Curves must be agree within the torelance of 1.5 * 1 MHz. - np.testing.assert_array_almost_equal( - test_yvals_pos, - ref_yvals_pos, - decimal=-6, - err_msg="Reconstructed phase polynominal on positive frequency shift side " - "doesn't match with the original curve.", - ) + # Check the fitted parameter can approximate the same polynominal. + # Note that coefficient values don't need to exactly match as long as + # frequency shift is predictable. + # Since the fit model is just an empirical polynomial, + # comparing coefficients don't physically sound. + # Curves must be agreed within the tolerance of 1.5 * 1 MHz. + fit_coeffs = exp_data.analysis_results("stark_coefficients").value + fit_yvals = util.convert_amp_to_freq(xvals, fit_coeffs) + np.testing.assert_array_almost_equal( - test_yvals_neg, - ref_yvals_neg, + yvals, + fit_yvals, decimal=-6, - err_msg="Reconstructed phase polynominal on negative frequency shift side " - "doesn't match with the original curve.", + err_msg="Reconstructed phase polynominal doesn't match with the actual phase shift.", ) diff --git a/test/library/driven_freq_tuning/test_utils.py b/test/library/driven_freq_tuning/test_utils.py new file mode 100644 index 0000000000..11a30460c3 --- /dev/null +++ b/test/library/driven_freq_tuning/test_utils.py @@ -0,0 +1,161 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2024. +# +# 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. + +"""Test for Stark coefficients utility.""" + +from test.base import QiskitExperimentsTestCase + +from ddt import ddt, named_data, data, unpack +import numpy as np + +from qiskit_experiments.library.driven_freq_tuning import coefficient_utils as util +from qiskit_experiments.test import FakeService + + +@ddt +class TestStarkUtil(QiskitExperimentsTestCase): + """Test cases for Stark coefficient utilities.""" + + def test_coefficients(self): + """Test getting group of coefficients.""" + coeffs = util.StarkCoefficients( + pos_coef_o1=1e6, + pos_coef_o2=2e6, + pos_coef_o3=3e6, + neg_coef_o1=-1e6, + neg_coef_o2=-2e6, + neg_coef_o3=-3e6, + offset=0, + ) + self.assertListEqual(coeffs.positive_coeffs(), [3e6, 2e6, 1e6]) + self.assertListEqual(coeffs.negative_coeffs(), [-3e6, -2e6, -1e6]) + + @named_data( + ["ordinary", 5e6, 200e6, -50e6, 5e6, -180e6, -40e6, 100e3], + ["asymmetric_inflection_1st_ord", 10e6, 200e6, -20e6, -50e6, -180e6, -20e6, -10e6], + ["inflection_3st_ord", 10e6, 200e6, -80e6, 80e6, -180e6, -200e6, 100e3], + ) + @unpack + def test_roundtrip_convert_freq_amp( + self, + pos_o1: float, + pos_o2: float, + pos_o3: float, + neg_o1: float, + neg_o2: float, + neg_o3: float, + offset: float, + ): + """Test round-trip conversion between frequency shift and Stark amplitude.""" + coeffs = util.StarkCoefficients( + pos_coef_o1=pos_o1, + pos_coef_o2=pos_o2, + pos_coef_o3=pos_o3, + neg_coef_o1=neg_o1, + neg_coef_o2=neg_o2, + neg_coef_o3=neg_o3, + offset=offset, + ) + target_freqs = np.linspace(-70e6, 70e6, 11) + test_amps = util.convert_freq_to_amp(target_freqs, coeffs) + test_freqs = util.convert_amp_to_freq(test_amps, coeffs) + + np.testing.assert_array_almost_equal(test_freqs, target_freqs, decimal=2) + + @data( + [-0.5, 0.5], + [-0.9, 0.9], + [0.25, 1.0], + ) + @unpack + def test_calculate_min_max_shift(self, min_amp, max_amp): + """Test estimating maximum frequency shift within given Stark amplitude budget.""" + + # These coefficients induce inflection points around ±0.75, for testing + coeffs = util.StarkCoefficients( + pos_coef_o1=10e6, + pos_coef_o2=100e6, + pos_coef_o3=-90e6, + neg_coef_o1=80e6, + neg_coef_o2=-180e6, + neg_coef_o3=-200e6, + offset=100e3, + ) + # This numerical solution is correct up to amp resolution of 0.001 + nop = int((max_amp - min_amp) / 0.001) + amps = np.linspace(min_amp, max_amp, nop) + freqs = util.convert_amp_to_freq(amps, coeffs) + + # This finds strict solution, unless it has a bug + min_freq, max_freq = util.find_min_max_frequency( + min_amp=min_amp, + max_amp=max_amp, + coeffs=coeffs, + ) + + # Allow 1kHz tolerance because ref is approximate value + self.assertAlmostEqual(min_freq, np.min(freqs), delta=1e3) + self.assertAlmostEqual(max_freq, np.max(freqs), delta=1e3) + + def test_get_coeffs_from_service(self): + """Test retrieve the saved Stark coefficients from the experiment service.""" + mock_experiment_id = "6453f3d1-04ef-4e3b-82c6-1a92e3e066eb" + mock_result_id = "d067ae34-96db-4e8e-adc8-030305d3d404" + mock_backend = "mock_backend" + + ref_coeffs = util.StarkCoefficients( + pos_coef_o1=1e6, + pos_coef_o2=2e6, + pos_coef_o3=3e6, + neg_coef_o1=-1e6, + neg_coef_o2=-2e6, + neg_coef_o3=-3e6, + offset=0, + ) + + service = FakeService() + service.create_experiment( + experiment_type="StarkRamseyXYAmpScan", + backend_name=mock_backend, + experiment_id=mock_experiment_id, + ) + service.create_analysis_result( + experiment_id=mock_experiment_id, + result_data={"value": ref_coeffs}, + result_type="stark_coefficients", + device_components=["Q0"], + tags=[], + quality="Good", + verified=False, + result_id=mock_result_id, + ) + + retrieved = util.retrieve_coefficients_from_service( + service=service, + backend_name=mock_backend, + qubit=0, + ) + + self.assertEqual(retrieved, ref_coeffs) + + def test_get_coeffs_no_data(self): + """Test raises when Stark coefficients don't exist in the result database.""" + mock_backend = "mock_backend" + + service = FakeService() + + with self.assertRaises(RuntimeError): + util.retrieve_coefficients_from_service( + service=service, + backend_name=mock_backend, + qubit=0, + )