# Learning a quantum channel from its steady-state

Yigal Ilin1<EMAIL_ADDRESS>Itai Arad2,3,4 1Andrew and Erna Viterbi
Department of Electrical and Computer Engineering, Technion – Israel Institute
of Technology, 3200003, Haifa, Israel  2Physics Department, Technion – Israel
Institute of Technology, 3200003, Haifa, Israel  3The Helen Diller Quantum
Center, Technion – Israel Institute of Technology, 3200003, Haifa, Israel
4Centre for Quantum Technologies, National University of Singapore, 117543
Singapore, Singapore

###### Abstract

We present a scalable method for learning local non-unital quantum channels
using only local expectation values measured on their steady state. Our method
can be seen as the channel equivalent of the algorithms for learning local
Hamiltonians from their ground states. Such channels can be efficiently
implemented on the quantum computer using mid-circuit measurements or RESET
gates. Using engineered dissipative dynamics, our method provides a simple way
to check how well a given noise model describes the actual noisy dynamics
within the hardware, when all qubits are actively used for a non-negligible
amount of time. Moreover, given a parameterized noise model, our method can be
applied to learn the underlying noise parameters for the entire system. We
demonstrate our method numerically and experimentally on an IBMQ machine, and
show that a full noise model can be verified and learned from Pauli
measurements on a _single_ state.

##  I Introduction

In recent years, with the advances in quantum information and quantum
computation, there have been a growing interest in methods for learning many-
body quantum systems. These methods include, for example, learning a many-body
Hamiltonian based on its dynamics [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11], its
steady states [12, 13, 14, 15], or using a trusted quantum simulator [16, 17,
18, 19, 20].

Among these methods, a promising line of research tries to learn the
underlying generators of dynamics using local measurements on a steady state
of the system [21, 22, 23, 24, 25, 26, 27]. It can be shown that if the
underlying Hamiltonians or Lindbladians are sufficiently generic, they can be
uniquely determined by the expectation values of local observables in these
steady states. Specifically, Ref. [25] showed that for a system governed by a
local Hamiltonian, there exists a set of linear constraints between
expectation values of local observables in a steady state of the system, from
which an unknown local Hamiltonian can be uniquely determined. Following that
work, Ref. [26] showed that a similar set of constraints can be obtained for
the case of open quantum systems, allowing reconstruction of unknown local
Lindbladians by performing local measurements on a steady state. Together,
Refs. [25, 26] imply that by measuring local expectation values from a single
state within a given patch of the system, we can infer the underlying
Hamiltonian or Lindbladian terms in that region 111Strictly speaking, this is
true up to an overall factor..

In some respect, at least in the case of local Hamiltonians, the method relies
on the non-commutativity of the local Hamiltonian terms, and is therefore a
‘purely quantum’ method, with no classical equivalent. It gives a surprisingly
simple and efficient method for learning the generators of the local dynamics,
and can also be used to test how well a particular model of the Hamiltonian or
Lindbladian agrees with the empirical measurements. All this is done _without
having to classically calculate local expectation values_ , which is crucial
for the scalablity of the approach.

The simplicity and efficiency of the methods developed in Refs. [25, 26]
raises the question of whether and how they can be applied to discrete quantum
channels, in particular those that describe, or engineered on quantum
computers. Is it possible to find constraints on local expectation values in
these states that can be used to efficiently learn the underlying quantum
channel?

A quantum computer has several key distinctions from local Hamiltonians and
Lindbladians. The latter evolve _continuously_ , whereas, at least in the
abstract level, a quantum computer evolves in a _discrete_ , gate-based
fashion. In addition, steady states of a unitary circuit are convex
combination of its eigenstates, and it is not clear how to prepare them
efficiently. Moreover, the presence of noise changes this picture
significantly by breaking unitarity and the purity of the underlying quantum
states.

In this work we propose a method that overcomes these obstacles and
generalizes the approach of Refs. [25, 26] to discrete channels that describe
gate-based quantum computers. Our main idea is to engineer a dissipative
channel that can be implemented on a quantum computer using non-unitary mid-
circuit measurement operations, such as the $\operatorname{RESET}$ gates
available on the IBMQ machines. By iteratively applying this channel, we
rapidly reach a fixed point, on which we can measure local expectation values.
These expectation values satisfy an extensive set of local linear constraints
among themselves, which can be efficiently verified _without having to
classically estimate these expectation values_. Therefore, these constraints
can be used to learn the underlying channel that describes the quantum
computer by examining only its fixed point represented by a single quantum
state.

Like the methods in Refs. [24, 25, 26, 27], our method is highly scalable, and
can be used on systems with a large number of qubits. To demonstrate it, we
introduce two possible dissipative maps. The first map is a stochastic map,
which can be described as a local Kraus map (at least in the absence of
noise). It can be implemented on a quantum computer by randomly applying a
gate from a set of gates with a prescribed set of probabilities. Consequently,
the steady state is probed by running many different random circuits, which
realize different trajectories of the dynamics. The second map is a
deterministic map that is based on a composed 2-local gate, which is both non-
unitary and entangling. Implementing this map on a quantum computer requires
only one circuit, which makes it easier to execute on currently available
hardware such as IBMQ. To test our method we used numerical simulations of
both maps on $5$–$11$ qubits, and in addition applied the deterministic map to
an IBMQ machine using $5$ qubits.

The structure of this paper is as follows. In Sec. II we provide brief
definitions and basic facts on quantum channels and gates, and the notation we
use throughout the paper. Then in Sec. III we introduce the theoretical
framework of learning a local quantum channel from its steady state, and
explain how such framework can be used, for example, to learn and validate the
noisy dynamics of a quantum computer. In Sec. IV we focus on the learning and
validation task by defining the noise model we will be using to characterize
and validate (benchmark) the noisy dynamics of IBMQ machines. The optimization
method used in the variational learning part is also described. In Sec. V we
describe the results of our numerical simulations, performed on the two types
of engineered dissipative dynamics. In Sec. VI we present the results of
applying our method on an actual IBMQ machine using 5 qubits, which include
benchmarking different noise models, as well as using our method to learn its
noise. Finally, in Sec. VII we present our conclusions and discuss possible
future research directions.

##  II Preliminaries

Throughout this work we consider a quantum computer (QC) made of $n$ qubits
with Hilbert space $(\mathbb{C}^{2})^{\otimes n}$. Quantum states will be
generally denoted by Greek letters $\rho$, $\sigma$, etc. The expectation
value of an observable $A$ with respect to the quantum state $\rho$ will be
generally denoted by ${\langle{A}\rangle}_{\rho}=\operatorname*{Tr}(\rho A)$.

Unitary gates will be denoted by $X,Y,Z,\operatorname{CX},R_{x},$ etc. We use
$\operatorname{CX}(i\to j)$ to denote the CNOT (Controlled-NOT) gate between
the control qubit $i$ and the target qubit $j$. We use $R_{\bm{x}}(\theta)$ to
denote a Bloch sphere rotation of angle $\theta$ around the axis $\bm{x}$.
Quantum channels, i.e., completely positive, trace-preserving (CPTP) maps will
usually be denoted by calligraphic letters such as
$\mathcal{E}(\cdot),\mathcal{N}(\cdot)$, etc. Given a unitary gate $U$, we
will denote its corresponding channel by $\mathcal{U}$, i.e.,
$\mathcal{U}(\rho)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}U\rho
U^{\dagger}$. For any quantum channel $\mathcal{E}$, its adjoint is the unique
superoperator $\mathcal{E}^{*}$ that satisfies
$(\mathcal{E}^{*}(A),B)=(A,\mathcal{E}(B))$ for any operators $A,B$, where
$(A,B)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\operatorname*{Tr}(A^{\dagger}B)$
is the Hilbert-Schmidt inner product between two operators. A map
$\mathcal{E}$ is trace preserving if and only if
$\mathcal{E}^{*}(\mathbb{I})=\mathbb{I}$. A channel $\mathcal{E}(\cdot)$ is
called _unital_ if it also holds that $\mathcal{E}(\mathbb{I})=\mathbb{I}$. In
other words, the maximally mixed state is a fixed point of the channel.

An important example of a non-unital channel that can be applied on a QC is
formed by the $\operatorname{RESET}$ gate, which performs an active mid-
circuit reset of a qubit to the ${|{0}\rangle}$ state. Ideally, it corresponds
to the channel

$\displaystyle\mathcal{E}_{\text{RESET}}(\rho)={|{0}\rangle\\!\,\langle{1}|}\rho{|{1}\rangle\\!\,\langle{0}|}+{|{0}\rangle\\!\,\langle{0}|}\rho{|{0}\rangle\\!\,\langle{0}|}.$

It can also be realized by measuring the qubit in the standard basis, and
applying $X$ gate whenever the result is ${|{1}\rangle}$.

##  III Learning quantum channels from their steady state

In this section we introduce our method for learning quantum channels using
local measurements on their steady state. We shall consider channels that can
be realized on a quantum computer. A summary of the main steps of our method
is given in Fig. 1. Two examples to keep in mind are (i) the application of
quantum gate drawn from a set of gates according to a fixed probability
distribution, and (ii) a low-depth quantum circuit. In both examples, we
assume that the overall channel is non-unital because of the noise and/or
because of actively using the RESET gate. Unlike unital channels, whose fixed
point is typically the maximally mixed state, the fixed points of non-unital
channels are typically non-trivial in the sense that two different channels
from a relevant set of channels will have different fixed points. We shall
further assume that $\mathcal{E}(\cdot)$ has a _unique_ fixed point
$\mathcal{E}({\rho_{\infty}})={\rho_{\infty}}$, which we call the _steady
state_. We shall also assume that ${\rho_{\infty}}$ can be approximately
reached by the QC using few applications of $\mathcal{E}$, starting from an
initial state $\rho_{0}={|{0}\rangle\\!\,\langle{0}|}^{\otimes n}$. For
brevity, we will denote expectation values with respect to ${\rho_{\infty}}$
by ${\langle{\cdot}\rangle}_{\infty}$ instead of
${\langle{\cdot}\rangle}_{\rho_{\infty}}$.

Figure 1: Block diagram summarizing our method. Given a description of the
non-unital channel $\mathcal{E}(\cdot)$ (e.g. its quantum gates, activation
probabilities, etc), we iteratively apply it on a QC to reach an approximate
steady state. The local Pauli expectation values are sampled from the steady
state, and used as input to the classical post-processing routine. They
satisfy a set of local constraints with coefficients that depend on the noise
parameters of the underlying quantum device. The local constraints can be
verified efficiently and their number grows linearly with the system size. The
classical post-processing routine is readily tailored to learn specific
channel properties such as rotation angles, coherence times, gate activation
probabilities and so on.

Given a channel $\mathcal{E}$, assume that the system is in its steady state
${\rho_{\infty}}$. Then the expectation value of any observable $A$ must
satisfy

$\displaystyle{\langle{A}\rangle}_{\infty}=\operatorname*{Tr}(A{\rho_{\infty}})=\operatorname*{Tr}(\mathcal{E}^{*}(A){\rho_{\infty}})={\langle{\mathcal{E}^{*}(A)}\rangle}_{\infty}.$
(1)

The above equation is the basis of our method. For the channels that will be
discussed below, when $A$ is a local observable, the observable
$\mathcal{E}^{*}(A)$ is also local, or at least can be estimated using local
measurements. In such case, we can use the equality between the LHS and RHS to
constrain, or even learn $\mathcal{E}$, without having to calculate the
expectation values in $\rho_{\infty}$ — which is in general exponentially
expensive in $n$. Moreover, as we are measuring the steady state of the
system, our results are independent from the initial state of the system, and
therefore insensitive to state preparation errors.

Below, we describe two methods to construct such $\mathcal{E}$, which we call
the stochastic method and the deterministic method. In Sec. III.1 and Sec.
III.2 we describe the first method, and in Sec. III.3 describe the second.

###  III.1 A stochastic map: the strictly local case

Let $\\{\mathcal{E}_{k}\\}$ be a set of local channels implementable by a QC.
We can think of every $\mathcal{E}_{k}$ as being the application of a quantum
gate or a small circuit acting on a few contiguous qubits. We assume that the
$\mathcal{E}_{k}$ channels are strictly local in the sense that they act
trivially outside the support of the small circuit that defines them. This
assumption breaks down in realistic quantum hardware, where additional noise
channels act _simulatiniously on all_ qubits. In the next subsection we will
show how to extend our method to many of these cases (where the noise is local
and Markovian), but for the sake of keeping the presentation clear, we start
by assuming that the $\mathcal{E}_{k}$ channels are strictly local, acting on
at most two neighboring qubits. In addition, we allow some of the
$\mathcal{E}_{k}$ to be non-unital by containing RESET gates.

Together with $\\{\mathcal{E}_{k}\\}$, we will use a probability distribution
$\\{p_{k}\\}$ and define our quantum channel by

$\displaystyle\mathcal{E}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\sum_{k}p_{k}\mathcal{E}_{k}.$
(2)

Physically, $\mathcal{E}$ can be implemented on a QC in a stochastic way: we
choose a $k$ with probability $p_{k}$ and apply $\mathcal{E}_{k}$.

For example, the three qubits non-unital quantum channel

$\displaystyle\mathcal{E}(\rho)$ $\displaystyle=0.2\cdot
X^{1/2}_{1}\rho{X^{{1/2}}_{1}}^{\dagger}+0.2\cdot\operatorname{CX}(1\to
2)\rho\operatorname{CX}(1\to 2)$ $\displaystyle+0.4\cdot\operatorname{CX}(2\to
3)\rho\operatorname{CX}(2\to 3)$
$\displaystyle+0.2\cdot\operatorname{RESET}_{3}(\rho)$

can be implemented on a QC by randomly applying one of the gates
$\\{X^{1/2}_{1},\operatorname{CX}(1\to
2),\operatorname{CX}(2\to~{}3),\operatorname{RESET}_{3}\\}$ according to the
probabilities $\\{0.2,0.2,0.4,0.2\\}$. A simple numerical check shows that the
steady state of the this channel is a non-product state which can be
approximated to within a trace distance of $10^{-5}$ using few tens of steps.

Plugging Eq. (2) in Eq. (1), we see that for any local observable $A$,

$\displaystyle{\langle{A}\rangle}_{\infty}=\sum_{k}p_{k}{\langle{\mathcal{E}_{k}^{*}(A)}\rangle}_{\infty}.$
(3)

Importantly, if $\mathcal{E}_{k}$ is a $2$-local channel, then
$\mathcal{E}_{k}^{*}$ also acts non-trivially only on the two-qubits on which
$\mathcal{E}_{k}$ is defined. Consequently, if $A$ is a $t$-local observable,
then $\mathcal{E}_{k}^{*}(A)$ is at most a $(t+2)$-local observable. In fact,
it is at most $(t+1)$ for the following reason. The only case where
$\mathcal{E}_{k}^{*}(A)$ might be $(t+2)$-local is when the supports of
$\mathcal{E}_{k}$ and $A$ are disjoint. But in such case,
$\mathcal{E}_{k}^{*}(A)=A\mathcal{E}_{k}^{*}(\mathbb{I})=A$, which is
$t$-local.

Using this observation, we may write Eq. (3) as

$\displaystyle{\langle{A}\rangle}_{\infty}=\sum_{k\in\operatorname{supp}(A)}p_{k}{\langle{\mathcal{E}_{k}^{*}(A)}\rangle}_{\infty}+\sum_{k\notin\operatorname{supp}(A)}p_{k}{\langle{A}\rangle}_{\infty},$

where we use the notation $k\in\operatorname{supp}(A)$ to denote enumeration
over all $k$ indices for which $\mathcal{E}_{k}$ acts non trivially on $A$.
Writing the LHS as
$\sum_{k\in\operatorname{supp}(A)}p_{k}{\langle{A}\rangle}_{\infty}+\sum_{k\notin\operatorname{supp}(A)}p_{k}{\langle{A}\rangle}_{\infty}$
and re-organizing the equation, we obtain the following equation, which holds
for every local observable $A$:

$\displaystyle\sum_{k\in\operatorname{supp}(A)}p_{k}\big{(}{\langle{\mathcal{E}_{k}^{*}(A)}\rangle}_{\infty}-{\langle{A}\rangle}_{\infty}\big{)}=0.$
(4)

The above equation can be used to _validate_ a given model of the local
channels $\\{\mathcal{E}_{k}\\}$. Indeed, given a model for the local
channels, together with a corresponding probability distribution
$\\{p_{k}\\}$, we define a _cost function_ $\Phi$ by taking the square of the
difference between the LHS and RHS of Eq. (4) for a set of independent local
observables $\\{A\\}$:

$\displaystyle\Phi\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\sum_{A}\Big{(}\sum_{k\in\operatorname{supp}(A)}p_{k}\big{(}{\langle{\mathcal{E}_{k}^{*}(A)}\rangle}_{\infty}-{\langle{A}\rangle}_{\infty}\big{)}\Big{)}^{2}.$
(5)

Since ${\langle{\mathcal{E}_{k}^{*}(A)}\rangle}_{\infty}$ is the expectation
value of a $(t+1)$-local observable, we can express it as a linear combination
of the expectation value of $(t+1)$-local Pauli strings on the same support.
The resulting $\Phi$ can then be written as a quadratic expression of these
Pauli expectation values, with coefficients $C_{\alpha\beta}$ that depend on
the details of underlying model:

$\displaystyle\Phi=\sum_{A}\sum_{\alpha,\beta}C^{(A)}_{\alpha\beta}{\langle{P_{\alpha}}\rangle}_{\infty}{\langle{P_{\beta}}\rangle}_{\infty}.$
(6)

Therefore, by measuring the $(t+1)$-local Pauli expectation values in the
steady state, we can calculate $\Phi$ and see how well the model fits the
actual quantum hardware.

This procedure can also be made local by considering local cost functions
$\Phi_{q}$ for qubit $q$. Denoting by $I_{q}$ the set of observables $A$ for
which $\mathcal{E}_{k}^{*}(A)$ acts non-trivially on $q$, the local cost
function is given by

$\displaystyle\Phi_{q}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\frac{1}{|I_{q}|}\sum_{A\in
I_{q}}\Big{(}\sum_{k\in\operatorname{supp}(A)}p_{k}{\langle{\mathcal{E}_{k}^{*}(A)}\rangle}_{\infty}-{\langle{A}\rangle}_{\infty}\Big{)}^{2},$
(7)

where the $|I_{q}|$ is the number of elements in $I_{q}$. In general,
$|I_{q}|$ can be different for different qubits, e.q., for qubits arranged on
a line the boundary qubits on each end of the system will have a smaller
$|I_{q}|$.

The local validation scheme allows us to identify the regions in the system
where the model performs well or badly in a setup in which possibly all the
device qubits are being actively used. Because it only requires $(t+1)$-local
expectation values, it is highly scalable and has a small computational cost.

The above idea can be pushed further; we can actually use Eq. (4) to _learn_
the local $\mathcal{E}_{k}$ channels. Suppose we have some parameterization of
$\mathcal{E}_{k}$, denoted by $\mathcal{E}_{k,{\bm{\theta}}}$, where
${\bm{\theta}}$ is a vector of parameters. We may define the
${\bm{\theta}}$-dependent cost function

$\displaystyle\Phi({\bm{\theta}})\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\sum_{A}\left(\sum_{k\in\operatorname{supp}(A)}p_{k}\big{(}{\langle{\mathcal{E}_{k,{\bm{\theta}}}^{*}(A)}\rangle}_{\infty}-{\langle{A}\rangle}_{\infty}\big{)}\right)^{2},$
(8)

and find the particular ${\bm{\theta}}$ that best describes the system by
_minimizing_ $\Phi({\bm{\theta}})$. As before, we can write
$\Phi({\bm{\theta}})$ as a quadratic expression of the expectation values of
$(t+1)$-local Paulis in the steady state with ${\bm{\theta}}$-dependent
coefficients:

$\displaystyle\Phi({\bm{\theta}})=\sum_{A}\sum_{\alpha,\beta}C^{(A)}_{\alpha\beta}({\bm{\theta}}){\langle{P_{\alpha}}\rangle}_{\infty}{\langle{P_{\beta}}\rangle}_{\infty}.$
(9)

Estimating these local Pauli expectations at the steady state allows us to
learn the best ${\bm{\theta}}$ that describes the quantum hardware.

The ${\bm{\theta}}$ parameterization can be full, in which case we would need
about $16\times 16$ parameters to describe all possible $2$-qubit channels
$\mathcal{E}_{k}$, or it could contain a much smaller number of parameters if
we use prior assumptions on the structure of $\mathcal{E}_{k}$.

For a system with $n$ qubits, there are exactly
$N_{t}=\sum_{k=1}^{t}3^{k}\left(n-k+1\right)$ $t$-local Pauli operators with
contiguous support, not counting the trivial identity operator. By selecting
$t$-local Pauli operators as the observables $A$ in Eq. (8), we effectively
minimize over $N_{t}$ constraints, and therefore it is desirable that the
number of free parameters will be smaller than $N_{t}$. For example, we may
model the channel that corresponds to an $R_{\bm{x}}(\alpha)$ rotation using
one parameter — $\alpha$, in order to see if the hardware gate under-rotates
or over-rotates with respect to the $\hat{x}$ axis.

###  III.2 A stochastic map: the non strictly-local case

In the previous section we derived our main equation, Eq. (4), under the
assumption that $\mathcal{E}_{k}$ are _strictly local_ channels, acting non
trivially on at most $O(1)$ qubits. This is assumption is no longer true in
realistic noise models, in which also idle qubits, far away from the support
of $\mathcal{E}_{k}$, experience noise such as dephasing and amplitude
damping. In this section we show how our method can be generalized to account
also for (some of) these cases. Specifically, we shall assume a local
Markovian noise model without cross talks, with the following properties:

  1. 1.

Each qubit $j$ has its own idle noise channel $\mathcal{N}_{j}$, which acts
independently of what happens to the other qubits. Thus the noise on the set
of idle qubits $I_{\text{idle}}$ factors into a tensor product of single-qubit
noise channels $\bigotimes_{j\in I_{\text{idle}}}\mathcal{N}_{j}$.

  2. 2.

When a unitary gate or a RESET gate acts on one or two qubits, the actual
channel applied by the QC is a noisy version of these gates that acts only on
the qubits in the support of the ideal gates (i.e., no spillover or cross-
talks to other qubits), _together_ with the idle noise channel on the rest of
the qubits.

Together, these two assumptions imply that every local channel
$\mathcal{E}_{k}$ in the strictly local case is replaced by the non-local
channel
$\tilde{\mathcal{E}}_{k}\otimes\bigotimes_{j\notin\operatorname{supp}(k)}\mathcal{N}_{j}$,
where $\tilde{\mathcal{E}}_{k}$ is the noisy version of $\mathcal{E}_{k}$ and
$\operatorname{supp}(k)$ denotes the qubits in the support of
$\mathcal{E}_{k}$ (which, under our assumptions, also define the support of
$\tilde{\mathcal{E}}_{k}$).

To proceed, we define superoperators $\mathcal{F}_{k}$

$\displaystyle\mathcal{F}_{k}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\Big{(}\bigotimes_{j\in\operatorname{supp}(k)}\mathcal{N}^{-1}_{j}\Big{)}\tilde{\mathcal{E}}_{k}.$
(10)

We will use $\mathcal{F}_{k}$ as a mathematical tool to simplify the equations
we derive below. Note that $\mathcal{F}_{k}$ is not necessarily a channel,
since $\mathcal{N}^{-1}_{j}$ by itself is not necessarily a channel.
Nevertheless it is a _local_ superoperator which satisfies
$\mathcal{F}_{k}(A)=A$ when $A$ is outside its support. The point of using
$\mathcal{F}_{k}$ is that it allows us to write the global action of the $k$th
local channel without the $j\notin\operatorname{supp}(k)$ condition:

$\displaystyle\tilde{\mathcal{E}}_{k}\otimes\bigotimes_{j\notin\operatorname{supp}(k)}\mathcal{N}_{j}=\Big{(}\bigotimes_{j}\mathcal{N}_{j}\Big{)}\cdot\mathcal{F}_{k}.$

Plugging this expression to Eq. (1), we obtain the non strictly-local version
of Eq. (3):

$\displaystyle{\langle{A}\rangle}_{\infty}=\sum_{k}p_{k}\big{\langle}\big{(}\mathcal{F}_{k}^{*}\cdot\bigotimes_{j}\mathcal{N}^{*}_{i}\big{)}A\big{\rangle}_{\infty}.$
(11)

By the same argument used in the previous section, we note that the local
noise channels act on $A$ non-trivially only if they are in its support, and
therefore we can define a “noisy version” of $A$ by

$\displaystyle\Big{(}\bigotimes_{j}\mathcal{N}^{*}_{j}\Big{)}(A)=\Big{(}\bigotimes_{j\in\operatorname{supp}(k)}\mathcal{N}^{*}_{j}\Big{)}(A)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\tilde{A}.$

Plugging this back into Eq. (11), we get the equation

$\displaystyle{\langle{A}\rangle}_{\infty}=\sum_{k}p_{k}{\langle{\mathcal{F}_{k}^{*}(\tilde{A})}\rangle}_{\infty}.$
(12)

This equation is almost identical to its strictly-local counterpart Eq. (3),
except that in the RHS $A$ is replaced by $\tilde{A}$ and
$\mathcal{E}_{k}^{*}$ is replaced by $\mathcal{F}_{k}^{*}$. Following the same
logic as in the previous case, we use the fact that as $\mathcal{F}_{k}$ is
local, ${\langle{F^{*}_{k}(A)}\rangle}_{\infty}={\langle{A}\rangle}_{\infty}$
whenever $A$ is outside its support. Therefore, subtracting
${\langle{\tilde{A}}\rangle}_{\infty}$ from both sides, we obtain the final
equation in which the summation in the RHS is only over channels that
intersect with $A$:

$\displaystyle{\langle{A-\tilde{A}}\rangle}_{\infty}=\sum_{k\in\operatorname{supp}(A)}p_{k}\big{(}{\langle{\mathcal{F}_{k}^{*}(\tilde{A})}\rangle}_{\infty}-{\langle{\tilde{A}}\rangle}_{\infty}\big{)}.$
(13)

As in the previous section, we can use the above set of constraints as a
validation tool for particular noise models by defining local cost functions,
as done in Eq. (7). We can also use it to learn the best noise model from a
family of noise models parameterized by ${\bm{\theta}}$, by using the global
cost function

$\displaystyle\Phi({\bm{\theta}})\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\sum_{A}\Big{(}\sum_{k\in\operatorname{supp}(A)}p_{k}\big{(}{\langle{\mathcal{F}_{k,{\bm{\theta}}}^{*}(\tilde{A}_{\bm{\theta}})}\rangle}_{\infty}-{\langle{\tilde{A}_{\bm{\theta}}}\rangle}_{\infty}\big{)}-{\langle{A-\tilde{A}_{\bm{\theta}}}\rangle}_{\infty}\Big{)}^{2}.$
(14)

As before, to evaluate $\Phi({\bm{\theta}})$ from local measurements, we
expand the operators $\tilde{A}_{\bm{\theta}}$,
$\mathcal{F}_{k,{\bm{\theta}}}^{*}(\tilde{A}_{\bm{\theta}})$ in terms of
$(t+1)$-local Pauli strings with ${\bm{\theta}}$-dependent coefficients and
write the cost function $\Phi(\theta)$ is a quadratic expression of these
expectations (see Eq. (9)).

###  III.3 Deterministic map

Realizing the method presented in Sec. III.2 on a QC requires the execution of
a large number of _different_ quantum circuits, or trajectories, in order to
properly sample the steady state. This makes it challenging to implement on
currently available devices, which are often limited by the number of distinct
circuits that can be executed in a single experimental batch.

To overcome this limitation, we propose an alternative way to construct a non-
unital channel, which relies on a _deterministic_ map that can be implemented
on a QC using a _single_ circuit. For simplicity, we shall describe our
construction in the 1D case with open boundary conditions, but generalization
to other geometries and higher dimensions is straightforward. In such case,
our map $\mathcal{E}$ is defined by a product of local channels
$\\{\mathcal{E}_{k}\\}$, where $\mathcal{E}_{k}$ acts non-trivially on the
neighboring qubits $k,k+1$. The $\mathcal{E}_{k}$ channels are organized in a
two layers brick-wall structure, as shown in Fig. 2.

Figure 2: Example of the deterministic channel implementation on a QC. Left:
Odd channel $\mathcal{E}_{I}$ applied on four qubits, the even layer
$\mathcal{E}_{\text{even}}$ is followed by the odd layer
$\mathcal{E}_{\text{odd}}$. During the action of the $\mathcal{E}_{2}$ on
qubits $2,3$ in the even layer, we assume that the idle noise channels
$\mathcal{N}_{1},\mathcal{N}_{4}$ are acting on the qubits $1,4$. Right: Inner
structure of the 2-local channel $\mathcal{E}_{3}$. The 2-local non-unital
channel $\mathcal{E}_{3}$ consists of applying a two-qubit unitary gate
$U_{3,1}$ on qubits $3,4$, followed by a RESET gate on qubit $3$ along with
idle noise channel $\mathcal{N}_{4}$ on qubit $4$, and ends with applying a
two-qubit unitary gate $U_{3,2}$. The two-qubit unitary gate $U_{3,1}$ has two
known single-qubit rotation gates $R_{\bm{\alpha}},R_{\bm{\beta}}$ about axes
$\hat{\alpha}$ and $\hat{\beta}$ together with a CX gate. The exact details of
our realization of the noisy non-unital channels $\\{\mathcal{E}_{k}\\}$ are
given in Sec. V.3.1.

As explained in Sec. III, to learn the parameters of the quantum channel, we
want its steady state to be non-trivial in the sense that it will not be the
steady state of other local channels. To that aim, we want the local
$\mathcal{E}_{k}$ channels to be non-unital and also entangling so that the
steady state will be a non trivially-entangled state. For example, such
$\mathcal{E}_{k}$ may be realized on a QC by some combination of CX, single-
qubit rotations, and RESET gates, as shown in Fig. 2 and discussed below.

Figure 3: A light cone originating from a 2-local observable $A_{j}$ (red
shaded area) has at most 4-local support. For the case of finite system with
open boundary conditions the light cone is truncated from 4-local to 3-local
at each end of the system.

To derive our constraints, we view $\mathcal{E}$ as the product of two layers.
We define the _odd layer_ by
$\mathcal{E}_{\text{odd}}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\bigotimes_{\text{odd
$k$}}\mathcal{E}_{k}$ and even layer by
$\mathcal{E}_{\text{even}}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\bigotimes_{\text{even
$k$}}\mathcal{E}_{k}$. There are two possible choices for the overall channel
$\mathcal{E}$, based on the order in which the layers are applied:
$\mathcal{E}_{I}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\mathcal{E}_{\text{odd}}\cdot\mathcal{E}_{\text{even}}$
and
$\mathcal{E}_{II}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\mathcal{E}_{\text{even}}\cdot\mathcal{E}_{\text{odd}}$.
An example of the $\mathcal{E}_{I}$ channel is given in Fig. 2. The steady
state of $\mathcal{E}_{I}$ is defined by
$\mathcal{E}_{I}(\rho^{I}_{\infty})=\rho^{I}_{\infty}$, and the steady state
of $\mathcal{E}_{II}$ by
$\mathcal{E}_{II}(\rho^{II}_{\infty})=\rho^{II}_{\infty}$. In general, they
will differ from each other. For brevity, expectation values calculated under
${\rho_{\infty}^{I}},{\rho_{\infty}^{II}}$ are denoted by
${\langle{\cdot}\rangle}_{I},{\langle{\cdot}\rangle}_{II}$ respectively.

Let us now apply Eq. (1) to the case of $\mathcal{E}_{I}$ and a local
observable $A$:

$\displaystyle{\langle{A}\rangle}_{I}={\langle{\mathcal{E}_{I}^{*}(A)}\rangle}_{I}={\langle{\mathcal{E}^{*}_{\text{even}}\big{(}\mathcal{E}_{\text{odd}}^{*}(A)\big{)}}\rangle}_{I}.$
(15)

To simplify the equation, we restrict our attention to observables $A_{j}$
acting on qubits $j,j+1$ for odd $j$. The support $A_{j}$ coincides with the
support of $\mathcal{E}_{j}$ from $\mathcal{E}_{\text{odd}}$. By the same
argument that was used in Sec. III.1, if the supports of $\mathcal{E}_{k}$ and
$A_{j}$ are disjoint, $\mathcal{E}_{k}^{*}(A_{j})=A_{j}$, and so
$\mathcal{E}_{\text{odd}}^{*}(A_{j})=\mathcal{E}^{*}_{j}(A_{j})$. Another
simplification comes from the fact that the light cone of $A_{j}$ has support
only on $\mathcal{E}^{*}_{j-1}$ and $\mathcal{E}^{*}_{j+1}$, as shown in Fig.
3. Therefore,
$\mathcal{E}^{*}_{\text{even}}\big{(}\mathcal{E}^{*}_{\text{odd}}(A_{j})\big{)}=(\mathcal{E}_{j-1}^{*}\otimes\mathcal{E}^{*}_{j+1})\big{(}\mathcal{E}^{*}_{j}(A_{j})\big{)}$.
All together, this allows us to rewrite Eq. (15) in terms of 4-local
expectation values on the RHS:

$\displaystyle{\langle{A_{j}}\rangle}_{I}={\langle{(\mathcal{E}^{*}_{j-1}\otimes\mathcal{E}^{*}_{j+1})\mathcal{E}^{*}_{j}(A_{j})}\rangle}_{I}.$
(16)

For observables $A_{j}$ with even $j$, a similar equation holds with the
steady state of $\mathcal{E}_{II}$ replacing the steady state of
$\mathcal{E}_{I}$.

Following our previous construction with the stochastic map, we can use the
constraints in Eq. (16) to validate a given model for the local channels, or
learn the model that best fits the data from a family of models characterized
by a set of parameters ${\bm{\theta}}$. For the learning task, we start by
defining the cost function $\Phi_{I}({\bm{\theta}})$ using expectation values
in the steady state ${\rho_{\infty}^{I}}$:

$\displaystyle\Phi_{I}({\bm{\theta}})\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\sum_{A_{j}}\left({\langle{[\mathcal{E}^{*}_{j-1,{\bm{\theta}}}\otimes\mathcal{E}^{*}_{j+1,{\bm{\theta}}}]\mathcal{E}^{*}_{j,{\bm{\theta}}}(A_{j})}\rangle}_{I}-{\langle{A_{j}}\rangle}_{I}\right)^{2}.$
(17)

Above, the sum is over a set of $2$-local observables $A_{j}$ defined on
$j,j+1$ for odd $j$, which we take to be Pauli operators. Similarly, we define
$\Phi_{II}({\bm{\theta}})$ from the expectation values at
${\rho_{\infty}^{II}}$, and set

$\displaystyle\Phi({\bm{\theta}})\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\Phi_{I}({\bm{\theta}})+\Phi_{II}({\bm{\theta}}).$
(18)

For a system with $n$ qubits, there are exactly $3n+9(n-1)=12n-9$ such
geometrically $2$-local Pauli operators, excluding the trivial identity
operator. This number provides a rough upper bound to the number of parameters
that can be learned using this method.

To validate a given model, we define local cost functions $\Phi_{q}$, which
for a given $q$ contain all terms in $\Phi_{I},\Phi_{II}$ that overlap $q$.
See Eq. (7) for the equivalent definition in the stochastic case.

As in the case of the stochastic map, we can expand
$\mathcal{E}^{*}_{I}(A),\mathcal{E}^{*}_{II}(A)$ in terms of $4$-local Pauli
operators, and write $\Phi({\bm{\theta}})$ as a quadratic expression of the
expectation values of these operators, as done in Eq. (9).

There are many ways to realize the non-unital channels $\mathcal{E}_{k}$ on a
QC. In this work we use a composed gate, which we call the
$\operatorname{RESU}$ gate ($\operatorname{RESU}$ = $\operatorname{RESET}$ \+
unitary). It is defined by two $2$-local unitaries $U_{k,1}$ and $U_{k,2}$
acting on qubits $k,k+1$ and a $\operatorname{RESET}$ gate on qubit $k$. Each
of the $2$-local unitaries is made of a $\operatorname{CNOT}$ and two general
1-qubit rotations, as shown in Fig. 2. Further details of how we model this
gate in the presence of noise are given in Sec. V.3.1.

In Sec. V.3 we show numerically that our protocol can be used to learn a given
noise model, and in Sec. VI we demonstrate its performance on actual quantum
hardware.

##  IV Noise models and classical optimization

In this section we describe the noise models and the numerical procedures we
used to find ${\bm{\theta}}$ by optimizing over the cost functions in Eqs.
(14, 18),

###  IV.1 Local Markovian noise models

We assume a local Markovian noise model, which can be generally treated using
the Lindbladian master equation [45, 46, 47]. To model the evolution of a
noisy quantum gate, we write

$\displaystyle\frac{d}{dt}\rho=-\frac{i}{T_{0}}[H_{k}({\bm{\theta}}_{k}),\rho]+\frac{1}{T_{0}}\sum_{m}\theta_{m}\mathcal{L}_{m}(\rho).$
(19)

Above, $T_{0}$ is a time scale, to be determined later, and the dimensionless
Hamiltonian $H_{k}({\bm{\theta}}_{k})$ defines the coherent evolution of the
gate $U_{k}$, i.e., $U_{k}=e^{-iH_{k}({\bm{\theta}}_{k})T/T_{0}}$, where $T$
is the gate time and $\\{{\bm{\theta}}_{k}\\}$ are variational gate
parameters. For example, the $\\{{\bm{\theta}}_{k}\\}$ parameters may describe
a source of coherent error in the two-qubit $\operatorname{CX}$ gate.

Under our local noise assumption, we take $\mathcal{L}_{m}$ to be _single-
qubit_ superoperators that represent different single-qubit noise processes.
They are given in terms of jump operators $L_{m}$,
$\mathcal{L}_{m}(\rho)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}L_{m}\rho
L^{\dagger}_{m}-\frac{1}{2}\\{L^{\dagger}_{m}L_{m},\rho\\}$, where
$\\{\cdot,\cdot\\}$ is an anti-commutator. The variational parameters
$\\{\theta_{m}\\}$ model the strength of the different noise processes. For
example, for a single-qubit dephasing noise, we use
$\mathcal{L}_{z}(\rho)=Z\rho Z-\rho$, and the corresponding $\theta/T_{0}$ is
the decoherence rate.

Given a set of channels $\\{\mathcal{E}_{k}\\}$ that define either the
stochastic or the deterministic maps, we take $T_{0}$ to be the maximum over
the running times of $\\{\mathcal{E}_{k}\\}$. Working with the IBMQ hardware,
this time scale was primarily due to the RESET gate time, which is typically
on the order of a microsecond [48, 49]. The $T_{0}$ and gate execution times
we used in our numerical simulations are given in Appendix B and Appendix C.

We allowed the coupling strength (variational) parameters to be different for
each qubit. For a qubit $j$, we use the notation
${\bm{\theta}}^{(j)}=\\{\theta^{(j)}_{m}\\}$ to denote its relevant noise
parameters.

Using the above assumptions, the noise channel
$\mathcal{N}_{{\bm{\theta}}^{(j)}}$ on an idle qubit $j$ for time $T$ is
represented by:

$\displaystyle\mathcal{N}_{{\bm{\theta}}^{(j)}}(\rho)=e^{\frac{T}{T_{0}}\sum_{m}\theta^{(j)}_{m}\mathcal{L}_{m}}(\rho).$
(20)

Similarly, the noisy version of a unitary gate $U_{k}$ described by a
Hamiltonian $H_{k}({\bm{\theta}})$ acting for time $T$ is given by:

$\displaystyle\tilde{\mathcal{E}}_{k,{\bm{\theta}}^{(j)}}(\rho)=e^{\frac{T}{T_{0}}\left(-i[H_{k}({\bm{\theta}}_{k}),\cdot]+\sum_{m}\theta^{(j)}_{m}\mathcal{L}_{m}\right)}(\rho).$
(21)

If $U_{k}$ is a two-qubits gate, the above equation will contain the
$\theta_{m}\mathcal{L}_{m}$ operators of both qubits.

Equation (21) in its entirety was in fact only used to model the noisy
$\operatorname{CX}$ gate when analyzing the experimental data from the IBMQ
hardware. For the single qubit gates, or for the $\operatorname{CX}$ gate in
the numerical simulations, no coherent errors were assumed, and
$H_{k}({\bm{\theta}}_{k})$ was taken to be the ideal Hamiltonian, which is
independent of ${\bm{\theta}}_{k}$. The full details of the
$\operatorname{CX}$ modeling we used are given in Appendix A.

Another simplification was the use of a first-order Trotter-Suzuki
approximation for the case of single-qubit gates. As the typical execution
time of these gates is much shorter than that of the $\operatorname{CX}$ gate
(in the IBMQ hardware their execution time is the range of tens of nanoseconds
[48, 49], which is about an order of magnitude shorter than the execution time
of $\operatorname{CX}$), we approximated Eq. (21) by

$\displaystyle\tilde{\mathcal{E}}_{k,{\bm{\theta}}^{(j)}}\simeq\mathcal{U}_{k}\cdot\mathcal{N}_{{\bm{\theta}}^{(j)}}.$
(22)

Above,
$\mathcal{U}_{k}(\cdot)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}e^{-i\frac{T}{T_{0}}[H_{k},\cdot]}$
is the ideal single-qubit gate channel. Notice that as the resultant error of
this approximation scales as
$O\left((T/T_{0})^{2}\cdot\big{\|}\big{[}[H_{k},\cdot],\sum_{m}\theta^{(j)}_{m}\mathcal{L}_{m}\big{]}\big{\|}\right)$,
and as ${\|H_{k}\|}$ and ${\|\theta_{m}\mathcal{L}_{m}\|}$ are all $O(1)$ (see
Appendix B and Appendix C for more details), we deduce that our Trotter-Suzuki
error for single-qubit gates is of the order of $\big{(}T/T_{0}\big{)}^{2}\sim
10^{-4}$, which is much smaller than all other error sources (e.g.,
statistical error) described below.

We conclude this section by describing how we modeled the noisy
$\operatorname{RESET}$ gate. Here we followed Ref. [32] and modeled the noisy
gate phenomenologically using a Kraus map representation. In this approach,
the noisy $\operatorname{RESET}$ gate on the qubit $j$ is described by two
variational parameters,
${\bm{\theta}}^{(j)}=(\theta^{(j)}_{0},\theta^{(j)}_{1})$. The
$\theta^{(j)}_{0}$ is the probability to measure 0 given that the qubit $j$
was in the state 0, and the $\theta^{(j)}_{1}$ is the probability to measure 1
given that the qubit was in the state 1. All together, the noisy
$\operatorname{RESET}$ gate is modeled by

$\displaystyle\operatorname{R\widetilde{ESE}T}(\rho)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\theta_{0}{|{0}\rangle\\!\,\langle{0}|}\,\rho\,{|{0}\rangle\\!\,\langle{0}|}+\theta_{1}{|{0}\rangle\\!\,\langle{1}|}\,\rho\,{|{1}\rangle\\!\,\langle{0}|}$
(23)
$\displaystyle+(1-\theta_{0}){|{1}\rangle\\!\,\langle{0}|}\,\rho\,{|{0}\rangle\\!\,\langle{1}|}+(1-\theta_{1}){|{1}\rangle\\!\,\langle{1}|}\,\rho\,{|{1}\rangle\\!\,\langle{1}|}.$

###  IV.2 Classical optimization

The minimization of the cost function $\Phi({\bm{\theta}})$ in Eqs. (14, 18)
was done using the stochastic optimization algorithm Adam [50], implemented by
PyTorch with an automatic differentiation engine [51]. We used the same
numerical optimization procedure for both the numerical simulations and the
experimental results on real quantum hardware. This is because the cost
function depends only on the local expectation values, which may come either
from numerical simulations or actual experiments.

Heuristically, we fixed the maximum number of optimization steps to be
$15,000$. We added another heuristic termination criterion for the
optimization process, wherein the optimization ceased if at step $i$, the
difference in the loss function from step $i-500$ was less than $0.25\%$.

In what follows, we denote by ${\bm{\theta}}_{est}$ the channel parameters
which achieve the lowest cost function value of $\Phi({\bm{\theta}}_{est})$.

##  V Numerical simulations

To test our method, we performed several numerical simulations of the
stochastic and deterministic maps, and optimized over the cost function to
learn the underlying noise model. Below, we briefly describe the technical
details of simulations, followed by their results.

###  V.1 Simulation of the dynamics and calculation of the expectation values

For both maps, we simulated systems of $n=5$ to $n=11$ qubits arranged on a
line, with the initial state $\rho_{0}={|{0}\rangle\\!\,\langle{0}|}^{\otimes
n}$. The simulations were done by evolving the full density-matrix on a
computer. An approximation for steady state ${\rho_{\infty}}$ of a quantum
channel $\mathcal{E}$ was obtained numerically by iteratively applying the
quantum channel until the convergence criteria
$D\big{(}\rho_{t},\mathcal{E}(\rho_{t})\big{)}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\frac{1}{2}{\|\rho_{t}-\mathcal{E}(\rho_{t})\|}_{1}<10^{-6}$
was reached. The approximate steady-state was then used to calculate the
expectation values of the local Pauli operators in the different cost
functions $\Phi({\bm{\theta}})$, as given in Eq. (9).

We modeled the statistical noise in the expectation value of the local
observables using a Gaussian random variable that approximates the binomial
distribution of measurement results. For an observable $A$ that is a product
of Pauli matrices, it is easy to verify that the statistical noise of $N$
shots is well approximated by a Gaussian random variable
$S\sim\mathcal{N}(0,\sigma)$ with
$\sigma=\sqrt{(1-{\langle{A}\rangle}^{2}_{\infty})/N}$.

###  V.2 Stochastic map simulations

For the stochastic map simulations, our noise model consisted of dephasing (in
the $\hat{z}$ direction) and amplitude damping, as well as the noisy
$\operatorname{RESET}$ gate from Eq. (23). This amounts to $4$ noise
parameters for each qubit: $2$ for the dephasing and amplitude damping and $2$
for the noisy $\operatorname{RESET}$ gate. In all unitary gates, we assumed no
coherent errors. The exact form of the Lindbladian dissipators is given in
Appendix B.2.

Figure 4: Structure of the $\operatorname{RESCX}(i\to j)$ gate acting two
neighboring qubits $j,i$. We first break the entanglement between $i,j$ by
resetting qubit $j$, and then recreate some entanglement between them by
applying $\operatorname{CX}(i\to j)$. During the action of the
$\operatorname{RESET}$ gate on qubit $j$, the idle noise channel
$\mathcal{N}_{i}$ is acting on qubit $i$.

We simulated 3 different stochastic maps, which used the same local channels
but with different probabilities $\\{p_{k}\\}$. We refer to them as
probability set 1, 2, and 3.

The set of gates that defined the stochastic map consisted of the single-
qubits gates $X^{1/2}$, $H$, $R_{\bm{\alpha}}(\psi)$, and a non-unitary
$2$-local gate between neighboring qubits $i,j$, which we call
$\operatorname{RESCX}(i\to j)$. The $R_{\bm{\alpha}}(\psi)$ gate is a single-
qubit rotation gate by an angle $\psi$ around the axis ${\bm{\alpha}}$, and
$\operatorname{RESCX}(i\to j)$ is discussed below.

Following the IBMQ hardware specifications, the $R_{\bm{\alpha}}(\psi)$ gate
was implemented as a combination of two $X^{1/2}$ and three
$R_{\bm{z}}(\phi_{i})$ gates, which are native gates in IBMQ:

$\displaystyle
R_{\bm{\alpha}}(\psi)=R_{\bm{z}}(\phi_{1})X^{1/2}R_{\bm{z}}(\phi_{2})X^{1/2}R_{\bm{z}}(\phi_{3}),$
(24)

In the formula above, the three rotation angles $\\{\phi_{i}\\}_{i=1}^{3}$
implicitly define the rotation angle $\psi$ and the rotation axis
${\bm{\alpha}}$.

The $\operatorname{RESCX}(i\to j)$ gate is a combination of the
$\operatorname{RESET}$ gate on $j$, followed by a $\operatorname{CX}(i\to j)$,
as shown in Fig. 4. The motivation behind this construction is that the
$\operatorname{RESET}$ gate makes the quantum channel non-unital, but also
breaks entanglement. This might lead to a trivial steady-state, from which it
is hard to learn the underlying noise model. Applying the
$\operatorname{CX}(i\to j)$ gate right after the $\operatorname{RESET}(j)$ was
applied regenerates some entanglement.

The full description of the 3 probability sets, together with the underlying
noise model, the gates and their execution times, are given in Appendix B.

Figure 5: The average diamond distance between the noisy $\operatorname{CX}$
and $\operatorname{RESET}$ gates, for a fixed system size $n=6$ qubits. Each
data point being the average of 10 statistical noise realizations in the
expectation values. The error bars indicate the uncertainty of two standard
deviations. (Top): The average diamond distance for the $\operatorname{CX}$
gates. (Bottom): The average diamond distance for the $\operatorname{RESET}$
gates. Each color represents a different probability set $\\{p_{k}\\}$. The
constant horizontal line (red color) represents the diamond distance between
the ideal case of the noiseless channel and the true noise channel parameters.
Figure 6: The average diamond distance between the noisy $\operatorname{CX}$
and $\operatorname{RESET}$ gates, for a fixed number of shots $N=10^{6}$,
using the probabilities $\\{p_{k}\\}$ from the set 2. Each data point being
the average of 10 statistical noise realizations in the expectation values.
The error bars indicate the uncertainty of two standard deviations.

Our three probability sets were chosen such that the probability of picking a
certain gate is independent of the qubit on which it acts. Taking the time
scale $T_{0}$ (see also Sec. IV.1) to be the maximum running time over the
gate set $\\{X^{1/2},H,R_{\bm{\alpha}}(\psi),\operatorname{RESCX}\\}$, the
implementation of the stochastic channel from Sec. III.2 on a quantum
computer, is described by the following algorithm:

  1. 1.

Pick a gate $1\leq k\leq 4$ from the set
$\\{X^{1/2},H,R_{\bm{\alpha}}(\psi),\operatorname{RESCX}\\}$ with probability
$p_{k}$.

  2. 2.

If the chosen gate is a single-qubit unitary, uniformly choose a qubit and
apply it.

  3. 3.

If the chosen gate is the $\operatorname{RESCX}$ gate, uniformly choose a
neighboring pair of control and target qubits and apply it.

  4. 4.

If the total running time $T$ of the chosen gate from steps 2 or 3 is smaller
than $T_{0}$, wait $\Delta T=T_{0}-T$.

Following the discussion in Sec. III.2, we chose the observables $A$ in the
cost function (14) to be all possible, geometrically 2-local Pauli operators.
Consequently, the optimization procedure used the expectation values of
3-local Paulis.

To evaluate the accuracy of the numerical optimization, we calculated the
diamond distance between the original noisy $\operatorname{CX}$ and
$\operatorname{RESET}$ gates and the ones obtained from the optimization
routine. We focused on these two gates because they have the longest gate
times, and therefore they are the most noisy gates in the simulation.

Fig. 5 and Fig. 6 present the average diamond distance between the noisy
$\operatorname{CX}$ and $\operatorname{RESET}$ gates. For a system of size
$n$, the average is taken over $n$ possible $\operatorname{RESET}$ gates and
$2(n-1)$ possible $\operatorname{CX}$ gates. Fig. 5 shows the average diamond
distance as a function of the number of shots per local Pauli observable, with
the system size fixed at $n=6$ qubits. The results are robust for the three
different probability sets.

With $N=10^{6}$ shots per Pauli observable, the average diamond distance for
the $\operatorname{CX}$ gate is approximately $2\times 10^{-5}$, while the
average diamond distance for the RESET gates is of the order $2\times
10^{-2}$. The difference between these two distances is probably because the
noise parameters that determine the noisy $\operatorname{CX}$ affect not only
the $\operatorname{RESCX}$ gate, but also the idle noise channels and the
noisy single qubit gates. In contrast, the noise parameters of the
$\operatorname{RESET}$ gate appear only when the $\operatorname{RESCX}$ gate
is applied.

The same behavior can be observed in Fig. 6, where we present the average
diamond distance as a function of the number of qubits in the system using
$10^{6}$ measurements per local Pauli operator for the probability set 2.

It would be interesting to find the _optimal_ choice of the gates and
probability distribution $\\{p_{k}\\}$ for learning a given model of noise. We
leave this for future works.

###  V.3 Deterministic maps simulations

####  V.3.1 The RESU gate

As mentioned in Sec. III.3, we realized the $2$-local non-unital channels
$\\{\mathcal{E}_{k}\\}$ of the deterministic maps using a composed gate called
the $\operatorname{RESU}$ gate. As shown in Fig. 2, the $\operatorname{RESU}$
gate on qubits $k,k+1$ consists of two $2$-local unitaries $U_{k,1}$ and
$U_{k,2}$ acting on qubits $k,k+1$ along with a $\operatorname{RESET}$ gate on
qubit $k$:

$\displaystyle\operatorname{RESU}_{k}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}U_{k,2}\cdot(\operatorname{RESET}_{k}\otimes\mathbb{I}_{k+1})\cdot
U_{k,1}.$ (25)

Each of the unitaries $U_{k,i=1,2}$ is made of a CX gate and two arbitrary
single-qubit rotation gates. Therefore, different $\operatorname{RESU}$ gates
are characterized by four single-qubit rotations $R_{\bm{\alpha}}(\psi)$,
which we denote by $\\{{\bm{\alpha}}_{k},\psi_{k}\\}_{k=1}^{4}$. As in Sec.
V.2, each single-qubit rotation $R_{\bm{\alpha}}(\psi)$ is implemented as a
product of two $X^{1/2}$ and three $R_{\bm{z}}(\phi_{i})$ gates, as described
in Eq. (24).

To model the noise in the deterministic map, we used the framework of Sec.
IV.1. For each unitary $U_{k,i=1,2}$ in we defined its noisy version
$\tilde{U}_{k,i=1,2}$ by replacing its noiseless composing gates by their
noisy versions. The noisy $\operatorname{CX}$ gates in $\tilde{U}_{k,i=1,2}$
were modeled by Eq. (21), and the noisy single-qubit rotations
$R_{\bm{\alpha}}(\psi)$ were modeled by Eq. (22). Finally, the noisy
$\operatorname{RESET}$ gate is modeled by Eq. (23). Together, the expression
of the noisy $\operatorname{RESU}$ gate becomes:

$\displaystyle\operatorname{RESU}_{k}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\tilde{U}_{k,2}\cdot(\operatorname{R\widetilde{ESE}T}_{k}\otimes\mathcal{N}_{k+1})\cdot\tilde{U}_{k,1},$
(26)

where $\mathcal{N}_{k+1}$ is an idle noise channel on qubit $k+1$ for the
duration of the noisy $\operatorname{RESET}$ gate.

For the time scale $T_{0}$, as defined in Sec. IV.1, we took the maximum
running time over all $\\{\operatorname{RESU}_{k}\\}$ gates. For any
$\operatorname{RESU}_{k}$ gate with the running time $T_{k}<T_{0}$, a two-
qubit idle noise channel $\mathcal{N}_{k}\otimes\mathcal{N}_{k+1}$ is applied
for the duration $\Delta T_{k}=T_{0}-T_{k}$. On a real QC the idle noise
channel is replaced by the _delay_ instruction.

####  V.3.2 Numerical results

We simulated three deterministic maps, which we refer to as map-1, map-2, and
map-3. Each map used two different $\operatorname{RESU}$ gates: one for the
even layer, and one for the odd layer. In each layer, the same
$\operatorname{RESU}$ gate was used on all qubits. We note that while the
general $\operatorname{RESU}$ gate is defined by four single-qubit rotations
$R_{\bm{\alpha}}(\psi)$, the actual gates used in the simulations and
experiments only used two single qubits rotations in each gate by setting
$R_{\bm{\alpha}}=R_{\bm{\beta}}$ and $R_{\bm{\gamma}}=R_{\bm{\delta}}$ (see
Fig. 2). This was due to a technical reason: limiting the number of parameters
helped us to find sets of angles for which the convergence to the steady state
was rapid (as required by the IBMQ hardware). All together, each deterministic
map was defined by a set of $4$ single-qubit rotations ($2$ for the even layer
$\operatorname{RESU}$ and $2$ for the odd layer $\operatorname{RESU}$). The
full details of the maps, including the rotation axes and angles, as well the
different gate times, are described in Appendix C.

Figure 7: The average diamond distance between the noisy $\operatorname{CX}$
and $\operatorname{RESET}$ gates, for a fixed system size $n=8$ qubits. Each
data point being the average of 10 statistical noise realizations in the
expectation values. The error bars indicate the uncertainty of two standard
deviations. (Top): The average diamond distance for the $\operatorname{CX}$
gates. (Bottom): The average diamond distance for the $\operatorname{RESET}$
gates. Each color represents a different maps defined by the rotation axes and
angles, $\\{{\bm{\alpha}}_{k},\psi_{k}\\}_{k=1}^{4}$. For a specific map, the
pairs $\\{{\bm{\alpha}}_{k},\psi_{k}\\}$ were randomly generated, and kept
constant. The constant horizontal line (red color) represents the diamond
distance between the ideal case of the noiseless channel and the true noise
channel parameters. Figure 8: The average diamond distance between the noisy
$\operatorname{CX}$ and $\operatorname{RESET}$ gates, for a fixed number of
shots $N=10^{6}$, using the rotations $\\{{\bm{\alpha}}_{k},\psi_{k}\\}$ from
the map 2. Each data point being the average of 10 statistical noise
realizations in the expectation values. The error bars indicate the
uncertainty of two standard deviations.

The noise model we used in the simulations of the deterministic maps consisted
of generalized amplitude damping and non-uniform depolarization. For each
qubit this amounted to 5 Lindblad operators: 3 depolarization operators (one
for each axis) and two operators for the generalized amplitude damping [52,
53, 54]. In addition, for each qubit we had two parameters modeling its
$\operatorname{RESET}$ gate according to Eq. (23). As in the stochastic case,
we assume no coherent errors in the unitary gates. The full description of the
noise model is given in Appendix C.2.

After running the simulations, we learned the noise model using the
optimization procedure of Sec. IV.2. To evaluate the accuracy of the noise
parameter estimation, we used the same procedure as in the stochastic maps
simulations, and compared the noisy $\operatorname{CX}$ and
$\operatorname{RESET}$ gates used in the simulations to the ones obtained from
the optimization routine.

In Fig. 7 and Fig. 8, the average diamond distance between the noisy
$\operatorname{CX}$ and $\operatorname{RESET}$ gates used in the simulations
and those obtained from the optimization is presented. For a system of size
$n$, the average is taken over the $n-1$ possible $\operatorname{RESET}$ gates
and $2(n-1)$ possible $\operatorname{CX}$ gates.

In Fig. 7, we present the average diamond distance as a function of the number
of shots per local Pauli observable for a fixed system size of $n=8$ qubits,
using three randomly generated sets of rotation axes and angles,
$\\{{\bm{\alpha}}_{k},\psi_{k}\\}_{k=1}^{4}$. Each of the three sets resulted
in a different steady state. The accuracy of the estimation is robust to the
different sets of rotation angles and axes (which are given in Appendix C.1).
For $N=10^{6}$ shots per local Pauli observable, the average diamond distance
for the $\operatorname{CX}$ gates is approximately $3\times 10^{-4}$, while
for the $\operatorname{RESET}$ gates it is approximately $1.4\times 10^{-2}$.

In Fig. 8 we kept the number of shots per Pauli fixed at $N=10^{6}$ and used
the rotation angles and axes from map 2. Increasing the system size from $n=5$
to $n=11$, we observe that the accuracy of the parameter estimation improves
slightly between $n=5$ and $n=7$, and remains robust up to $n=11$. The initial
increase in accuracy for smaller systems in Fig. 8 may be explained by
examining the light cone shown in Fig. 3 in Sec. III.3. For a smaller system
size with open boundary conditions, the light cone is truncated at each end of
the system, causing the expectation values at the edges of the system to be
3-local instead of 4-local. As a result, the optimization routine uses fewer
statistics for the qubits located at the ends of the system. However, for
larger system sizes, this boundary effect becomes negligible as the relative
fraction of boundary qubits decreases.

As in the case of the stochastic map in Sec. V.2, the average diamond distance
for the $\operatorname{CX}$ gates is much smaller than the
$\operatorname{RESET}$ gates. This is probably because the noise channel
parameters within the $\operatorname{CX}$ gates appear in multiple
combinations within a single $\operatorname{RESU}$ gate, while the channel
parameters of the $\operatorname{RESET}$ gate only appear once per
$\operatorname{RESU}$ gate.

##  VI Quantum device results

In this section we describe the experiments we performed on actual quantum
hardware to test our method. All experiments were done on an IBMQ machine via
the IBMQ cloud.

###  VI.1 Details of the experiment

For the quantum hardware demonstration we chose the _ibm_lagos_ v1.0.32, which
is one of the IBM Quantum Falcon Processors [55] with $7$ qubits in the shape
of a rotated H. Out of the $7$ qubits, we chose $n=5$ qubits arranged in a
line, which we labeled by $0,1,2,3,4$. The IBMQ labels of these qubits were
$6,5,3,1,0$ in corresponding order.

We used the deterministic map-1 in Sec. III.3 to learn the noise of the
quantum hardware. See Appendix C.1 for an exact description of the map. Our
noise model was the same model used in the deterministic map simulations,
which consisted of non-uniform depolarization errors, generalized amplitude
damping, and $\operatorname{RESET}$ errors (total of 7 free parameters per
qubit). In addition to these parameters we added a modeling of coherent
$\operatorname{CX}$ errors as free parameters in its Hamiltonian — see
Appendices A, C.2 for a full description of the model.

After learning the hardware noise model using map-1, we ran deterministic
map-2 and measured $2$-local Pauli expectation values on its steady state. We
used these measurements to test the quality of the noise model we learned from
map-1.

To reduce SPAM errors, we used a readout-error mitigation scheme, which was
partially based on Refs. [56, 57]. The scheme relied on a preliminary step in
which we estimated the conditional probability of measuring a $4$-local bit
string in the computational basis, assuming the qubits were prepared in a
different computational basis. See Appendix D.1 for an full description of our
scheme.

Overall, our experiment consisted of three steps: (i) preliminary readout-
error mitigation measurements, (ii) running deterministic map-1 and measuring
the expectation values of $4$-local Paulis, and (iii) running deterministic
map-2 and measuring the expectation values of $2$-local Paulis.

In all three parts of the experiment we needed to estimate the expectation
values of local Pauli strings. To that aim, we used the overlapping local
tomography method, described by Zubida _et al._ in Ref. [8], which
simultaneously estimates all geometrically $k$-local Pauli strings using a set
of $3^{k}$ different measurement circuits. This means that using a total
budget of $M$ shots, each Pauli string was estimated using $M/3^{k}$ shots.

We now describe the 3 parts of the experiment in more details.

_(i) Readout error mitigation measurements_ : Here we performed preliminary
readout-error mitigation measurements for the protocol that is described in
Appendix D.1. Our measurements were used to mitigate the readout errors of
steps (ii),(iii). We used the overlapping local tomography method (see Zubida
_et al._ in Ref. [8]), to obtain $4$-local readout statistics from a set of
$2^{4}$ different measurement circuits of unit depth. On each circuit we
performed $2\times 10^{5}$ measurements (about $M=3\times 10^{6}$ shots in
total).

_(ii) $4$-local channel measurements_: In this step we ran the deterministic
map-1 on the quantum computer and prepared many copies of the steady states of
$\mathcal{E}_{I}$ and $\mathcal{E}_{II}$ channels. Following the numerical
simulations, we assumed that the steady state of each channel was
approximately reached after 10 steps. We then used the overlapping local
tomography to estimate the expectation values of the $4$-local Pauli operators
using $1.25\times 10^{5}$ measurements for every $4$-local expectation value
(and about $M=2\times 10^{7}$ shots in total).

_(iii) $2$-local channel measurements_: In this step we ran the deterministic
map-2 on the quantum computer and performed $2$-local Pauli measurements using
the overlapping local tomography. As explained in Sec. V.3, map-2 produces a
different steady state from the map-1. We used these $2$-local Pauli
measurements to asses the accuracy of our noise model we learned from the
$4$-local expectation values of step (ii). As in part (ii), we assumed that
the steady states of the $\mathcal{E}_{I}$ and $\mathcal{E}_{II}$ channels
were approximately reached after 10 steps. The $2$-local Pauli expectation
values were obtained using overlapping local tomography that used $1.25\times
10^{5}$ measurements for every $2$-local expectation value (and about
$M=2\times 10^{6}$ shots in total).

Upon finishing the measurements in (ii), (iii), we used the results of readout
measurements (i) to correct these measurement statistics (see Appendix D.1 for
details).

###  VI.2 Results I: Noise model validation

The main application of our method is to _validate_ (benchmark) a noise model
of a QC by calculating its local cost functions $\Phi_{q}$ given a noise model
that defines the noisy channels $\\{\tilde{\mathcal{E}}_{k}\\}$, as explained
in Sec. III.3. The local cost function $\Phi_{q}$ for a qubit $q$ is the sum
of all the $A_{j}$ entries in $\Phi_{I}$ and $\Phi_{II}$ (see Eq. (17)) that
overlap qubit $q$, divided by the number of such $A_{j}$. This happens when
the light cone of $A_{j}$ includes qubit $q$, i.e., when
$q\in\\{j-1,\ldots,j+2\\}$ (see Fig. 3 for an illustration).

Figure 9: The cost map for the _ibm_lagos_ device for three noise models: the
ideal noiseless case (top row), the learned noise model from IBMQ Qiskit
backend (middle row), and the learned noise model using our method (bottom
row).

We used our validation protocol to compare three noise models: (i) the ideal
noiseless model, (ii) our estimate of the IBMQ Qiskit [58] noise model,
imported from the device’s backend properties at the time of our experiment,
and (iii) a noise model that we learned from the IBMQ hardware using our
optimization method (see Sec. IV).

To calculate model (ii), we needed to know
$\mathcal{E}^{*}_{odd},\mathcal{E}^{*}_{even}$ of the IBMQ Qiskit noise model,
from which the local cost function is calculated (see Eqs. (15, 16, 17)).
However, as far as we could tell, IBMQ does not fully publish the exact
formulas it uses to model the local noise. Instead, it publishes general noise
parameters such as $T_{1},T_{2}$, and provides code that simulates noisy
circuits. Therefore, in order to calculate
$\mathcal{E}^{*}_{odd},\mathcal{E}^{*}_{even}$, we first _learned_ IBMQ Qiskit
noise model using the deterministic map-3, and used the resultant model as an
approximation to IBMQ Qiskit noise model. In more details, we used IBMQ Qiskit
simulator to simulate the steady state of map-3, and calculated local Pauli
expectations with vanishing statistical noise. We then used our optimization
routine to learn Qiskit’s noise model. We used the parameters described in
Appendix C.2, but without any coherent errors, and with depolarization only
along the $\hat{z}$ axis. As a sanity check, we used the model we learned to
simulate the 5 qubits steady states of maps-1,2,3 and compared them to
Qiskit’s steady states. We noticed that in all cases the trace distance
between our steady state and Qiskit’s steady state was less than $1\%$.

Fig. 9 shows the single-qubit cost function, as a heat map. The different rows
represent $\log_{10}\Phi_{q}$ for the 3 noise models. Unsurprisingly, it
suggests that IBMQ Qiskit’s noise model provides a better description of the
hardware than the ideal, noiseless channel. It also suggests that our
estimated noise model performs significantly better than the two other models.
This might not be surprising, given the fact that we learned our noise model
using the same steady state that we used to validate it. In other words, the
learned noise model parameters are the ones that minimize the global cost
function, which is roughly the sum of the local cost functions that appear in
Fig. 9. In the machine learning jargon, this might imply that our noise model
estimation procedure overfitted the particular steady state from which the
4-local measurements were sampled. To rule this out (at least to some extent),
we show in the next subsection that our estimated noise model has a better
predictive power than the other two: it is able to better approximate the
2-local expectation values of a steady state of map-2, which is different from
the steady state of map-1 with which we learned the noise model.

###  VI.3 Results II: Noise model characterization

In this section we present the results of applying our optimization algorithm
to _characterize_ the noise within a quantum device. The full learned noise
model parameter details are presented in Appendix D.2.

To evaluate the quality of our noise characterization, we took the noise model
that we learned and used it to numerically simulate the steady states of map-2
(of $\mathcal{E}_{I}$ and $\mathcal{E}_{II}$). On these steady states we
calculated the expectation values of $2$-local Paulis in the bulk of the 5
qubits line (we used qubits $1,2$ for $\mathcal{E}_{I}$ and $2,3$ for
$\mathcal{E}_{II}$). We then compared these simulated expectation values to
the actual results we measured on the quantum hardware when running map-2
(part (iii) of our experiment described in Sec. VI.1). We used the same
procedure to test the two other noise models from the previous section, i.e.,
the ideal (noiseless) model, and IBMQ Qiskit’s model. In contrast to the
previous section, to sample the $2$-local Paulis from the IBMQ Qiskit noise
model we don’t need to know $\mathcal{E}^{*}_{odd},\mathcal{E}^{*}_{even}$.
Hence, we calculated the $2$-local Pauli expectation values directly from the
IBMQ Qiskit simulator using the noise model imported from the device’s backend
properties at the time of our experiment.

Our results are presented as histograms in Fig. 10, where we show the
expectation values of the $9$ possible $2$-local Pauli expectation values for
the steady states of $\mathcal{E}_{I},\mathcal{E}_{II}$ of map-2. The
expectation values from the quantum hardware are shown in red. As we used more
than $10^{5}$ measurements for each data point, the statistical error is
smaller than $0.005$, which is negligible comparing to the differences with
respect to other models. For the 3 analytical models, we preset the exact
expectation values, without any statistical error. Our model prediction is
shown in blue, Qiskit’s model in yellow and the noiseless model in green. As
can be seen in the figure, while the predictions of all the models are in
decent agreement with the empirical results, our model is in better agreement
than the other two models.

Figure 10: Expectation values of the 2-local Pauli operators obtained from:
the _ibm_lagos_ device (red), the learned noise model using our method (blue),
the IBMQ Qiskit simulator (yellow) and the ideal noiseless case (green). We
focus on qubits $1,2$ for $\mathcal{E}_{I}$ (top) and $2,3$ for
$\mathcal{E}_{II}$ (bottom) to avoid including the boundary qubits $0$ and
$4$, for which the measurement statistics become $3$-local instead of
$4$-local (see Sec. V.3). See Table 1 for the corresponding trace distances
between the measured RDM and the RDMs of the three models.

To compare our estimate in a quantitative way, we used the 2-local Paulies
expectation values of the IBMQ hardware to perform quantum state tomography
and estimate the underlying 2-local reduced density matrices (RDMs) in the
steady states of $\mathcal{E}_{I},\mathcal{E}_{II}$ of map-2. We compared it
to the RDMs produced by simulations of our model, the Qiskit model and the
ideal model. For each case, we calculated the trace distance
$D(\rho,\sigma)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\frac{1}{2}{\|\rho-\sigma\|}_{1}$
between the simulation and the estimated RDM of the IBMQ hardware. The results
are given in Table 1, showing that indeed the trace distance between IBMQ RDMs
and our model RDMs are $50\%-70\%$ of the distances from ideal model or the
Qiskit model.

It is worth noting that the 2-local Pauli operators used for the reliability
assessment were collected after the 4-local Pauli expectation values used for
noise model characterization were obtained. As a result, the noise model at
the time of the 2-local Pauli measurements may have undergone some drifting.
This may account for some of the variations seen in the 2-local expectation
values depicted in Fig. 10 [59, 60, 61].

| $\mathcal{E}_{I}$ | $\mathcal{E}_{II}$
---|---|---
Ideal model | 0.19 | 0.21
Qiskit model | 0.174 | 0.19
Learned model | 0.09 | 0.13
Table 1: The trace distance of the 2-local reduced density matrices obtained
from the IBMQ hardware using quantum state tomography and the corresponding
density matrices taken from simulations using the ideal model, Qiskit noise
model, and our learned model. We calculated the distances on the steady states
of $\mathcal{E}_{I}$ and $\mathcal{E}_{II}$ of map-2.

##  VII Summary and outlook

We have presented and demonstrated numerically and experimentally a new
framework for learning a non-unital quanutum channel from its steady state,
which can be used to learn the noise channel of a quantum computer. It is a
natural generalization of the recent results of learning a local Hamiltonian
or Lindbladian from their steady states [24, 25, 26, 27].

We have shown that if the underlying noise channel is Markovian and local,
then it can be learned using our method by driving the system into an
engineered a dissipative steady state. We have shown that there exist a set of
linear constraints between the expectation values of local observables in the
steady state, which can be used to learn or validate the underlying noise
channel. This steady state can be a correlated, entangled state, and might
even be classically inaccessible. Nevertheless, as we showed, the constraints
between its local expectation values are local and can be efficiently verified
once we measure these expectation values. We can therefore apply our method to
a system with a large number of qubits and test how well a given noise model
_globally_ describes the system, which is in a state that might be classically
inaccessible.

For generic dissipative maps, the steady state is unique, and consequently,
our method is independent of the initial state and is insensitive to state
preparation errors. In addition, the fact that we measure a steady state that
was created using several iterations allows us to test possible deviations
from Markovianity. Finally, by optimizing over a parametrized set of local
Markovian noise models, our method allows learning of the optimal noise model.

We have suggested two generic ways to engineer a dissipative map on a quantum
computer, both of which rely on the non-unital $\operatorname{RESET}$ gate.
The first is a stochastic map that is implemented on a quantum computer by
applying a set of gates according to a prescribed probability distribution.
The second, uses a brick-wall type circuit, whose building blocks are composed
$2$-local non-unital gates and can therefore be applied using a single
circuit. We have demonstrated numerically both methods using simulations of
$5-11$ qubits, and showed that it can successfully learn an underlying noise
model. In our numerical tests we managed to learn the underlying
$\operatorname{CNOT}$ and $\operatorname{RESET}$ gates up to diamond distances
of $10^{-3}-10^{-5}$ and $10^{-1}-10^{-2}$ respectively, using $10^{4}-10^{6}$
shots per observable.

We have also tested our method experimentally on on the IBMQ ibm_lagos machine
using $5$ qubits. We have shown that given a noise model, we can construct a
‘heat map’ view of the qubits showing how well the model describes each qubit.
We have also showed how our method can learn an underlying noise model for
that machine, and that the model it learned better predicted the outcomes of
another experiment than the noiseless model or IBMQ Qiskit’s model.

Our work leaves several theoretical and practical questions. On the
theoretical side, it would be interesting to characterize the complexity of
the steady states of dissipative maps like the one we used. Can they be
sampled classically, or can they encode some BQP hard problems? More general
engineered dissipative maps are known to be BQP-complete [37], but our maps
are much simpler. They also differ from the much studied class of randomly
generated circuits, for which there are several complexity bounds (see, for
example, the recent result in Ref. [62] and references within), in several
aspects: our maps are not completely random (we apply the same map over and
over), they are strongly dissipative and non-unital, and can be applied for
$\operatorname{poly}(n)$ time.

It would also be interesting to find a systematical way of engineering the
dissipative maps (of both types) to have an optimal estimate of the error
parameters. This might be possible to do by studying the quantum Fisher
information [63] of the steady states with respect to the free parameters of
the map. Related to that, it would be interesting to understand how the
quantum Fisher information, or more general our ability to use the steady
state for learning, depends on the deviation of the channel from unitarity. Is
there a phase transition behavior, as seen in related measurement induced
phase transitions literature (see Ref. [64] and references within), or is it a
smooth transition?

On the more practical direction, it would be interesting to see if there is a
way to make our method completely insensitive to SPAM errors, by removing its
dependence of measurement errors. The readout error mitigation procedure that
we employed arguably removes a large part of these errors, but not in a
controlled, systematic way, as done in methods like RB [38, 39, 40, 41] and
GST [42, 43, 44].

## Acknowledgements

We thank Eyal Bairey and Netanel Lindner for enlightening discussions. IA
acknowledges the support of the Israel Science Foundation (ISF) under the
Research Grants in Quantum Technologies and Science No. 2074/19, and the joint
NRF-ISF Research Grant No. 3528/20. This research was supported by the
National Research Foundation, Singapore, under its Quantum Engineering
Programme. We are very grateful for the support of the National University of
Singapore (NUS) and the Centre of Quantum Technologies (CQT) for their help in
running the IBMQ experiments. We acknowledge the use of IBM Quantum services
for this work and to advanced services provided by the IBM Quantum Researchers
Program. The views expressed are those of the authors, and do not reflect the
official policy or position of IBM or the IBM Quantum team. This research
project was partially supported by the Helen Diller Quantum Center at the
Technion – Israel Institute of Technology.

## References

  * Burgarth  _et al._ [2009] D. Burgarth, K. Maruyama, and F. Nori, Phys. Rev. A 79, 10.1103/PhysRevA.79.020305 (2009), arXiv:0810.2866 .
  * Di Franco  _et al._ [2009] C. Di Franco, M. Paternostro, and M. S. Kim, Phys. Rev. Lett. 102, 10.1103/PhysRevLett.102.187203 (2009), arXiv:0812.3510 .
  * Shabani  _et al._ [2011] A. Shabani, M. Mohseni, S. Lloyd, R. L. Kosut, and H. Rabitz, Phys. Rev. A 84, 10.1103/PhysRevA.84.012107 (2011), arXiv:1002.1330 .
  * Da Silva  _et al._ [2011] M. P. Da Silva, O. Landon-Cardinal, and D. Poulin, Phys. Rev. Lett. 107, 10.1103/PhysRevLett.107.210404 (2011), arXiv:1104.3835 .
  * Zhang and Sarovar [2014] J. Zhang and M. Sarovar, Phys. Rev. Lett. 113, 80401 (2014).
  * De Clercq  _et al._ [2016] L. E. De Clercq, R. Oswald, C. Flühmann, B. Keitch, D. Kienzler, H. Y. Lo, M. Marinelli, D. Nadlinger, V. Negnevitsky, and J. P. Home, Nat. Comm. 7, 10.1038/ncomms11218 (2016), arXiv:1011.1669v3 .
  * Wang  _et al._ [2018] Y. Wang, D. Dong, B. Qi, J. Zhang, I. R. Petersen, and H. Yonezawa, IEEE Transactions on Automatic Control 63, 1388 (2018), arXiv:1610.08841 .
  * Zubida  _et al._ [2021] A. Zubida, E. Yitzhaki, N. H. Lindner, and E. Bairey, arXiv e-prints , arXiv:2108.08824 (2021), arXiv:2108.08824 [quant-ph] .
  * Kokail  _et al._ [2021] C. Kokail, R. van Bijnen, A. Elben, B. Vermersch, and P. Zoller, Nature Physics 17, 936 (2021).
  * Yu  _et al._ [2022] W. Yu, J. Sun, Z. Han, and X. Yuan, arXiv e-prints , arXiv:2201.00190 (2022), arXiv:2201.00190 [quant-ph] .
  * Wilde  _et al._ [2022] F. Wilde, A. Kshetrimayum, I. Roth, D. Hangleiter, R. Sweke, and J. Eisert, arXiv , 1 (2022), arXiv:2209.14328 .
  * Kappen [2020] H. J. Kappen, Journal of Physics A: Mathematical and Theoretical 53, 214001 (2020).
  * Anshu  _et al._ [2021] A. Anshu, S. Arunachalam, T. Kuwahara, and M. Soleimanifar, Nature Physics 17, 931 (2021).
  * Haah  _et al._ [2021] J. Haah, R. Kothari, and E. Tang, arXiv e-prints , arXiv:2108.04842 (2021).
  * Lifshitz  _et al._ [2021] Y. Y. Lifshitz, E. Bairey, E. Arbel, G. Aleksandrowicz, H. Landa, and I. Arad, arXiv e-prints , arXiv:2112.10418 (2021), arXiv:2112.10418 [quant-ph] .
  * Granade  _et al._ [2012] C. E. Granade, C. Ferrie, N. Wiebe, and D. G. Cory, New Journal of Physics 14, 10.1088/1367-2630/14/10/103013 (2012), arXiv:1207.1655v1 .
  * Wiebe  _et al._ [2014a] N. Wiebe, C. Granade, C. Ferrie, and D. G. Cory, Phys. Rev. Lett. 112, 10.1103/PhysRevLett.112.190501 (2014a), arXiv:1309.0876 .
  * Wiebe  _et al._ [2014b] N. Wiebe, C. Granade, C. Ferrie, and D. Cory, Phys. Rev. A 89, 10.1103/PhysRevA.89.042314 (2014b), arXiv:1311.5269 .
  * Wiebe  _et al._ [2015] N. Wiebe, C. Granade, and D. G. Cory, New Journal of Physics 17, 10.1088/1367-2630/17/2/022005 (2015), arXiv:1409.1524 .
  * Wang  _et al._ [2017] J. Wang, S. Paesani, R. Santagati, S. Knauer, A. A. Gentile, N. Wiebe, M. Petruzzella, J. L. O’brien, J. G. Rarity, A. Laing, and M. G. Thompson, Nat. Phys. 13, 551 (2017), arXiv:1703.05402 .
  * Sone and Cappellaro [2017] A. Sone and P. Cappellaro, Physical Review A 95, 1 (2017), arXiv:1609.09446 .
  * Greiter  _et al._ [2018] M. Greiter, V. Schnells, and R. Thomale, Phys. Rev. B 98, 081113 (2018).
  * Chertkov and Clark [2018] E. Chertkov and B. K. Clark, Phys. Rev. X 8, 031029 (2018).
  * Qi and Ranard [2019] X.-L. Qi and D. Ranard, Quantum 3, 159 (2019).
  * Bairey  _et al._ [2019] E. Bairey, I. Arad, and N. H. Lindner, Phys. Rev. Lett. 122, 020504 (2019).
  * Bairey  _et al._ [2020] E. Bairey, C. Guo, D. Poletti, N. H. Lindner, and I. Arad, New Journal of Physics 22, 10.1088/1367-2630/ab73cd (2020), arXiv:1907.11154 .
  * Evans  _et al._ [2019] T. J. Evans, R. Harper, and S. T. Flammia, arXiv e-prints , arXiv:1912.07636 (2019), arXiv:1912.07636 [quant-ph] .
  * Note [1] Strictly speaking, this is true up to an overall factor.
  * Del Re  _et al._ [2020] L. Del Re, B. Rost, A. F. Kemper, and J. K. Freericks, Phys. Rev. B 102, 125112 (2020).
  * García-Pérez  _et al._ [2020] G. García-Pérez, M. A. C. Rossi, and S. Maniscalco, npj Quantum Information 6, 1 (2020).
  * Sun  _et al._ [2021] S. Sun, L.-C. Shih, and Y.-C. Cheng, arXiv e-prints , arXiv:2106.12882 (2021), arXiv:2106.12882 [quant-ph] .
  * Rost  _et al._ [2021] B. Rost, L. Del Re, N. Earnest, A. F. Kemper, B. Jones, and J. K. Freericks, arXiv e-prints , arXiv:2108.01183 (2021), arXiv:2108.01183 [quant-ph] .
  * Sewell and Jordan [2021] T. J. Sewell and S. P. Jordan, arXiv e-prints , arXiv:2109.09787 (2021), arXiv:2109.09787 [quant-ph] .
  * Tornow  _et al._ [2022] S. Tornow, W. Gehrke, and U. Helmbrecht, Journal of Physics A: Mathematical and Theoretical 55, 245302 (2022).
  * Cattaneo  _et al._ [2022] M. Cattaneo, M. A. C. Rossi, G. García-Pérez, R. Zambrini, and S. Maniscalco, arXiv e-prints , arXiv:2201.11597 (2022), arXiv:2201.11597 [quant-ph] .
  * Harrington  _et al._ [2022] P. M. Harrington, E. J. Mueller, and K. W. Murch, Nature Reviews Physics 4, 660 (2022).
  * Verstraete  _et al._ [2009] F. Verstraete, M. M. Wolf, and J. Ignacio Cirac, Nature Physics 5, 633 (2009).
  * Emerson  _et al._ [2005] J. Emerson, R. Alicki, and K. Życzkowski, Journal of Optics B: Quantum and Semiclassical Optics 7, S347 (2005).
  * Emerson  _et al._ [2007] J. Emerson, M. Silva, O. Moussa, C. Ryan, M. Laforest, J. Baugh, D. G. Cory, and R. Laflamme, Science 317, 1893 (2007), https://www.science.org/doi/pdf/10.1126/science.1145699 .
  * Knill  _et al._ [2008] E. Knill, D. Leibfried, R. Reichle, J. Britton, R. B. Blakestad, J. D. Jost, C. Langer, R. Ozeri, S. Seidelin, and D. J. Wineland, Phys. Rev. A 77, 012307 (2008).
  * Magesan  _et al._ [2012] E. Magesan, J. M. Gambetta, and J. Emerson, Phys. Rev. A 85, 042311 (2012).
  * Blume-Kohout  _et al._ [2013] R. Blume-Kohout, J. K. Gamble, E. Nielsen, J. Mizrahi, J. D. Sterk, and P. Maunz, arXiv 10.48550/ARXIV.1310.4492 (2013).
  * Kim  _et al._ [2015] D. Kim, D. R. Ward, C. B. Simmons, J. K. Gamble, R. Blume-Kohout, E. Nielsen, D. E. Savage, M. G. Lagally, M. Friesen, S. N. Coppersmith, and M. A. Eriksson, Nature Nanotechnology 10, 243 (2015).
  * Nielsen  _et al._ [2021] E. Nielsen, J. K. Gamble, K. Rudinger, T. Scholten, K. Young, and R. Blume-Kohout, Quantum 5, 557 (2021).
  * Head-Marsden  _et al._ [2021] K. Head-Marsden, S. Krastanov, D. A. Mazziotti, and P. Narang, Physical Review Research 3, 13182 (2021), arXiv:2005.00029 .
  * Fratus  _et al._ [2022] K. R. Fratus, K. Bark, N. Vogt, J. Leppakangas, S. Zanker, M. Marthaler, and J.-M. Reiner, arXiv (2022), arXiv:2210.11371 .
  * van den Berg  _et al._ [2022] E. van den Berg, Z. K. Minev, A. Kandala, and K. Temme, arXiv (2022), arXiv:2201.09866 .
  * Kjaergaard  _et al._ [2020] M. Kjaergaard, M. E. Schwartz, J. Braumuller, P. Krantz, J. I.-J. Wang, S. Gustavsson, and W. D. Oliver, Annual Review of Condensed Matter Physics 11, 369 (2020), https://doi.org/10.1146/annurev-conmatphys-031119-050605 .
  * Jurcevic  _et al._ [2021] P. Jurcevic, A. Javadi-Abhari, L. S. Bishop, I. Lauer, D. F. Bogorin, M. Brink, L. Capelluto, O. Günlük, T. Itoko, and N. Kanazawa, Quantum Science and Technology 6, 025020 (2021).
  * Kingma and Ba [2015] D. P. Kingma and J. L. Ba, 3rd International Conference on Learning Representations, ICLR 2015 - Conference Track Proceedings , 1 (2015), arXiv:1412.6980 .
  * Paszke  _et al._ [2019] A. Paszke, S. Gross, F. Massa, A. Lerer, J. Bradbury, G. Chanan, T. Killeen, Z. Lin, N. Gimelshein, L. Antiga, A. Desmaison, A. Kopf, E. Yang, Z. DeVito, M. Raison, A. Tejani, S. Chilamkurthy, B. Steiner, L. Fang, J. Bai, and S. Chintala, in  _Advances in Neural Information Processing Systems 32_ (Curran Associates, Inc., 2019) pp. 8024–8035.
  * Srikanth and Banerjee [2008] R. Srikanth and S. Banerjee, Physical Review A - Atomic, Molecular, and Optical Physics 77, 10.1103/PhysRevA.77.012318 (2008), arXiv:0707.0059 .
  * Cafaro and Van Loock [2014] C. Cafaro and P. Van Loock, Physical Review A - Atomic, Molecular, and Optical Physics 89, 1 (2014), arXiv:1308.4582 .
  * Krantz  _et al._ [2019] P. Krantz, M. Kjaergaard, F. Yan, T. P. Orlando, S. Gustavsson, and W. D. Oliver, Applied Physics Reviews 6, 10.1063/1.5089550 (2019), arXiv:1904.06560 .
  * ref [2022a] Ibm quantum. https://quantum-computing.ibm.com/, 2022 (2022a).
  * Geller [2020] M. R. Geller, Quantum Science and Technology 5, 03LT01 (2020).
  * Geller and Sun [2021] M. R. Geller and M. Sun, Quantum Science and Technology 6, 025009 (2021).
  * ref [2022b] Qiskit: An open-source framework for quantum computing (2022b).
  * Chow  _et al._ [2009] J. M. Chow, J. M. Gambetta, L. Tornberg, J. Koch, L. S. Bishop, A. A. Houck, B. R. Johnson, L. Frunzio, S. M. Girvin, and R. J. Schoelkopf, Physical Review Letters 102, 1 (2009), arXiv:0811.4387 .
  * Fogarty  _et al._ [2015] M. A. Fogarty, M. Veldhorst, R. Harper, C. H. Yang, S. D. Bartlett, S. T. Flammia, and A. S. Dzurak, Physical Review A - Atomic, Molecular, and Optical Physics 92, 1 (2015), arXiv:1502.05119 .
  * Klimov  _et al._ [2018] P. V. Klimov, J. Kelly, Z. Chen, M. Neeley, A. Megrant, B. Burkett, R. Barends, K. Arya, B. Chiaro, and Y. Chen, Physical Review Letters 121, 90502 (2018).
  * Aharonov  _et al._ [2022] D. Aharonov, X. Gao, Z. Landau, Y. Liu, and U. Vazirani, arXiv e-prints , arXiv:2211.03999 (2022), arXiv:2211.03999 [quant-ph] .
  * Meyer [2021] J. J. Meyer, Quantum 5, 539 (2021).
  * Li  _et al._ [2022] Y. Li, Y. Zou, P. Glorioso, E. Altman, and M. P. A. Fisher, arXiv e-prints , arXiv:2209.00609 (2022), arXiv:2209.00609 [quant-ph] .
  * Steffen  _et al._ [2011] M. Steffen, D. P. DiVincenzo, J. M. Chow, T. N. Theis, and M. B. Ketchen, IBM Journal of Research and Development 55, 13:1 (2011).
  * Rigetti and Devoret [2010] C. Rigetti and M. Devoret, Physical Review B - Condensed Matter and Materials Physics 81, 1 (2010).
  * Malekakhlagh  _et al._ [2020] M. Malekakhlagh, E. Magesan, and D. C. McKay, Physical Review A 102, 1 (2020), arXiv:2005.00133 .
  * Magesan and Gambetta [2020] E. Magesan and J. M. Gambetta, Physical Review A 101, 10.1103/PhysRevA.101.052308 (2020), arXiv:1804.04073 .
  * Kandala  _et al._ [2021] A. Kandala, K. X. Wei, S. Srinivasan, E. Magesan, S. Carnevale, G. A. Keefe, D. Klaus, O. Dial, and D. C. McKay, Physical Review Letters 127, 130501 (2021), arXiv:2011.07050 .

##  Appendix A Modeling the noisy CX

We modeled the noisy $\operatorname{CX}$ gate to approximate the actual
$\operatorname{CX}$ gate available on the IBMQ hardware, which is based on the
supercoducting qubits [65]. In IBMQ machines, a two-qubit
$\operatorname{CX}(i\to j)$ gate is realized as a microwave-activated two-
qubit cross-resonance ($\operatorname{CR}$) gate that consists of driving the
control $i$ qubit at the frequency of the target qubit $j$ [66].

Based on the discussions in Refs. [67, 68, 69], we modeled the effective
Hamiltonian $\operatorname{CX}(i\to j)$ by:

$\displaystyle H_{\operatorname{CX}(i\to
j)}\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\theta^{(ij)}_{1}Z_{i}\otimes
X_{j}-\theta^{(ij)}_{2}\mathbb{I}_{i}\otimes
X_{j}-\theta^{(ij)}_{3}Z_{i}\otimes\mathbb{I}_{j},$ (27)

where $\theta^{(ij)}_{m=1,2,3}$ are variational parameters that represent the
strength for each of the interaction terms. The ideal value for all these
parameters is $1$, which leads to a $\operatorname{CX}$ gate when using
dimensionless execution time $T/T_{0}=\pi/4$. In reality, these parameters can
be different from $1$, leading to coherent $\operatorname{CX}$ errors. We
refer the reader to the Refs. [66, 67, 68, 69] for more details on the
$\operatorname{CR}$ gate implementation on transmons.

##  Appendix B Details of the stochastic maps

###  B.1 Gates and probabilities

As discussed in Sec. V.2, our stochastic maps were defined using the gates
$\operatorname{RESCX}$, $X^{1/2}$, $H$, $R_{\bm{\alpha}}(\phi)$. We used 3
stochastic maps, defined by 3 probability distributions. In each map, the
probability of acting with certain gate was independent of the qubit number.
The 3 probability distributions are given in the following table. The
parameters of the single-qubit rotations $R_{\bm{\alpha}}(\phi)$ were chosen
randomly.

Gate | $RESCX$ | $X^{1/2}$ | $H$ | $R_{\bm{\alpha}}(\phi)$
---|---|---|---|---
Probability set 1 | 0.3 | 0.2 | 0.2 | 0.2
Probability set 2 | 0.5 | 0.1 | 0.2 | 0.2
Probability set 3 | 0.7 | 0.1 | 0.1 | 0.1
Table 2: Three different probability distributions $\\{p_{k}\\}$ used in the
numerical simulations. The rotation axis ${\bm{\alpha}}$ and the rotation
angle $\psi$ were defined by three rotation angles $\\{\phi_{i}\\}_{i=1}^{3}$
which we randomly generated.

###  B.2 Noise model in the stochastic case

In the stochastic case, we used 4 parameters for each qubit to model its
noise. We used $\theta_{z},\theta_{amp}$ to model dephasing and amplitude
damping, which are given by the Lindbladians

$\displaystyle\mathcal{L}_{z}(\rho)$
$\displaystyle\stackrel{{\scriptstyle\mathrm{def}}}{{=}}Z\rho Z-\rho,$ (28)
$\displaystyle\mathcal{L}_{amp}(\rho)$
$\displaystyle\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\sigma_{+}\rho\sigma_{-}-\frac{1}{2}\\{\sigma_{-}\sigma_{+},\rho\\},$
(29)

where

$\displaystyle\sigma_{\pm}$
$\displaystyle\stackrel{{\scriptstyle\mathrm{def}}}{{=}}(X\pm iY)/2.$ (30)

In addition, we used $\theta_{0},\theta_{1}$ to model the noisy RESET gate, as
defined in Eq. (23).

We note that in all the simulations of the stochastic maps, we assumed that
the $\operatorname{CX}$ gate had no coherent errors; all its errors came from
he single-qubit errors of its two-qubits.

For every qubit, the dimensionless noise parameters $\theta_{z},\theta_{amp}$
that we used in our simulations were randomly picked between the values of
$0.01$ and $0.07$. In terms of the coherence times (given by
$T_{0}/\\{\theta_{m}\\}$), our choices correspond approximately to a range of
$20\,[\mu sec]$ to $150\,[\mu sec]$.

The parameters for the noisy $\operatorname{RESET}$ were randomly generated
between the values $0.87$ and $0.98$. Following the discussion in Ref. [32],
we required $\theta_{0}\geq\theta_{1}$ to match the typical values on the IBMQ
hardware. However, we note that the experimentally determined noise parameters
of the $\operatorname{RESET}$ gates did not follow this assumption.

###  B.3 Gate durations

The noisy versions of all single-qubit unitary gates were simulated using
Suzuki-Trotter expansion, as outlined in Eq. (22). For the noisy version of
the $\operatorname{CX}$ we used Eq. (21) with the Hamiltonian defined in
Appendix A with vanishing coherent errors.

To define the gate channels in the stochastic map, we used the execution times
that resemble the execution times of actual quantum gates in superconducting
qubits [48, 49]. For the $\operatorname{RESET}$ gates, duration times were
generated between $0.8\,[\mu sec]$ and $1.1\,[\mu sec]$. For the
$\operatorname{CX}$ gates, durations were generated between $0.3\,[\mu sec]$
and $0.45\,[\mu sec]$. The duration time of the $X^{1/2}$ gate was set to
$15\,[nsec]$, and for the $H$ gate the duration was set to $20\,[nsec]$. As
stated in Sec. V.2, the $R_{\bm{\alpha}}(\psi)$ gate is a product of two
$X^{1/2}$ and three $R_{\bm{z}}(\phi_{i})$ gates. The $R_{\bm{z}}(\phi_{i})$
gates are implemented on the superconducting platforms, such as IBMQ, with
zero error and duration [48]. Therefore, we took the duration of the general
single-qubit rotation gate $R_{\bm{\alpha}}(\psi)$, defined in Eq. (24), to be
twice the duration of the $X^{1/2}$ gate.

In our simulations, the time scale $T_{0}$ was dictated by duration of the
$\operatorname{RESCX}$ gate, and was roughly $1.5\,[\mu sec]$.

##  Appendix C Details of the deterministic map

###  C.1 Rotation angles in the RESU gates

We used the deterministic map _both_ in numerical simulations and real quantum
hardware. Consequently, it was important that the maps we defined will
converge to the fixed point in at most 10 steps, which was the maximal number
of steps that we could run of these types of circuits on the IBMQ hardware. To
that aim, we first ran few numerical simulations to find a combination of
rotation angles that will guarantee fast convergence.

Recall from Fig. 2 that the $\operatorname{RESU}$ gate is defined by 4 single-
qubit rotations, which amounts to 12 rotation parameters per
$\operatorname{RESET}$ gate (3 parameters for every single-qubit rotation). To
facilitate the search for optimal rotation angles (which will give rapid
convergence), we reduced this number by half by using the same parameters in
the two single-qubit rotations in $U_{1}$ (the 2-local unitary before the
$\operatorname{RESET}$ gate), and similarly in $U_{2}$ (the 2-local unitary
after the $\operatorname{RESET}$ gate). In other words, we set
$R_{\bm{\alpha}}=R_{\bm{\beta}}$ and $R_{\bm{\delta}}=R_{\bm{\gamma}}$.

The parameters we used for map-1,2,3 are presented in Table 3 and Table 4 for
the even and odd layers, respectively.

| map-1 | map-2 | map-3
---|---|---|---
$U_{1}$, $\phi_{1}$ | -2.63 | 2.13 | 0.77
$U_{1}$, $\phi_{2}$ | 1.71 | -1.31 | -2.46
$U_{1}$, $\phi_{3}$ | 0.21 | -0.71 | 2.59
$U_{2}$, $\phi_{1}$ | -1.07 | -2.07 | -0.79
$U_{2}$, $\phi_{2}$ | -2.12 | -1.12 | 1.09
$U_{2}$, $\phi_{3}$ | 2.73 | -2.73 | -1.80
Table 3: Rotation angles used to define the $\operatorname{RESU}$ gates in the three simulated deterministic maps, for the even layer. The angles are given in radians. | map-1 | map-2 | map-3
---|---|---|---
$U_{1}$, $\phi_{1}$ | 1.42 | -2.32 | 2.57
$U_{1}$, $\phi_{2}$ | -2.92 | 1.92 | 3.04
$U_{1}$, $\phi_{3}$ | -0.33 | -0.73 | 2.67
$U_{2}$, $\phi_{1}$ | 2.83 | -2.83 | -2.72
$U_{2}$, $\phi_{2}$ | 1.36 | 1.87 | -1.79
$U_{2}$, $\phi_{3}$ | -1.08 | -1.98 | 0.21
Table 4: Rotation angles used to define the $\operatorname{RESU}$ gates in the
three simulated deterministic maps, for the odd layer. The angles are given in
radians.

###  C.2 Noise model in the deterministic case

In the deterministic case, we used a more expressive noise model than was used
in the stochastic case, since we wanted the same model to describe also the
real quantum hardware experiments of Sec. VI. For each qubit, we used 7 noise
parameters: $\theta_{x},\theta_{y},\theta_{z}$, to model depolarization in the
$\hat{x},\hat{y},\hat{z}$ directions, $\theta_{amp},\theta_{ex}$ to model
general amplitude damping noise, and $\theta_{0},\theta_{1}$ to model the
noisy RESET gate (as was done in the stochastic case).

The depolarization noise channels in the different directions were modeled by
the Lindbladians

$\displaystyle\mathcal{L}_{\sigma}(\rho)$
$\displaystyle\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\sigma\rho\sigma-\rho,$
$\displaystyle\sigma$ $\displaystyle\in\\{X,Y,Z\\},$ (31)

For the generalized amplitude damping channel [52, 53], we used the
Lindbladian

$\displaystyle\theta_{amp}\big{[}(1-\theta_{ex})\mathcal{L}_{amp,+}+\theta_{ex}\mathcal{L}_{amp,-}\big{]},$

where

$\displaystyle\mathcal{L}_{amp,\pm}(\rho)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}\sigma_{\pm}\rho\sigma_{\mp}-\frac{1}{2}\\{\sigma_{\mp}\sigma_{\pm},\rho\\}.$
(32)

Physically, the $\theta_{amp}$ models the strength of the general amplitude
damping process, and $\theta_{ex}$ models the occupation of the qubit’s state
${|{1}\rangle}$ at thermal equilibrium. The simple amplitude damping model is
described by $\mathcal{L}_{amp,+}$, which corresponds $\theta_{ex}=0$, where
the equilibrium state of the qubit is assumed to be
${|{0}\rangle\\!\,\langle{0}|}$ (which correspond to a qubit at zero
temperature).

In the numerical simulations we chose the dimensionless noise parameters
$\\{\theta_{x},\theta_{y},\theta_{z},\theta_{amp},\theta_{ex}\\}$ randomly
between the values of $0.01$ and $0.1$. The $\theta_{0},\theta_{1}$ parameters
for the noisy $\operatorname{RESET}$ were randomly generated between the
values $0.87$ and $0.98$, as in the stochastic case.

In our numerical simulations we have used several randomly generated gate
durations using the same bounds as for the stochastic map, given in Appendix
B.3. The resulting time scale $T_{0}$ was approximately $2\,[\mu sec]$.

##  Appendix D Details of the quantum device noise model characterization
results

###  D.1 Details of the readout error mitigation scheme for k-local Paulis

For quantum readout error mitigation (QREM), we followed the approach of Ref.
[56], but adapted it to $k$-local readout statistics obtained using the
overlapping local tomography from Ref. [8]. The method empirically evaluates
the conditional probability distribution $P(\underline{x}|\underline{y})$ of
measuring $\underline{x}=(x_{1},\ldots,x_{k})$ classical bit string in
standard basis after the system is prepared the quantum state
${|{\underline{y}}\rangle\\!\,\langle{\underline{y}}|}$ for
$\underline{y}=(y_{1},\ldots,y_{k})$. In our case, we used $k=4$ for learning
the noise channel (part (ii) in Sec. VI.1) and $k=2$ for testing our
characterization (part (iii) in Sec. VI.1). In the ideal case,
$P(\cdot|\cdot)$ is a Kronecker delta-function. Upon obtaining the
$P(\cdot|\cdot)$, we used to ‘fix’ an empirical probability
$P_{noisy}(\underline{x})$ of measuring $\underline{x}$ in the computational
basis, as follows. Assuming that
$P_{noisy}(\underline{x})=\sum_{\underline{y}}P(\underline{x}|\underline{y})\cdot
P_{ideal}(\underline{y})$, we extracted $P_{ideal}(\underline{x})$ by

$\displaystyle
P_{ideal}(\underline{x})=\sum_{\underline{y}}P^{-1}(\underline{x}|\underline{y})\cdot
P_{noisy}(\underline{y}),$

where $P^{-1}(\underline{x}|\underline{y})$ is the matrix inverse of
$P(\underline{x}|\underline{y})$.

###  D.2 Estimated noise parameters of the ibm_lagos v1.0.32 device

The experimental results presented in Sec. VI.2 and Sec. VI.3 were obtained
from the _ibm_lagos_ v1.0.32 quantum device on the IBMQ platform on September
11, 2022. The time scale based on the data from the device backend was
$T_{0}=1.792\,[\mu sec]$.

In Table 5, we give the estimated values of the dimensionless noise (coupling
strength) parameters
$\\{\theta_{x},\theta_{y},\theta_{z},\theta_{amp},\theta_{ex}\\}$, which are
described in Appendix C.2. The values of $10^{-5}$ and $0.1$ are the lower and
the upper bounds on the parameters in the numerical optimization, at which the
gradient with respect to the specific parameter is forced to zero as an
overall stability measure.

| qubit 0 | qubit 1 | qubit 2 | qubit 3 | qubit 4
---|---|---|---|---|---
$\theta_{x}$ | $0.1$ | $10^{-5}$ | $10^{-5}$ | $10^{-5}$ | $10^{-5}$
$\theta_{y}$ | $10^{-5}$ | $0.1$ | $0.1$ | $10^{-5}$ | $10^{-5}$
$\theta_{z}$ | $10^{-5}$ | $0.1$ | $0.1$ | $10^{-5}$ | $10^{-5}$
$\theta_{amp}$ | $0.0598$ | $0.0377$ | $0.1$ | $10^{-5}$ | $0.0749$
$\theta_{ex}$ | $10^{-5}$ | $0.0961$ | $0.0705$ | $0.000279$ | $10^{-5}$
Table 5: Estimated values of the dimensionless noise parameters.

The estimated values of the $\theta_{0},\theta_{1}$ parameters of the noisy
$\operatorname{RESET}$ gates are given in Table 6 (see Eq. (23) in Sec. IV.1).
In the optimization, we used $0.99$ as an upper bound to
$\theta_{0},\theta_{1}$. We note that for the $\operatorname{RESET}$ gate on
the qubit $4$ the optimization was not performed. This is because in the
definition the $\operatorname{RESU}$ gate (see Fig. 2) the
$\operatorname{RESET}$ gate acts on the first qubit in the pair, and so we did
not use $\operatorname{RESET}$ on qubit $4$.

| qubit 0 | qubit 1 | qubit 2 | qubit 3 | qubit 4
---|---|---|---|---|---
$\theta_{0}$ | $0.926$ | $0.941$ | $0.905$ | $0.99$ | $0.96$
$\theta_{1}$ | $0.99$ | $0.927$ | $0.929$ | $0.895$ | $0.96$
Table 6: Estimated values of the $\operatorname{RESET}$ gate parameters.

In Table 7, we give the estimated values of $\operatorname{CX}$ gate noise
parameters, representing the coherent errors in the two-qubit gate (see Eq.
(27) in Appendix A). We note that the $\theta_{3}$ parameter in the gates
$\operatorname{CX}(0\to 1),\operatorname{CX}(1\to 2),\operatorname{CX}(2\to
3),\operatorname{CX}(3\to 4)$ was not affected by the optimization (i.e., the
gradient with respect to $\theta_{3}$ vanished for these gates), and remained
in its ideal inital value of $1.0$. In a hindsight, this can be explained by
the fact that $\theta_{3}$ is equivalent to adding a $\hat{z}$ rotation after
$\operatorname{CX}$. In the gates $\operatorname{CX}(0\to
1),\operatorname{CX}(1\to 2),\operatorname{CX}(2\to 3),\operatorname{CX}(3\to
4)$, this rotation is followed by a $\operatorname{RESET}$ gate, which is
insensitive to phases in the standard basis, and as a result the $\theta_{3}$
parameter does not affect $\operatorname{RESU}$ gate (and therefore it does
not affect the steady state).

| $\theta_{1}$ | $\theta_{2}$ | $\theta_{3}$
---|---|---|---
$\operatorname{CX}(0\to 1)$ | $0.976$ | $1.352$ | $1.0$
$\operatorname{CX}(1\to 0)$ | $1.053$ | $1.029$ | $0.587$
$\operatorname{CX}(1\to 2)$ | $0.791$ | $0.711$ | $1.0$
$\operatorname{CX}(2\to 1)$ | $0.987$ | $1.020$ | $0.936$
$\operatorname{CX}(2\to 3)$ | $1.082$ | $0.973$ | $1.0$
$\operatorname{CX}(3\to 2)$ | $0.970$ | $1.011$ | $0.627$
$\operatorname{CX}(3\to 4)$ | $0.958$ | $0.930$ | $1.0$
$\operatorname{CX}(4\to 3)$ | $1.062$ | $1.008$ | $1.030$
Table 7: Estimated values of the $\operatorname{CX}$ gate parameters. The
notation $\operatorname{CX}(i\to j)$ corresponds to the gate with control
qubit $i$ and target qubit $j$.