Spectral Quantum Tomography

TL;DR

Spectral quantum tomography is a resource-efficient method for extracting eigenvalues of noisy quantum gates, providing detailed diagnostics beyond traditional fidelity measures, and is applicable to real superconducting qubit systems.

Contribution

We introduce spectral quantum tomography, a novel SPAM-resistant technique for eigenvalue extraction of quantum gates using minimal resources.

Findings

Successfully applied to superconducting qubits on cloud platforms.

Revealed effects of cross-talk, leakage, and non-Markovian errors.

Provides detailed gate diagnostics beyond fidelity measures.

Abstract

We introduce spectral quantum tomography, a simple method to extract the eigenvalues of a noisy few-qubit gate, represented by a trace-preserving superoperator, in a SPAM-resistant fashion, using low resources in terms of gate sequence length. The eigenvalues provide detailed gate information, supplementary to known gate-quality measures such as the gate fidelity, and can be used as a gate diagnostic tool. We apply our method to one- and two-qubit gates on two different superconducting systems available in the cloud, namely the QuTech Quantum Infinity and the IBM Quantum Experience. We discuss how cross-talk, leakage and non-Markovian errors affect the eigenvalue data.

Equations

$$S_{\mu\nu}={\rm Tr}\big{[}P_{\mu}{\cal S}(P_{\nu})\big{]},\;\;\mu,\nu=0,\ldots,N.$$

$${\cal S}\leftrightarrow S=\left(\begin{array}[]{cc}1&0\\ {\bf s}&T^{\cal S}\\ \end{array}\right),$$

$${\cal F}_{\rm ent}(\mathcal{S},U)={\rm Tr}\big{[}I\otimes U\ket{\Psi}\!\bra{\Psi}I\otimes U^{\dagger}(I\otimes{\cal S})(\ket{\Psi}\!\bra{\Psi})\big{]},$$

$${\cal F}_{\rm ent}(\mathcal{S},U)=\frac{1}{d^{2}}\sum_{\mu}{\rm Tr}\big{[}UP_{\mu}U^{\dagger}{\cal S}(P_{\mu})\big{]}=\frac{1}{d^{2}}\big{(}1+\operatorname{Tr}\big{[}T^{U^{\dagger}}T^{\cal S}\big{]}\big{)}.$$

$${\mathcal{F}}_{\mathrm{ent}}(\mathcal{S},U)\le \frac{1}{d^2}\left[1+(d^2-1)\left(\sqrt{1-\frac{{\sum }_j|{\lambda }_j{|}^2}{d^2-1}}+{\xi }_{\max }\right)\right],$$

$$u({\cal S})=\frac{d}{d-1}\int d\phi\;{\rm Tr}\big{[}[{\cal S}^{\prime}(\ket{\phi}\!\bra{\phi})]^{\dagger}{\cal S}^{\prime}(\ket{\phi}\!\bra{\phi})\big{]},$$

$$u({\cal S})=\frac{1}{d^{2}-1}{\rm Tr}\big{[}{T^{\mathcal{S}}}^{\dagger}{T^{\mathcal{S}}}\big{]}=\frac{1}{d^{2}-1}\sum_{i}\sigma_{i}(T^{\mathcal{S}})^{2},$$

$$u(\mathcal{S})\ge \frac{1+{\sum }_{i=1}^{d^2-1}{\left|{\lambda }_i\right|}^2-d}{d(d-1)},$$

$${\left|\left|\mathbf{s}\right|\right|}^2\le 1-{\left|{\lambda }_1\right|}^2-{\left|{\lambda }_2\right|}^2-{\left|{\lambda }_3\right|}^2+2{\lambda }_1{\lambda }_2{\lambda }_3.$$

$$\dot{\rho }=\mathcal{L}(\rho ).$$

$$\mathcal{L}(\rho)=-i[H,\rho]+\frac{1}{T_{1}}\Bigl{(}\sigma_{-}\rho\sigma_{+}-\frac{1}{2}\{\sigma_{+}\sigma_{-},\rho\}\Bigr{)}+\frac{1}{2T_{\phi}}(Z\rho Z-\rho).$$

$$L^{\mathcal{L}}=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& -{\Gamma }_2& h_z& h_y\\ 0& -h_z& -{\Gamma }_2& h_x\\ {\Gamma }_1& -h_y& -h_x& -{\Gamma }_1\end{array}\right).$$

$${\Omega }_j\in \{-{\Gamma }_2,-{\Gamma }_2,-{\Gamma }_1\},$$

$${\lambda }_j\in \{e^{-{\Gamma }_2\tau },e^{-{\Gamma }_2\tau },e^{-{\Gamma }_1\tau }\},$$

$${\Omega }_j\in \{-{\Gamma }_2+i{h}_z,-{\Gamma }_2-i{h}_z,-{\Gamma }_1\},$$

$${\lambda }_j\in \{e^{-{\Gamma }_2\tau }{e}^{i{h}_z\tau },e^{-{\Gamma }_2\tau }{e}^{-i{h}_z\tau },e^{-{\Gamma }_1\tau }\},$$

$$\Omega_{2,3}(h_{x})=-\frac{1}{2}\Bigl{(}\Gamma_{1}+\Gamma_{2}\pm\sqrt{(\Gamma_{1}-\Gamma_{2})^{2}-4h_{x}^{2}}\Bigr{)}.$$

$$g(k)=\sum_{\mu=1}^{N}{\rm Tr}\big{[}P_{\mu}{\cal N}_{\rm meas}\circ{\cal S}^{k}\circ{\cal N}_{\rm prep}(P_{\mu})\big{]}.$$

$$g^{\mbox{\tiny NO SPAM}}(k)=\sum_{\mu=1}^{N}(T^{k})_{\mu\mu}={\rm Tr}\big{[}T^{k}\big{]}=\sum_{j=1}^{N}\lambda_{j}^{k},$$

$$g(k)={\rm Tr}\big{[}T_{\rm meas}T^{k}T_{\rm prep}\big{]}={\rm Tr}\big{[}A_{\rm SPAM}D^{k}\big{]}=\sum_{j=1}^{N}A_{j}\lambda_{j}^{k},$$

$$Y$$

$$\displaystyle=\sum_{j=1}^{N}A_{j}\left(\begin{array}[]{cccc}1&\lambda_{j}&\ldots&\lambda_{j}^{L}\\ \lambda_{j}&\lambda_{j}^{2}&\ldots&\lambda_{j}^{L+1}\\ \lambda_{j}^{2}&\vdots&&\vdots\\ \vdots&&&\vdots\\ \lambda_{j}^{K-L}&\ldots&\ldots&\lambda_{j}^{K}\end{array}\right).$$

$$\mathfrak{T}\left(\begin{array}[]{c}1\\ \lambda_{j}\\ \vdots\\ \lambda_{j}^{M}\end{array}\right)=\lambda_{j}\left(\begin{array}[]{c}1\\ \lambda_{j}\\ \vdots\\ \lambda_{j}^{M}\end{array}\right).$$

$${ϵ}_N^{\mathrm{rms}}={\left(\frac{1}{K+1}\sum _{k=0}^K|g(k)-g_N^{\mathrm{est}}(k){|}^2\right)}^{1/2}.$$

$$T^{\rm GM}=\left(\begin{array}[]{cc}T_{\rm comp}&T_{\rm seep}\\ T_{\rm leak}&T_{\rm beyond}\end{array}\right),$$

$$g^{\text{mbox}NOSPAM,LEAK}(k)=\sum _{j=1}^8{\tilde{A}}_j{\lambda }_j^k,{\tilde{A}}_j=⟨{\sigma }_j|V^{-1}{\Pi }_{\mathrm{comp}}V|{\sigma }_j⟩.$$

$$\mathcal{S}_{k}(\rho)=\int_{-\infty}^{+\infty}d\phi\,\mathbb{P}_{\sigma}(\phi)\bigl{(}R_{x}(\phi)\,\mathrm{CZ}\bigr{)}^{k}\rho\bigl{(}\mathrm{CZ}\,R_{x}(\phi)^{\dagger}\bigr{)}^{k}.$$

$$P_e(t)=\frac{1}{2}+\frac{1}{2}\sum _{n=0}^{\infty }{p}_{\alpha }(n)\cos (\Omega t\sqrt{n+1}),$$

$$J_{i}=\left(\begin{array}[]{ccccc}\lambda_{i}&1&&\\ &\lambda_{i}&\ddots&\\ &&\ddots&1\\ &&&\lambda_{i}\end{array}\right),$$

$$S=\left(\begin{array}[]{cccc}1&0&0&0\\ 0&1-\epsilon x&-\epsilon y&0\\ 0&\epsilon y&1-\epsilon(x+2y)&0\\ 2\epsilon x&0&0&1-2\epsilon(x+y)\end{array}\right)+O(\epsilon^{2}).$$

$$g^{\mbox{\tiny NO SPAM}}(k)={\rm Tr}\big{[}T^{k}\big{]}={\rm Tr}\big{[}(D+E)^{k}\big{]}={\rm Tr}\big{[}D^{k}\big{]},$$

$$g(k)={\rm Tr}\big{[}W^{\dagger}T_{\rm prep}T_{\rm meas}W(D+E)^{k}\big{]},$$

Full text

I Introduction

A central challenge on the path towards large-scale quantum computing is the engineering of high-quality quantum gates. To achieve this goal, many methods which accurately and reliably characterize quantum gates have been developed. Some of these methods are scalable, meaning that they require an effort which scales polynomially in the number of qubits on which the gates act. Scalable protocols, such as randomized benchmarking knill+:RB ; MGE:RB ; magesan+:IRB ; WG:leak ; wallman+:unitarity ; DHW:URB ; erhard+:cycle ; OWE:RB necessarily give a partial characterization of the gate quality, for example an average gate fidelity. Other protocols such as robust tomography kimmel+:tomo or gate-set tomography BK:GST ; greenbaum:GST trade scalability for a more detailed characterization of the gate. A desirable feature of all the above protocols is that they are resistant to state-preparation and measurement (SPAM) errors. The price of using SPAM-resistant (scalable) methods is that these protocols have significant experimental complexity and/or require assumptions on the underlying hardware to properly interpret their results.

In this work we present spectral quantum tomography, a simple non-scalable method that extracts spectral information from noisy gates in a SPAM-resistant manner. To process the tomographic data and obtain the spectrum of the noisy gate, we rely on the matrix-pencil technique, a well-known classical signal processing method. This technique has been advocated in OWE:RB in the context of randomized benchmarking and has also been used in OTT:QPE for processing data in the algorithm of quantum phase estimation. It has also been used, under the phrase ‘linear systems identification’, in BL:system to predict the time evolution of quantum systems. While the matrix pencil technique leads to explicitly useful estimates of eigenvalues and their amplitudes, we note that the same underlying idea is used in the method of “delayed vectors” which has been proposed in WP:dim to assess the dimensionality of a quantum system from its dynamics. This “delayed vectors” approach has been applied to assess leakage in superconducting devices in SDK:leak-witness .

The spectral information of a noisy gate $\mathcal{S}$, which approximates some target unitary $U$, is given by the eigenvalues of the so-called Pauli transfer matrix representing ${\cal S}$. These eigenvalues, which are of the form $\lambda=\exp(-\gamma)\exp(i\phi)$, contain information about the quality of the implemented gate. Intuitively, the parameter $\gamma$ captures how much the noisy gate deviates from unitarity due to entanglement with an environment, while the angle $\phi$ can be compared to the rotation angles of the targeted gate $U$. Hence $\phi$ gives information about how much one over- or under-rotates. The spectrum of $\mathcal{S}$ can also be related to familiar gate-quality measures such as the average gate fidelity and the unitarity. Moreover, in the case of a noisy process modeled by a Lindblad equation, the spectrum can be easily related to the more familiar notions of relaxation and dephasing times.

The main advantage of spectral quantum tomography is its simplicity, requiring only the (repeated) application of a single noisy gate $\mathcal{S}$, as opposed to the application of a large set of gates as in randomized benchmarking, gate-set tomography and robust tomography. Naturally, simplicity and low-cost come with some drawback, namely the method does not give information about the eigenvectors of the noisy gate, such as the axis around which one is rotating. However, information about the eigenvectors is intrinsically hard to extract in a SPAM-resistant fashion since SPAM errors can lead to additional rotations RPZ:gauge . Another feature of spectral quantum tomography is that it can be used to extract signatures of non-Markovianity, namely the phenomenon where the noisy gate $\mathcal{S}$ depends on the context in which it is applied (e.g. time of application, whether any gates have been applied before it). As we show in this paper, our method can be used to detect various types of non-Markovian effects such as coherent revivals, parameter drifts, and Gaussian-distributed time-correlated noise. It is also possible to distinguish non-Markovian effects from qubit leakage. For these reasons we believe that spectral quantum tomography adds a useful new tool to the gate-characterization toolkit. The method could also have future applications in assessing the performance of logical gates in a manner which is free of logical state preparation and measurement errors, see the Discussion Section.

II Results

II.1 Eigenvalues of Trace-Preserving Completely-Positive (TPCP) maps

Take a unitary gate $U$ on a $d$-dimensional space with $U\ket{\psi_{j}}=e^{i\phi_{j}}\ket{\psi_{j}}$. The corresponding TPCP map ${\cal S}_{U}(\rho)=U\rho U^{\dagger}$ has one trace-full eigenvector, namely $I$ with eigenvalue 1, as well as $d^{2}-1$ traceless eigenvectors. In particular, there are $d^{2}-d$ traceless eigenvectors of the form $\ket{\psi_{j}}\bra{\psi_{l}}$ for $j\neq l$ with eigenvalues $\exp(i(\phi_{j}-\phi_{l}))$, and $d-1$ traceless eigenvectors of the form $\ket{\psi_{1}}\bra{\psi_{1}}-\ket{\psi_{j}}\bra{\psi_{j}}$ for $j=2,\ldots,d$ with eigenvalue $1$.

For general TPCP maps it is convenient to use the Pauli transfer matrix formalism. For an $n$-qubit system ($d=2^{n}$) consider the normalized set of Pauli matrices $P_{\mu}$ for $\mu=0,\ldots,N$ with $N+1=4^{n}=d^{2}$, where $P_{0}=I/\sqrt{2^{n}}$ and the normalization is chosen such that {\rm Tr}\big{[}P_{\mu}P_{\nu}\big{]}=\delta_{\mu\nu}. For a TPCP map ${\cal S}$ acting on $n$ qubits, the Pauli transfer matrix is then defined as

[TABLE]

The form of the Pauli transfer matrix $S$ is RSW:TCP

[TABLE]

where $T^{\cal S}$ is a real $N\times N$ matrix and ${\bf s}$ is a $N$-dimensional column vector. The 1 and 0’s in the top row of the Pauli transfer matrix are due to the fact that ${\cal S}$ is trace-preserving. For a unital ${\cal S}$ which obeys ${\cal S}(I)=I$, the vector ${\bf s}=0$.

A few properties are known of the eigenvalue-eigenvector pairs of $S$, i.e. the pairs $(\lambda,\vec{\bf v})$ with $S{\bf v}=\lambda{\bf v}$:

• •

The eigenvalues of $S$ are $1$ and the eigenvalues of $T^{\cal S}$ since the solutions of the equation ${\rm det}(S-\lambda I)=0$ are the solutions of the equation $(1-\lambda){\rm det}(T^{\cal S}-\lambda I)=0$.

• •

The eigenvalues of $S$, and thus the eigenvalues of $T^{\cal S}$, come in complex-conjugate pairs. This is true because $T^{\cal S}$ is a real matrix.

• •

The eigenvalues of $T^{\cal S}$ (or $S$ for that matter) have modulus less than 1, i.e. $|\lambda|\leq 1$ (see e.g. Proposition 6.1 in wolf ).

If $T^{\cal S}$ is diagonalizable as a matrix, it holds that $T^{\cal S}=VDV^{-1}$ where $D$ is a diagonal matrix and $V$ a similarity transformation. Generically, $T^{\cal S}$ will be diagonalizable, in which case there are $N$ eigenvalue-eigenvector pairs for $T$. A sufficient condition for diagonizability is, for example, that all the eigenvalues of $T^{\cal S}$ are distinct. In Appendix IV.1 we give examples and discuss what it means if $T^{\cal S}$ is not diagonalizable.

For some simple single-qubit channels we can explicitly compute the spectrum. For instance, for a single-qubit depolarizing channel with depolarizing probability $p$, the eigenvalues of the sub-matrix $T^{\cal S}$ of the Pauli transfer matrix are $\{1-p,1-p,1-p\}$. For a single qubit amplitude-damping channel with damping rate $p$ they are $\{\sqrt{1-p},\sqrt{1-p},1-p\}$ greenbaum:GST .

II.1.1 Relation to gate-quality measures

The eigenvalues of the Pauli transfer matrix of a noisy gate ${\cal S}$ can be related to several other known measures of gate quality such as the average gate fidelity ${\cal F}(\mathcal{S},U)$, the gate unitarity $u({\cal S})$ and, for a single qubit ($n=1$), the gate unitality.

The average gate fidelity is defined as ${\cal F}(\mathcal{S},U)=\int d\phi\bra{\phi}U^{\dagger}{\cal S}(\ket{\phi}\!\bra{\phi})U\ket{\phi}$. This fidelity relates directly to the entanglement fidelity ${\cal F}_{\rm ent}(\mathcal{S},U)$ via ${\cal F}=\frac{{\cal F}_{\rm ent}d+1}{d+1}$ HHH:singlet , where the entanglement fidelity is defined as

[TABLE]

where $\ket{\Psi}=\frac{1}{\sqrt{d}}\sum_{i=1}^{d}\ket{i,i}$ is a maximally entangled state. Using that $\ket{\Psi}\bra{\Psi}=\frac{1}{d}\sum_{\mu=0}^{N}P_{\mu}\!\otimes\!P_{\mu}$ and $UP_{\mu}U^{\dagger}=\sum_{\kappa}T_{\mu\kappa}^{U^{\dagger}}P_{\kappa}$ we can write

[TABLE]

Thus for the (entanglement) fidelity of a noisy gate ${\cal S}$ with respect to the identity channel $U=I$, one has ${\cal F}_{\rm ent}(\mathcal{S},I)=\frac{1}{d^{2}}(1+\sum_{i}\lambda_{i})$, implying a direct relation to the spectrum $\{\lambda_{i}\}$ of $T^{\cal S}$. A more interesting relation is how the eigenvalues of $T^{\cal S}$ bound the fidelity with respect to a targeted gate $U$. In Appendix IV.2 we prove that the entanglement fidelity can be upper bounded as

[TABLE]

where $\xi_{\rm max}=\frac{1}{d^{2}-1}|\sum_{j}\lambda^{\rm ideal}_{j}\lambda_{j}^{*}|$ with $\lambda^{\rm ideal}_{j}$ the eigenvalues of $T^{U}$ with $U$ the targeted unitary, ordered such that the sum $|\sum_{j}\lambda^{\rm ideal}_{j}\lambda_{j}^{*}|$ is maximal.

This upper bound is not particularly tight, but for the case of a single qubit we can make a much stronger numerical statement, see Section IV.2.

Another measure of gate quality, namely the unitarity or the coherence of a channel wallman+:unitarity on a $d$-dimensional system, is defined as

[TABLE]

where $\mathcal{S}^{\prime}(\rho)\coloneqq\mathcal{S}(\rho)-\operatorname{Tr}[\mathcal{S}(\rho)]I/\sqrt{d}$. A more convenient but equivalent definition is

[TABLE]

where $\{\sigma_{i}\}$ are the singular values of the matrix $T^{\mathcal{S}}$.

The unitarity captures how close the channel is to a unitary gate. A lower bound on the unitarity is given by Proposition 2 in RPZ:gauge :

[TABLE]

where $\{\lambda_{i}\}$ are the eigenvalues of $T^{\mathcal{S}}$. For a single qubit, an upper bound on the unitarity can also be given in terms of a non-convex optimization problem, see Section IV.2.

The unitality of a TPCP map is defined as $1-\left|\left|{\mathbf{s}}\right|\right|^{2}$ with ${\bf s}$ in Eq. (2). Specifically, for single-qubit channels one can derive the bound RPZ:gauge

[TABLE]

II.1.2 Relation to relaxation and dephasing times

We consider the eigenvalues of a superoperator induced by a simple Lindblad equation modeling relaxation and decoherence of a driven qubit, as an example. We have a Lindblad equation with time-independent Lindbladian $\mathcal{L}$:

[TABLE]

The formal solution of Eq. 8 is given by $\rho(t)=e^{t\mathcal{L}}(\rho(t=0))$, where $e^{t\mathcal{L}}$ is a TPCP map for every $t$. We are interested in the total evolution after a certain gate time $\tau$ and set $\mathcal{S}_{\tau}=e^{\tau\mathcal{L}}$. We assume a simple model in which a qubit evolves according to a Hamiltonian $H=(h_{x}X+h_{y}Y+h_{z}Z)/2$ and is subject to relaxation and pure dephasing processes, according to the Lindbladian:

[TABLE]

We define the relaxation respectively dephasing rates $\Gamma_{1}=1/T_{1}$ and $\Gamma_{2}=1/T_{2}=1/(2T_{1})+1/T_{\phi}$. The Pauli transfer matrix $L^{\mathcal{L}}$ of $\mathcal{L}$ then takes the form

[TABLE]

We will denote the eigenvalues of $L^{\mathcal{L}}$ by $\Omega_{j}$ for $j\in\{0,\ldots,3\}$ and the eigenvalues of $\mathcal{S}_{\tau}$ by $\lambda_{j}$ for $j\in\{0,\ldots,3\}$. As expected, $\Omega_{0}=0$ implying that $\lambda_{0}=e^{0}=1$ is an eigenvalue of $\mathcal{S}_{\tau}$. The other three eigenvalues of $L^{\mathcal{L}}$ can be found from the $3\times 3$ sub-matrix in the lower-right corner. Here we consider some simple cases.

Case 1: $h_{x}=h_{y}=h_{z}=0$. In this case, for $j=1,2,3$ the three eigenvalues of ${\mathcal{L}}$ and ${\cal S}_{\tau}$ are clearly

[TABLE]

thus relating directly to the relaxation and dephasing rates.

Case 2: $h_{x}=h_{y}=0$. In this case we have

[TABLE]

where we have separated the decaying part of the $\lambda_{j}$ (corresponding to the real part of the $\Omega_{j}$) and their phases (corresponding to the imaginary part). If we work in the rotating frame of the qubit, $h_{z}$ can be understood as an over-rotation along the $Z$-axis, which would appear in the spectrum as an extra phase imparted to two of the eigenvalues. Again we see that the decaying part of the eigenvalues directly relates to the relaxation and dephasing rates.

Case 3: $h_{y}=h_{z}=0$. This case shows that over-rotations can also modify the decay strength of the eigenvalues. We analyze the eigenvalues as a function of $h_{x}$. From $L^{\mathcal{L}}$ in Eq. (9) we see that $\Omega_{1}(h_{x})=-\Gamma_{2}$ for all $h_{x}$. For the other eigenvalues we have

[TABLE]

We see that if $|h_{x}|<|\Gamma_{1}-\Gamma_{2}|/2\equiv h_{x}^{\mathrm{cr}}$, only the moduli of $\lambda_{2}$ and $\lambda_{3}$ are affected as compared to Case 1, in other words, $\lambda_{2}$ and $\lambda_{3}$ only decay with no extra phases. On the contrary, the phases of these eigenvalues becomes non-zero when the driving is sufficiently strong: $|h_{x}|>h_{x}^{\mathrm{cr}}$. It implies that if we look at the dynamics induced by the Lindblad equation, real oscillations, not only decay, will be present as a function of $\tau$. Hence these two scenarios represent respectively the overdamped ($|h_{x}|<h_{x}^{\mathrm{cr}}$) and underdamped regime ($|h_{x}|>h_{x}^{\mathrm{cr}}$), similar to the dynamics of a vacuum-damped qubit-oscillator system, see e.g. book:HR . At $|h_{x}|=h_{x}^{\mathrm{cr}}$, the system is critically damped and $L^{\mathcal{L}}$ does not have $4$ linearly-independent eigenvectors, meaning that the Pauli transfer matrix of ${\cal S}_{\tau}$ is not diagonalizable. In this case the dynamics also has a linear dependence on $t$ besides the exponential decay with $t$, see the discussion in Appendix IV.1.

II.2 Spectral tomography

In this section we describe the spectral tomography method, which estimates the eigenvalues of ${\cal S}$, where ${\cal S}$ is a TPCP implementation of a targeted unitary gate.

We model state-preparation errors as a perfect preparation step followed by an unknown TPCP map ${\cal N}_{\rm prep}$. Similarly, measurement errors are modeled by a perfect measurement preceded by an unknown TPCP map ${\cal N}_{\rm meas}$. We assume that when we apply the targeted gate $k$ times, an accurate model of the resulting noisy dynamics is ${\cal S}^{k}$.

The spectral tomography method can be applied without this assumption but the interpretation of its results is more difficult, see Section II.4 for a discussion.

The method works by constructing the following signal function, for $k=0,1,\ldots,K$ for some fixed $K$:

[TABLE]

Gathering the data to estimate $g(k)$ requires (1) picking a traceless $n$-qubit Pauli $P_{\mu}$, (2) preparing an $n$-qubit input state in one of the $2^{n}$ basis states corresponding to this chosen Pauli, (3) applying the gate $k$ times and measuring in the same chosen Pauli basis, and (4) repeating (1-3) over different Pauli’s, basis states and experiments to get good statistics. As in standard process tomography Nie00 , one takes linear combinations of the estimated probabilities for the outcomes to construct an estimator of a Pauli operator on a Pauli input. This gives an estimate of $g(k)$ for a fixed $k$. Repeating this process for $k\in\{0,\ldots,K\}$ we reconstruct the entire signal function. In Section II.2.2 we discuss the cost of doing these experiments as compared to randomized benchmarking.

Let us now examine how $g(k)$ depends on the eigenvalues of the matrix $T$. When there are no SPAM errors, that is, ${\cal N}_{\rm meas}$ and ${\cal N}_{\rm prep}$ are identity channels, we have

[TABLE]

where $\{\lambda_{j}\}$ are the eigenvalues of $T$. The last step in this equality follows directly when $T$ is diagonalizable, but it can alternatively be proved using the so-called Schur triangular form of $T$ (we give this proof in Appendix IV.1).

When ${\cal N}_{\rm meas}$ and ${\cal N}_{\rm prep}$ are not identity channels, we have

[TABLE]

where $T_{\rm meas}$ and $T_{\rm prep}$ are respectively the $T$-submatrices of the Pauli transfer matrix of ${\cal N}_{\rm meas}$ and ${\cal N}_{\rm prep}$. Here we assume that $T=VDV^{-1}$ is diagonalizable and the matrix $A_{\rm SPAM}=V^{-1}T_{\rm prep}T_{\rm meas}V$ captures the SPAM errors. One may expect that $A_{\rm SPAM}$ is close to the identity matrix in the typical case of low SPAM errors, in particular one may expect that $A_{j}\neq 0$ for all $j$ so that all eigenvalues of $T$ are present in the signal $g(k)$.

In principle, one could take more tomographic data and consider a full matrix-valued signal c_{\mu\nu}(k)={\rm Tr}\big{[}P_{\mu}{\cal N}_{\rm meas}\circ{\cal S}^{k}\circ{\cal N}_{\rm prep}(P_{\nu})\big{]} instead of only Eq. (11). This requires doing many more experiments and there is no clear advantage in terms of the ability to determine the spectrum.

II.2.1 Signal analysis or matrix-pencil method for extracting eigenvalues

In this section we review the classical signal-processing method which reconstructs, from the (noisy) signal $g(k)=\sum_{j=1}^{N}A_{j}\lambda_{j}^{k}$ for $k=0,\ldots,K$, an estimate for the eigenvalues $\lambda_{j}$ and the amplitudes $A_{j}$. Note that we have $g(k)\in\mathbb{R}$ due to Eq. (11). Not surprisingly, this signal-processing method has been employed and reinvented in a variety of scientific fields. We implement the so-called ESPRIT analysis described in SP:pencil , but see also pt:prony . In the context of spectral tomography we know that the signal $g(k)$ will in principle contain $N$ eigenvalues (which are possibly degenerate). However, we can vary the number of eigenvalues we use to fit the signal to see whether a different choice than $N$ gives a significantly better fit. This is relevant in particular when the implemented gate contains leakage or non-Markovian dynamics, see Section II.4.

We require at least $K\geq 2N-2$ in order to determine the eigenvalues accurately. This implies that for a single-qubit gate with $N=3$ we need at least $K=4$ and for a two-qubit gate with $N=15$ we need at least $K=28$. However, the signal $g(k)$ has sampling noise due to a bounded $N_{\rm samples}$ and in practice it is good to choose $K$ larger than strictly necessary to make the reconstruction more robust against noise. We study the effect of varying $K$ in Fig. 1 (left panel).

The method goes as follows and relies on picking a so-called pencil parameter $L$.

Let us assume for now that each $g(k)$ is learned without sampling noise. One constructs a $(K-L+1)\times(L+1)$-dimensional data matrix $Y$ as

[TABLE]

Note that ${\rm rank}(Y)\leq N$ since $Y$ is a sum of at most $N$ rank-1 matrices when there are $N$ eigenvalues. Consider two submatrices of $Y$: the matrix $G_{0}$ is obtained from $Y$ by deleting the last column of $Y$, while the matrix $G_{1}$ is obtained by deleting the first column of $Y$. When $L=\frac{K}{2}$, the matrices $G_{0}$ and $G_{1}$ are square matrices of dimension $M=\frac{K}{2}+1$. For this choice of $L$, the smallest value of $K$ so that $M=N$ is $2N-2$. We seek a time-shift matrix $\mathfrak{T}$ such that $\mathfrak{T}G_{0}=G_{1}$. When $M\geq N$, there certainly exists a matrix $\mathfrak{T}$ such that for all $j\in\{1,\ldots,N\}$:

[TABLE]

Furthermore, if $G_{0}^{-1}$ exists, which is the case when ${\rm rank}(G_{0})=M$, this matrix $\mathfrak{T}$ will be uniquely given as $G_{1}G_{0}^{-1}$. Hence, in this case there is a unique matrix $\mathfrak{T}$, obtained by constructing $G_{1}G_{0}^{-1}$ from the data, which is guaranteed to have $\{\lambda_{j}\}$ as eigenvalues. When the pencil parameter $L>\frac{K}{2}$, one needs to ensure that there are at least $N$ rows of the matrix $Y$: if not, $\mathfrak{T}$ would be of dimension less than $N$, not giving $N$ eigenvalues. This implies $K\geq N+L-1$.

The general method for a non-square $Y$ which includes an additional sampling-noise reduction step then goes as follows. The choice for $N$ in the procedure can be varied from its minimal value equal to $d^{2}-1$ to a larger value, depending on a goodness-of-fit.

Construct a singular-value decomposition of the matrix $Y$, i.e. $Y=R_{1}\Sigma R_{2}^{T}$ and replace the diagonal matrix $\Sigma$ by a diagonal matrix $\Sigma_{\rm clean}$ with only the largest $N$ singular values. Let $Y_{\rm clean}=R_{1}\Sigma_{\rm clean}R_{2}^{T}$. This step is to reduce sampling noise. 2. 2.

Take the submatrices $G_{0}$ and $G_{1}$ of $Y_{\rm clean}$. 3. 3.

Compute $\mathfrak{T}=G_{1}G_{0}^{+}$ where $G_{0}^{+}$ is the Moore-Penrose pseudo-inverse of the matrix $G_{0}$ so that $\mathfrak{T}$ is a matrix with at most $N$ non-zero eigenvalues. 4. 4.

Compute the eigenvalues of $\mathfrak{T}$: these will be the estimates $\lambda_{j}^{\rm est}$ of $\lambda_{j}$ for all $j\in\{1,\ldots N\}$. Formally, the linear matrix pencil is $G_{0}-\lambda G_{1}$ and the eigenvalues of this matrix pencil, i.e. the values where ${\rm det}(G_{0}-\lambda G_{1})=0$, are the $\lambda_{j}^{\rm est}$.

We have first applied this method on the signal $g(k)$ of a randomly chosen single-qubit channel: by varying $K$ and $L$ we want to understand the role of the matrix-pencil parameter $L$ and the choice for a larger $K$. The results are shown in Fig. 1 (left panel). Note that the chosen $K$’s are quite far above the bound $K\geq N+L-1$ to effectively suppress sampling noise. For each $K$ there is a flat region in $L$ where $\Delta^{2}$ is roughly constant. In the remainder we will choose $L=K/2$, putting ourselves in the middle of this region. Fig. 1 (right panel) shows how increasing $N_{\rm samples}$ lowers the total variance of the estimated eigenvalues.

An additional processing step is the determination of the (complex) amplitudes $\{A_{j}\}$. Viewing $g(k)$ as a set of $K+1$ inner products between the vector $(A_{1},\ldots,A_{N})$ and the linearly-independent vectors $(\lambda_{1}^{k},\lambda_{2}^{k},\ldots,\lambda_{N}^{k})$, it is clear that, given perfect knowledge of $g(k)$, the $\{A_{j}\}$ are uniquely determined when $K+1\geq N$. Since $g(k)$ is known with sampling noise, the $\{A_{j}\}$ can be found by solving the least-squares minimization problem $\min_{A_{j}}\sum_{k}|g(k)-\sum_{j}A_{j}(\lambda_{j}^{\rm est})^{k}|^{2}$. The optimal values in this minimization $A_{j}^{\rm est}$ and $\lambda_{j}^{\rm est}$ together form the reconstructed signal $g^{\rm est}(k)$ and the error is given by

[TABLE]

II.2.2 Resources

It is interesting to consider the amount of experiments that must be done to perform spectral quantum tomography. One must estimate the function $g(k)$ defined in Eq. 13. This reconstruction process requires running $2^{n}\times N\times(K+1)$ different experiments and repeating each experiment $N_{\rm samples}$ times. For a single-qubit gate we need $6(K+1)$ experiments, while for a two-qubit gate we need $60(K+1)$, see Sec. II.2.2 for a comparison with randomized benchmarking. Note that while the number of experiments scales exponentially with qubit number (not surprising for a tomographic protocol), the number of experiments needed for performing spectral tomography on single and two-qubit gates is comparable to the number of experiments that must be performed in randomized benchmarking on one or two qubits (which provides only average gate information). In randomized benchmarking one must sample $M$ random sequences for each sequence length $k\in[0:K]$, yielding $M\times(K+1)$ experiments. This $M$ is independent of the number of qubits helsen2017multi . In experiments $M$ is often chosen between $M\approx 40$ barends2014superconducting ; xue2019benchmarking at the low end and $M\approx 150$ at the higher end ballance2016high . Values of $K$ reported in randomized benchmarking experiments are also comparable to (or even higher than, see barends2014superconducting where $K\approx 300$ is considered) the values of $K$ used for single and two qubit spectral tomography (see Section II.3).

II.3 Spectral tomography on two superconducting chips

We have executed the spectral tomography method on a single-qubit $\pi/4$ rotation around the $X$-axis: $R_{x}(\pi/4)=\exp(-i\pi X/8)$. For this gate the ideal matrix $T^{R_{x}(\pi/4)}$ should have eigenvalues $1$, $\exp(i\pi/4)$ and $\exp(-i\pi/4)$. We execute this gate on two different systems available in the cloud: the two-qubit Quantum Infinity provided by the DiCarlo group at QuTech (for internal QuTech use) and the ibmqx4 (IBM Q5 Tenerife) available at https://quantumexperience.ng.bluemix.net/qx/editor. The results of this experiment are shown in Fig. 2 (left panel) in a polar plot which we refer to as the ‘spectral footprint’ of the gate. For clarity, in Fig. 2 (right panel) we have also plotted the phase deviation from ideal for the implemented gates.

On the two-qubit ($q_{0}$ and $q_{1}$) Quantum Infinity chip, we perform the single-qubit gate experiment on $q_{0}$ twice to study cross-talk: in one case the undriven qubit $q_{1}$ on the chip is in state $\ket{0}$, in the other case $q_{1}$ is in state $\ket{1}$. Since the residual off-resonant qubit coupling, mediated through a common resonator, is non-zero, we observe a small difference between these two scenarios, see Fig. 2. For the Quantum Infinity chip, when $q_{1}$ is $\ket{0}$ we estimate $\lambda_{j}^{\rm est}\in\{0.691+0.719i,0.691-0.719i,0.997\}$, while $\lambda_{j}^{\rm est}\in\{0.687+0.7239i,0.687-0.724i,0.998\}$ when $q_{1}$ is $\ket{1}$. Using the single-qubit fidelity bound given in Section IV.2, we can compute that the fidelity with respect to the targeted gate $R_{x}(\pi/4)$ can be no more than $0.999$ regardless of the state of $q_{1}$. We can also compute upper and lower bounds on the unitarity (see Section II.1 and Section IV.2) which yields $0.994\leq u\leq 0.996$ regardless of the state of $q_{1}$.

Regarding the ibmqx4 chip, the data are taken when all other qubits are in state $\ket{0}$. The reconstructed eigenvalues $\lambda_{j}^{\rm est}\in\{0.735+0.671i,0.735-0.671i,0.996\}$ turn out to be lower in magnitude. From these numbers we can conclude that the fidelity to the target gate is no higher than $0.998$ and the unitarity lies between $0.988$ and $0.991$.

For all these numbers a two-way $95\%$ confidence interval (for both real and imaginary parts) deviates by less than $0.005$ from the quoted values. The confidence intervals are obtained through a Wild resampling bootstrap with Gaussian kernel wu1986jackknife .

We have considered whether the data can be better fitted with more than $N=3$ eigenvalues. For each experiment we fit the data using $N$ eigenvalues with $N\in\{4,\ldots 15\}$ and we test whether there is a significant increase in goodness-of-fit using a standard F-test (sebernonlinear, , Section 2.1.5). For no experiment and value of $N$ does the resultant $p$-value drop below $0.05$, leading us to conclude that increasing the number of eigenvalues does not significantly increase the accuracy of the fit.

We have also executed a two-qubit CNOT gate on ibmqx4 (Fig. 3). The $T$ matrix of the ideal CNOT gate has 15 eigenvalues and a very degenerate spectrum: $6$ eigenvalues are equal to $-1$ and $9$ eigenvalues are equal to $1$, but our data, taking $K=50$, shows that a best fit is obtained using $4$ instead of $2$ eigenvalues. Using an F-test shows that the goodness-of-fit is significantly improved using 4 eigenvalues rather than 2 or 3, whereas adding more eigenvalues beyond $4$ does not significantly improve the goodness-of-fit $(p>0.05)$.

We have not tried using larger $K$ (which may lead to a resolution of more eigenvalues) since this would break the requirement that our experiments are executed as a single job performed in a short amount of time on the IBM Quantum Experience. The eigenvalues are $\lambda_{j}^{\rm est}\in\{0.939+0.059i,0.938-0.059i,-0.961+0.067i,-0.961-0.067i\}$, all with a $95\%$ confidence interval smaller than $\pm 3\times 10^{-3}$ for both real and imaginary parts. It is important to note that these 4 eigenvalues, coming in 2 complex-conjugate pairs, cannot be the spectrum of a two-qubit TPCP map $\mathcal{S}$, for the following reasons. As observed in Section II.1, the submatrix $T^{\mathcal{S}}$ of the Pauli transfer matrix of $\mathcal{S}$ is a real matrix of odd $(4^{2}-1=15)$ dimension. Since any complex eigenvalues of a real matrix come in conjugate pairs, $T^{\mathcal{S}}$ must have at least one real eigenvalue. Moreover, the data cannot be explained by allowing for leakage, as any eigenvalues associated to a small amount of leakage must have small associated amplitude, as we discuss in Section II.4. This is not the case for the eigenvalues plotted in Fig. 3 as all their amplitudes have comparable magnitude $A^{\rm est}\in\{3.34-1.70i,3.34+1.70i,1.57+0.91i,1.57-0.91i\}$.

In Section IV.3 we propose a simple model based on a frame mismatch accumulation that qualitatively reproduces these eigenvalues. This model is not stochastic but coherent, and it violates the assumption that the applied CNOT gate can be fully modeled as a TPCP map. A possible physical mechanism producing a frame mismatch accumulation can be a drift in an experimental parameter.

We do not compute bounds on the fidelity or unitarity of the CNOT gate since the bounds in Section II.1.1 do not apply when the evolution is non-Markovian.

II.4 Leakage and non-Markovian noise

In this section we consider how spectral tomography behaves under error models that violate the assumptions that go into Eq. 13.

Three common mechanisms for gate inaccuracy are (1) cross-talk, meaning the gate depends on or affects the state of other “spectator” qubits, (2) leakage, meaning that the dynamics of the gate acts outside of the computational qubit subspace and (3) non-Markovian dynamics, meaning that the assumption that $k$ applications of the noisy gate are equal to ${\cal S}^{k}$ for some TPCP map $\mathcal{S}$ is incorrect. Characterizing gates with respect to these features is important for assessing their use in multi-gate/multi-qubit devices for the purpose of quantum error correction or plainly reliable quantum computing WG:leak ; rol+:netzero .

One can see that all three scenarios are due to the dynamics taking place in a larger Hilbert space than the targeted computational qubit space. In the case of leakage, the larger space is an extension of the computational space, while in the other cases the larger space is the tensor product of the computational space with the state space of spectator qubits (1), as explored in Section II.3, or other quantum or classical degrees of freedom in the environment (3).

II.4.1 Leakage

Let us consider how gate leakage affects the signal $g(k)$, making the analysis for one or two qutrits. One can choose an operator basis for the qutrit space such as the basis of the 8 traceless (normalized) Gell-Mann matrices $\sigma_{\mu}^{\rm GM}$ for $\mu=1,\ldots,8$, together with the normalized identity $\sigma_{0}^{\rm GM}=\frac{1}{\sqrt{3}}I_{3}$. For a single qutrit, we can consider the ‘Pauli’ transfer matrix in this Gell-Mann basis, i.e. $S^{\rm GM}_{\mu\nu}={\rm Tr}[\sigma_{\mu}^{\rm GM}{\cal S}(\sigma_{\nu}^{\rm GM})]$ and its submatrix $T^{\rm GM}$.

For a single qutrit, the signal $g^{\mbox{\tiny NO SPAM}}(k)$ in Eq. (12) then equals ${\rm Tr}_{\rm comp}[(T^{\rm GM})^{k}]$ where ${\rm Tr}_{\rm comp}[A]$ represents the trace over a $3\times 3$ submatrix of $A$, corresponding to the Gell-Mann matrices which act like X, Y, and Z in the two-dimensional computational space. In other words, we can see the matrix $T^{\rm GM}$ as being composed of blocks:

[TABLE]

where the upper-left block is the sub-matrix whose trace we take in $g^{\mbox{\tiny NO SPAM}}(k)$. In the absence of other noise sources, $T^{\rm GM}$ corresponds to the evolution of a unitary gate and (assuming it is diagonalizable) it can be diagonalized by a rotation $V$ as $T^{\rm GM}=VDV^{-1}$, where $D$ is a diagonal matrix with all the eigenvalues $\{\lambda_{j}\}$. If we assume that leakage is low, meaning that $T_{\rm leak}$ and $T_{\rm seep}$ have small norm of $O(\epsilon)$, then at lowest order in $\epsilon$ the diagonalizing transformation $V$ will be block-diagonal, i.e. $V\approx V_{\rm comp}\oplus V_{\rm beyond}$. This means that g^{\mbox{\tiny NO SPAM}}(k)={\rm Tr}_{\rm comp}\big{[}(T^{\rm GM})^{k}\big{]}={\rm Tr}_{\rm comp}\big{[}VD^{k}V^{-1}\big{]}\approx\sum_{j=1}^{3}\lambda_{j}^{k}+O(\epsilon). Thus, at lowest order, the signal will have large amplitude on 3 relevant eigenvalues of the spectrum of $T^{\rm GM}$ and these eigenvalues could have been perturbatively shifted from their ideal location by low leakage. If the leakage is stronger, we can more generally write

[TABLE]

Here $\ket{\sigma_{j}}$ is a vector notation for one of the 8 Gell-Mann matrices $\sigma_{j}$ and $\Pi_{\rm comp}$ is the projector onto the basis spanned by the 3 Gell-Mann matrices which are the Paulis in the computational space. From this expression we see that the effect of leakage is the contribution of more eigenvalues to the signal $g(k)$. For low leakage we may expect three dominant eigenvalues with relatively large amplitude $\tilde{A}_{j}$ and five eigenvalues with small amplitude.

For a gate on two qutrits, identical remarks apply, except that an additional basis transformation is required from the orthogonal Gell-Mann basis to the computational qubit Pauli basis in order to keep the same division of $T^{\rm GM}$ as in Eq. 28. If we have two qutrits, the 80-dimensional traceless subspace is spanned by the matrices $\sigma_{\mu}^{\rm GM}\otimes\sigma_{\nu}^{\rm GM}$ for $\mu,\nu=0,\ldots,8$ except $\mu=\nu=0$. The issue is related to terms such as $\sigma_{0}^{\rm GM}\otimes\sigma_{\nu\neq 0}^{\rm GM}$ since the normalization of the qutrit identity ($\sigma_{0}^{\rm GM}=\frac{1}{\sqrt{3}}I_{3}$) is different from the normalization of the qubit identity ($P_{0}=\frac{1}{\sqrt{2}}I_{2}$). This suggests that for two qutrits it is better to write $T^{\rm GM}$ in a basis which includes the Pauli matrices in the computational subspace ($P_{\mu}\otimes P_{\nu}$ for $\mu,\nu=0,\ldots,3$ except $\mu=\nu=0$) as a sub-basis. For two qutrits, the signal may then contain up to 80 eigenvalues of which all but 15 are expected to have small amplitude in case of low leakage.

II.4.2 Non-Markovianity: time-correlated noise

Non-Markovian behavior of a gate can be due to temporal correlations in the classical or quantum environment of the driven qubit(s). Abstractly, we can include the environment in the gate action so that the evolution for each gate application is a unitary given by some $U_{\rm total}$ acting on system and environment. We can expand the Pauli transfer matrix of $U_{\rm total}$ in a Pauli basis for system and environment and view $T_{\rm comp}$ as a sub-block of $T_{\rm total}$, similar as in the case of leakage. Diagonalizing $T_{\rm total}$ and taking the trace over the computational space will result in an expression such as Eq. (29). For example, an additional spectator or environment qubit can lead to a signal $g(k)$ of a single-qubit gate having contributions from 15 eigenvalues. Choosing a sufficiently large $K$ may allow one to resolve these eigenvalues, even those with small amplitude.

A more malicious, but physically reasonable rol+:netzero , form of classical non-Markovian noise makes gate-parameters temporally correlated. In order to numerically study the effect of non-Markovian noise, we consider a toy example in which a perfect CZ gate is followed by a rotation around the $X$ axis on one qubit. For a series of $k$ repetitions of a perfect CZ gate, we assume that each one is followed by the same rotation $R_{x}(\phi)$ acting always on the same qubit. We assume that the angle $\phi$ is Gaussian-distributed with mean 0 and standard deviation $\sigma$: $\mathbb{P}_{\sigma}(\phi)=\exp(-\phi^{2}/2\sigma^{2})/\sqrt{2\pi}\sigma$. The time evolution for $k$ repetitions is then given by

[TABLE]

Note that $\mathcal{S}_{k}\neq(\mathcal{S}_{1})^{k}$ since this noise is correlated across multiple repetitions of the gate. Furthermore, we assume perfect state-preparation and measurement. In this case, one can represent the noisy gate by some unitary $U_{\rm total}$ acting on the two qubits and on a classical state in a Gaussian stochastic mixture of angles $\phi$. The continuous nature of this classical environment state leads to a lack of a hard cut-off on the number of eigenvalues in $g(k)$.

We apply the matrix-pencil method to the corresponding signal $g^{\mbox{\tiny NO SPAM}}(k)$ and we use an F-test to determine the optimal number of eigenvalues for each $\sigma$ (Fig. 4). For $\sigma=22.9^{\circ}$ and $K=50$ we find eigenvalues with modulus clearly larger than 1. Those are unphysical but not excluded by the matrix-pencil method. We expect that such $|\lambda^{\rm est}|>1$ disappear when considering a longer signal, since $g(k)$ does not increase exponentially in $k$. In other words, this is a sign that the signal contains more spectral content than can be resolved from the time scale set by $K$. Indeed, for $\sigma=22.9^{\circ}$ we have made the same analysis for larger $K$’s up to $K=200$ and we find that those eigenvalues get closer and closer to 1. If instead we fix $K=50$ and consider different $\sigma$’s, we find that for a low $\sigma$ (e.g. $5.7^{\circ}$) unphysical eigenvalues are not present (Fig. 4), whereas for $\sigma>22.9^{\circ}$ (e.g. $40.1^{\circ}$) they get again closer and closer to 1. This latter fact can be understood by noting that increasing $\sigma$ is analogous to enlarging the time scale set by $K$, as the characteristic time scale of dephasing gets shorter for a fixed $K$. Based on these observations, we conclude that there is a certain intermediate time scale at which eigenvalues larger than 1 are extrapolated from the data in the presence of sufficiently-strong non-Markovian noise of the kind described in this section. Section IV.3 discusses a model with a different kind of time-correlation leading to a spectral footprint which is incompatible with that of a TPCP map.

II.4.3 Non-Markovianity: coherent revivals

In order to better understand the occurrence of eigenvalue estimates $|\lambda^{\rm est}|>1$, we apply the matrix-pencil method on a signal (of a somewhat different physical origin), which has a revival over the time period set by $K$.

It is well-known that in the exchange of energy between a two-level atom with a bosonic mode, the Rabi oscillations of the two-level atom are subject to temporal revivals. These revivals are due to the fact that the bosonic driving field is not purely classical, but rather gets entangled with the state of the qubit via the Jaynes-Cummings interaction. In particular, for a coherent driving field with coherent amplitude $\alpha$ with average photon number $\bar{n}=|\alpha|^{2}$, the probability for the atom to be excited equals (see Section 3.4.3 in book:HR ):

[TABLE]

with $p_{\alpha}(n)=\exp(-|\alpha|^{2})\frac{|\alpha|^{2n}}{n!}$. We consider $\bar{n}=5$ and sample the damped oscillatory function $P_{e}(t)-\frac{1}{2}$ at regular intervals $k\Omega\delta t$ with $\Omega\delta t=0.05$ and $k=0,\ldots,K=900$. The signal function $g(t)=P_{e}(t)-\frac{1}{2}$ contains eigenvalues equal to $\lambda_{n}=\exp(\pm i\,0.05\sqrt{n+1})$ with amplitudes according to the Poisson distribution $p_{\alpha}(n)$ with mean photon number $\bar{n}$.

We observe that the matrix-pencil method finds eigenvalues larger than $1$, see Fig. 5, which contribute significantly ($p<0.01$ via F-test) to the reconstructed signal. We can understand this feature of eigenvalues exceeding 1 as a way in which the matrix-pencil method handles revivals: the signal has more spectral content than what can be resolved from the window of time given by $K$, in particular there is no hard cut-off on the number of eigenvalues which contribute. We have observed that an analysis of the signal over a longer period of time, that is, a larger $K$ up to $K=5000$, gives eigenvalues whose norm converges to at most 1.

III Discussion

We have introduced spectral quantum tomography, a simple method that uses tomographic data of the repeated application of a noisy quantum gate to reconstruct the spectrum of this quantum gate in a manner resistant to SPAM errors. We have experimentally validated our method on one- and two-qubit gates and have also numerically investigated its behavior in the presence of temporally-correlated non-trivial error models.

The effective upshot of leakage and non-Markovian noise is that the signal will have more spectral content than what can be resolved given a chosen sequence length $K$, leading to unphysical features in the spectrum such as an eigenvalue estimate larger than 1, or the absence of a real eigenvalue. Even though we have seen in our examples that a physical spectrum can be regained by going to larger $K$, depending on the noise model, this convergence may be very slow requiring much data-taking time. Hence these unphysical features are useful markers for deviations from our model of repeated TPCP qubit maps ${\cal S}^{k}$. We view it as an open question how well one can reliably distinguish different sources of deviations.

III.1 Logical Spectral Quantum Tomography

An interesting application of the spectral tomography method could be the assessment of logical gates on encoded quantum information in a SPAM-resistant fashion. In this logical scenario (for, say, a single logical qubit), one first prepares the eigenstates of the logical Pauli operators $\overline{X},\overline{Y}$ and $\overline{Z}$. One then applies a unit of error-correction $k=0,\ldots,K$ times: a single unit could be, say, the repeated error correction for $L$ rounds of a distance-$L$ surface code. Or a unit is the application of a fault-tolerant logical gate, e.g. by means of code-deformed error correction or a transversal logical gate followed by a unit of error correction. After $k$ units one measures the logical Pauli operators fault-tolerantly, and repeats experiments to obtain the logical signal $\overline{g}(k)$. Studying the spectral features of such logical channel will give information about the efficacy of the quantum error correction unit and/or the applied logical gate while departures from the code space or a need to time-correlate syndrome data beyond the given QEC unit can show up as leakage and non-Markovian errors.

IV Methods

IV.1 Single-qubit case with non-diagonalizable matrix $T$

In general, a matrix $T$ can be brought to Jordan normal form by a similarity transformation, i.e. $T=VJV^{-1}$ with $J=\oplus_{i}J_{i}$ where each Jordan block $J_{i}$ is of the form

[TABLE]

see e.g. Theorem 3.1.11 in book:HJ . $T$ is diagonalizable when each Jordan block is fully diagonal.

An example of a non-diagonalizable Lindblad superoperator on a single qubit has been constructed in SL:nondiag . Using this, one can easily get a single-qubit superoperator ${\cal S}$ for which the traceless block of the Pauli transfer matrix is a non-diagonalizable matrix $T$ as follows. Let ${\cal S}(\rho)=\exp({\cal L}\epsilon)(\rho)\approx\rho+\epsilon{\cal L}(\rho)+O(\epsilon^{2})$ with ${\cal L}(\rho)=-i[\frac{yZ}{2},\rho]+{\cal D}[(2x)^{1/2}\sigma_{-}](\rho)+{\cal D}[y^{1/2}X](\rho)$ with ${\cal D}[A](\rho)=A\rho A^{\dagger}-\frac{1}{2}\{A^{\dagger}A,\rho\}$ and real parameters $x,y\geq 0$. This implies that ${\cal S}$ has the $4\times 4$ Pauli transfer matrix

[TABLE]

Taking some small $\epsilon$ and $x\neq 0$, one can check that the submatrix $T$ does not have 3 eigenvectors and it has a pair of degenerate eigenvalues, so $T$ is not diagonalizable. When we take $x=0$, ${\cal S}$ is unital, that is ${\cal S}(I)=I$, and the submatrix $T$ is not diagonalizable either.

Even though a matrix $T$ is not always diagonalizable, there still exists the so-called Schur triangular form for any matrix $T$ book:HJ . This form says that $T=W(D+E)W^{\dagger}$ with $W$ a unitary matrix, $D$ a diagonal matrix with the eigenvalues of $T$, and $E$ a strictly upper-triangular “nilpotent” matrix with non-zero entries only above the diagonal. Since the $N\times N$ matrix $E$ is strictly upper-triangular, one has {\rm Tr}\big{[}D^{i}E^{j}\big{]}=0 for all $j\neq 0$. If we use this form in Eq. (12), one obtains for any $k$

[TABLE]

since any product of the form $D^{l_{1}}E^{l_{2}}D^{l_{3}}\ldots E^{l_{m}}$ with some non-zero $l_{i}>0$ is a matrix with zeros on the diagonal. In case of SPAM errors and non-diagonalizable $T$ we consider

[TABLE]

where $W^{\dagger}T_{\rm prep}T_{\rm meas}W$ is not the identity matrix due SPAM errors, implying that $g(k)$ can depend on $E$ and have a non-exponential dependence on $k$. Thus, in the special case of a non-diagonalizable matrix $T$, the signal $g(k)$ would not have the dependence on the eigenvalues as in Eq. (13).

In particular, we can examine the physically-interesting non-diagonalizable Case 3 in Section II.1.2 in this light, taking $h_{y}=h_{z}=0$ and a critical $h_{x}^{\mathrm{cr}}=\frac{\Gamma_{1}-\Gamma_{2}}{2}$. The dynamics of the Lindblad equation after time $t$ induces a superoperator ${\cal S}_{t}$ which will have the following action on the Pauli operators:

[TABLE]

Here we can note the linear dependence on $t$ due to the system being critically damped. If we consider the signal $g(t)=\sum_{\mu}{\rm Tr}[P_{\mu}{\cal S}_{t}(P_{\mu})]$ we see that this linear dependence on $t$ drops out in accordance with Eq. (33), i.e. this trace only depends on the eigenvalues and has an exponential dependence on $t$. In the presence of SPAM errors, some of the linear dependence could still be observable for such critically-damped system. In addition, coefficients such as $c_{\mu\nu}(t)={\rm Tr}[P_{\mu}{\cal S}_{t}(P_{\nu})]$ can depend linearly on $t$, making such tomographic data less suitable to extract eigenvalue information.

IV.2 Upper bound on the entanglement fidelity with the targeted gate

In this section we show how to relate the eigenvalues of the Pauli transfer matrix of a TPCP map ${\cal S}$ to an upper bound on the entanglement fidelity (and hence the average gate fidelity) with the targeted unitary gate $U$. Naturally, one can only expect to obtain an upper bound on the gate fidelity, since the eigenvalues do not provide information about the eigenvectors of ${\cal S}$. If the actual eigenvectors deviate a lot from ideal, the actual gate fidelity could be very low, so one can certainly not derive a lower bound on the fidelity based on the eigenvalues.

Lemma 1**.**

Let the eigenvalues of the $N\times N$ matrix $T^{\cal S}$ be $\{\lambda_{i}\}_{i=1}^{N}$ with $N=d^{2}-1$ for a $d$-dimensional system. Let $U$ be the targeted gate with eigenvalues $\{\lambda_{i}^{\rm ideal}\}_{i=1}^{N}$ and let there be permutation $\pi$ of $i$-th eigenvalue $\lambda_{i}$ which maximizes $|\sum_{i}\lambda_{\pi(i)}^{*}\lambda_{i}^{\rm ideal}|$ so that $0\leq\xi_{\rm max}=\max_{\pi}\frac{1}{N}|\sum_{i}\lambda_{\pi(i)}^{*}\lambda_{i}^{\rm ideal}|\leq 1$. The entanglement fidelity {\cal F}_{\rm ent}(U,{\cal S})=\frac{1}{N+1}(1+{\rm Tr}\big{[}T^{U^{\dagger}}T^{\cal S}\big{]}) is upper bounded as

[TABLE]

where $u(\mathcal{S})$ is the unitarity of $\mathcal{S}$.

Proof.

We write $T^{\cal S}$ in Schur triangular form as $T^{{\cal S}}=W(D^{\cal S}+E)W^{\dagger}$ with $W$ a unitary matrix, $D^{\cal S}$ a diagonal matrix with the eigenvalues of $T^{\cal S}$, and $E$ a strictly upper-triangular “error” matrix with non-zero entries only above the diagonal book:HJ . Using the Cauchy-Schwartz inequality one has

[TABLE]

Note that for a unitary gate $U$, $T^{U^{\dagger}}=(T^{U})^{T}=(T^{U})^{\dagger}$ and $T^{U^{\dagger}}T^{U}=I$ implying that $T$ is an orthogonal matrix with unit singular values. We thus have ({\rm Tr}\big{[}{T^{U}}^{\dagger}T^{U}\big{]})^{1/2}=\sqrt{N}. One has {\rm Tr}\big{[}{T^{\cal S}}^{\dagger}T^{\cal S}\big{]}={\rm Tr}\big{[}(D^{\cal S}+E)^{\dagger}(D^{\cal S}+E)\big{]}={\rm Tr}\big{[}({D^{\cal S}}^{\dagger}D^{\cal S}+E^{\dagger}E)\big{]}, using the strict upper-triangularity of $E$. In other words, {\rm Tr}\big{[}E^{\dagger}E\big{]}={\rm Tr}\big{[}{T^{\cal S}}^{\dagger}T^{\cal S}\big{]}-\sum_{i}|\lambda_{i}|^{2} where $\lambda_{i}$ are the eigenvalues of $T^{\cal S}$.

Recognizing that \frac{1}{N}{\rm Tr}\big{[}{T^{\cal S}}^{\dagger}T^{\cal S}\big{]}=u(\mathcal{S}), we obtain an upper bound on the second term in Eq. (36).

Now let’s upper bound the first term in Eq. (36) for unknown unitary $W$. Assume w.l.o.g. that $T^{U}$ and $D^{\cal S}$ are diagonal in the same basis (the additional rotation between these eigenbases can be absorbed into $W$). Let $T^{U}=\sum_{i}\lambda_{i}^{\rm ideal}P_{i}$ and $D^{\cal S}=\sum_{i}\lambda_{i}P_{i}$ with orthogonal projectors $P_{i}$ and $\sum_{i}P_{i}=I$. Define the matrix $M$ with entries M_{ij}={\rm Tr}\big{[}P_{i}WP_{j}W^{\dagger}\big{]}. The matrix $M$ is doubly-stochastic, since $\sum_{i}M_{ij}=1=\sum_{j}M_{ij}$ which implies that $M=\sum_{m}q_{m}\pi_{m}$ with $q_{m}\geq 0,\sum_{m}q_{m}=1$ (Birkhoff-von Neumann theorem book:HJ ) with permutation matrix $\pi_{m}$. With these facts and the convention $\bra{i}\lambda^{\cal S}\rangle=\lambda_{i}$ we can bound

[TABLE]

These bounds together then lead to Eq. (35). ∎

An immediate corollary of Lemma 1 is

[TABLE]

since $u(\mathcal{S})\leq 1$ for TPCP maps. However, this is in general not a very strong upper bound on the fidelity.

We can do better in the single-qubit case by realizing that there are strong relations between the singular values $\sigma_{i}$ of $T^{\mathcal{S}}$ and the absolute values of the eigenvalues $|\lambda_{i}|$ of $T^{\mathcal{S}}$. Ordering both the singular values and the eigenvalue magnitudes in descending order, we have the following (weak Majorization) inequalities for arbitrary matrices

[TABLE]

For single-qubit channels we can also impose TPCP constraints to the singular values of the channel. In particular we have (wolf2010inverse, , Eq. (4))

[TABLE]

Using these relations we can upper bound the unitarity of a single-qubit channel $\mathcal{S}$, given its eigenvalues, using the optimization:

[TABLE]

This is a non-convex optimization problem in three variables, for which a global minimum can be numerically computed given $\lambda_{1},\lambda_{2},\lambda_{3}$. This gives an upper bound on the unitarity of $\mathcal{S}$ and hence on the entanglement fidelity of $\mathcal{S}$ to the target unitary $U$. In the main text we use this optimization to give non-trivial upper bounds on the fidelities of single-qubit gates realized on superconducting chips and analyzed using the spectral tomography method.

IV.3 Frame Mismatch Accumulation

In Section II.3 we noted that the data gathered for the CNOT gate cannot be explained by a model of a noisy TPCP map $\mathcal{S}$ repeated $k$ times. Here we propose a simple coherent model that qualitatively reproduces the features observed in Fig. 3 and we call this the frame mismatch accumulation model. Let $\mathcal{S}_{0}$ be a TPCP map that is a good approximation of the targeted gate applied exactly once (in the main text this was the CNOT) and let $V$ be some unitary. In the frame mismatch accumulation model we assume that $k$ consecutive applications of the gate are equal to:

[TABLE]

Intuitively, this can be interpreted as an increasing mismatch between the frame in which $\mathcal{S}_{0}$ was defined and the frame in which the gate was implemented at the $i$-th repetition, up to $i=k$.

We apply this model to a CNOT gate, choosing $\mathcal{S}_{0}$ to be an ideal CNOT gate and choosing $V=\exp(-i\frac{\theta}{2}I~{}\otimes~{}Y)$ with $\theta=0.05$ deg. In the case of the cross-resonance CNOT gate performed on ibmqx4, this may correspond to an imperfect cancellation of the $I\otimes Y$ term Sheldon16 . In Fig. 6 we see that this example closely reproduces the eigenvalues shown in Fig. 3. At the same time, we note that the qualitative features observed in Fig. 6 do not depend on the choice of the rotation axis of $V$ (for either qubit), as long as the rotation does not commute with $\mathcal{S}_{0}$ (which would leave the gate unaffected by the frame mismatch).

V Data Availability

Experimental data gathered for Figs. 2 and 3, as well as an implementation of the matrix pencil algorithm can be found online at https://doi.org/10.5281/zenodo.2613856.

VI Competing interests

The authors declare no competing interests. The views expressed in this manuscript are those of the authors and do not reflect the official policy or position of IBM or the IBM Quantum Experience Team

VII Author contributions

All authors contributed to the development of the theoretical concepts presented. The experiments on the IBM QE and the QuTech QI were performed and analysed by JH and the simulations on non-Markovianity were performed by JH and FB, under supervision of BMT. All authors contributed to the writing of the manuscript.

VIII Acknowledgements

The work by FB and BMT was supported by ERC grant EQEC No. 682726. JH is funded by STW Netherlands, NWO VIDI, an ERC Starting Grant and by the NWO Zwaartekracht QSC grant. We thank Andrew Cross for generous access to the IBM Quantum Experience for a TU Delft MSc project which indirectly led to this work, we thank Jarn de Jong for discussions on quantum tomography and Adriaan Rol for assistance with the QuTech Quantum Infinity. JH would also like to thank Arnaud Carignan-Dugas for interesting discussions.

IX Additional information

The authors declare no competing interests. The views expressed in this manuscript are those of the authors and do not reflect the official policy or position of IBM or the IBM Quantum Experience Team.