Learning the structure of open quantum systems Laura Lewis∗
Ewin Tang∗
John Wright∗
arXiv:2606.30358v1 [quant-ph] 29 Jun 2026
Abstract We design an algorithm for learning the coefficients of an n-qubit constant-local Lindbladian to ε error with O(gd2 log(n)/ε2 ) total evolution time, where g is the single-site energy and d is the (approximate) degree of the interaction graph. Though Lindbladians present new challenges not present in the special case of Hamiltonians, our algorithm achieves the suite of desiderata attained by state-of-the-art Hamiltonian learning algorithms: (1) it uses non-adaptive, ancilla-free randomized Pauli measurement circuits with a time resolution of only Θ(1/g); (2) it works without knowledge of the structure of the unknown Lindbladian; (3) it depends on a smooth form of degree, thereby supporting the learning of quasi-local and power-law Lindbladians. Our algorithm is a simple iterative method, where the objective function consists of Fourier coefficients of the Lindbladian restricted to few-site regions. Its analysis identifies the difficulty unique to open systems, which we call “confusing” terms. For settings where the “confusion” is limited, the performance of the algorithm improves. We demonstrate this for the case of structure learning of Hamiltonians from access to real-time evolution, where we obtain a new algorithm that is significantly simpler than previous work. In addition, using the same iterative method, we design the first efficient algorithm for structure learning Hamiltonians from high-temperature Gibbs states.
Contents 1 Introduction 1.1 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Technical overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Related work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2 2 5 7 9
2 Preliminaries 2.1 Lindbladians . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2 Fourier analysis of quantum channels . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2.1 Local Fourier coefficients . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2.2 Estimating the local Fourier coefficients . . . . . . . . . . . . . . . . . . . . . . . . . .
9 10 11 12 13
3 Local Fourier coefficients of the time evolution operator
15
4 Series expansions 20 4.1 Cluster expansions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 4.2 Bounds on derivatives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 5 Algorithm 26 5.1 Overview of the algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 5.2 The algorithm and guarantee . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 5.3 Proof of correctness . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 5.4 Sample and time complexity analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 5.5 Applications to specific settings . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 5.6 Structure learning Hamiltonians . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 ∗ UC Berkeley. {lllewis,ewin,jswright}@berkeley.edu
1
5.6.1 5.6.2
1
Real-time evolution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . High-temperature Gibbs states . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
41 45
Introduction
Much of modern physics has been built by probing quantum mechanical systems. With the rise of controllable quantum systems as a means to run experiments at large scale [Aru+19; Zho+20; Wu+21; Zhu+23; Aba+25], we must demand a more precise understanding of how to learn the behavior of quantum systems. This drives the field of quantum learning theory, wherein one of the most active topics is Hamiltonian learning from real-time evolution: given some form of access to the dynamics e−iHt of an unknown Hamiltonian, estimate H. However, this assumes that the underlying evolution is performed in a closed system. Far less is known about the analogous question for open systems, which evolve according to a Lindbladian eLt (assuming that the evolution is time-independent and Markovian). This problem is a natural extension of Hamiltonian learning because Hamiltonians form a subclass of Lindbladians. Moreover, this question is particularly timely with the current influx of research interest in Lindbladians, due to recent advances in designing Lindbladians for efficiently preparing Gibbs states [CKBG25; CKG23; DLL24; SA26; RSA26; BLMT26; BLMT24a; BCL24; BC26]. In this work, we present an algorithm for structure learning an unknown geometrically local Lindbladian given access to its real-time dynamics1 . As with Hamiltonian learning, one may hope to optimize many different figures of merit associated to a Lindbladian learning algorithm, including the simplicity of its circuits, total evolution time, time resolution, and generality. Our algorithm is able to achieve performance matching state-of-the-art Hamiltonian learning results in all of these figures of merit. Moreover, to showcase the generality of our framework, a simplification of our main result yields new algorithms for structure learning Hamiltonians from both real-time evolution and high-temperature Gibbs states. To our knowledge, this is the first result for structure learning Hamiltonians from the Gibbs state access model. We first define the Lindbladian learning problem formally. Consider a quantum system consisting of n qubits. The Lindbladian L of this system describes its evolution under Markovian dynamics via the master equation [Lin76; GKS76]: X 1 1X αP [P, ρ] + DP1 ,P2 P1 ρP2 − {P2 P1 , ρ} , (1.1) L(ρ) = 2 2 P ̸=I
P1 ,P2 ̸=I
where P , P1 , and P2 are n-qubit Pauli operators. Throughout, we consider k-local Lindbladians, which means that each term only acts on at most k qubits. Given access to the real-time evolution eLt , the goal is to learn the coefficients of the Lindbladian to error ε > 0 in ∞-norm. We refer to the first sum in Equation (1.1) as the coherent or Hamiltonian part, while the second sum is the dissipative part. The dissipative part models interactions between the system and the environment. In particular, if DP1 ,P2 = 0 for all P1 , P2 , then eLt is simply a Hamiltonian evolution. While evolution under a Hamiltonian e−iHt is always unitary, including dissipative terms in the Lindbladian evolution eLt results in more complex, non-unitary dynamics.
1.1
Results
We focus on learning a physically-motivated class of Lindbladians, where we require two main assumptions. First, we assume that the total interaction strength of terms that act on any given qubit is bounded: ! X X ∥L∥B1 ≜ max |αP | + |DP1 ,P2 | ≤ g. (1.2) i∈[n]
P :supp(P )∋i
P1 ,P2 supp(P1 )∪supp(P2 )∋i
1 Hamiltonian and Lindbladian learning fall into two broad regimes. The first is the (geometrically) local regime, where the terms respect an underlying locality graph, and the total time evolution scales only logarithmically with the system size. The second is the sparse regime, where this condition is not imposed, and in exchange the time evolution scales polynomially with the system size. This work focuses on the first.
2
Such a bound implies that any qubit can only be involved in few terms with large interaction strengths. This norm is sometimes referred to as the one-spin energy [AKL16; Alh23], and an analogous condition has been considered in the Hamiltonian learning literature [BLMT24c]. Second, we require the notion of approximate degree. Recall the typical notion of degree is the maximum number of Lindbladian terms acting on any site. The approximate degree of L can be viewed as a smooth version of the degree: it measures the degree of a modified Lindbladian, obtained by disregarding sufficiently small interactions from L. For a parameter η > 0, it is defined as degη (L) ≜
min
Lbig :∥L−Lbig ∥B1 <η
deg(Lbig ).
(1.3)
To our knowledge, this parameter, although a natural generalization of the degree, has not appeared in any prior Lindbladian or Hamiltonian learning paper. However, similar parameters, e.g. the “effective sparsity” parameter in [BLMT24c], have been considered previously. We defer to Section 2.1 for further discussion of these definitions. With this, we can state our main result, which obtains an algorithm for structure learning local Lindbladians. A more detailed statement can be found in Theorem 5.1, and the full algorithm is presented in Algorithm 2. Theorem 1.1 (Learning local Lindbladians from real-time evolution). Let ε > 0, and let k = O(1). Let L be a k-local Lindbladian with unknown structure and known bounds on the local one-norm ∥L∥B1 ≤ g and approximate degree, degε/(100·16k ) (L) ≤ d. Then there exists a quantum algorithm A that outputs estimates b with the following guarantees: (b α, D) b 1. (Accuracy) With probability 0.99, ∥Lb − L∥B1 ≤ ε, where Lb is the Lindbladian with coefficients α b, D. b This implies that ∥b α − α∥∞ ≤ ε and ∥D − D∥∞ ≤ ε. 2. (Evolution time) A applies eLt with a total evolution time of ttotal = O(gd2 log(n)/ε2 ). 3. (Time resolution) A only applies eLt with t ≥ tmin = Θ(1/g). 4. (Quantum measurements) A performs O(g 2 d2 log(n)/ε2 ) quantum experiments of the following form: (i) prepare a Pauli eigenstate, (ii) apply eLt , (iii) measure in a Pauli eigenbasis. 5. (Classical overhead) A has classical runtime O(nk d log d + (4d)Ck log(dg/ε) + g 2 d2 nk log(n)/ε2 ). We highlight that our learning algorithm only utilizes extremely simple quantum experiments: prepare a Pauli eigenstate, apply time evolution under the unknown Lindbladian, and measure in a Pauli eigenbasis (see Figure 1). We also remark that for g, k, d = O(1), our classical time complexity is O(nk poly(1/ε)). Moreover, the assumptions in our main result are fairly general and encompass a wide range of natural, physically motivated settings (see Section 5.5). Namely, for suitable choices of g and d, we can specialize to the four following cases. We emphasize that all of these results apply to the problem of structure learning: the algorithm does not know a priori which of the unknown Lindbladian coefficients are large. Corollary 1.2 (Geometrically local Lindbladians; Informal version of Corollary 5.6). Let ε > 0. Let L be a k-local Lindbladian with bounded coefficients |αP |, |DP1 ,P2 | ≤ 1 such that each qubit only interacts with O(1) b such that ∥Lb − L∥B ≤ ε with nonzero terms. Then, there exists an algorithm which finds estimates α b, D 1 2 probability at least 0.99 using ttotal = O(log(n)/ε ) total time evolution and time resolution tmin = Θ(1), b where Lb is the Lindbladian with coefficients α b, D. Corollary 1.3 (Local Lindbladians; Informal version of Corollary 5.10). Let ε > 0. Let L be a k-local b such that Lindbladian with ∥L∥B1 = O(1). Then, there exists an algorithm which finds estimates α b, D 2k−2 2 b e ∥L − L∥B1 ≤ ε with probability at least 0.99 using ttotal = O(n log(n)/ε ) total time evolution and time b resolution tmin = Θ(1), where Lb is the Lindbladian with coefficients α b, D. Corollary 1.4 (Informal version of Corollary 5.8). Let ε > 0. Let L be a k-local, quasi-local Lindbladian b such on a p-dimensional lattice. If p, k = O(1), then there exists an algorithm which finds estimates α b, D 2pk 2 that ∥Lb − L∥B1 ≤ ε with probability at least 0.99 using ttotal = O(log(n)(log(1/ε)) /ε ) total time evolution, b where Lb is the Lindbladian with coefficients α b, D. 3
Corollary 1.5 (Informal version of Corollary 5.9). Let ε > 0. Let L be a k-local Lindbladian on a pdimensional lattice with γ-power-law decay for γ − p > 0. Let κ=
2pk . γ−p
(1.4)
b such that ∥Lb − L∥B ≤ ε with probability at least Then, there exists an algorithm which finds estimates α b, D 1 γκ 2+κ b 0.99 using ttotal = O(2 log(n)/ε ) total time evolution, where Lb is the Lindbladian with coefficients α b, D. The last two applications are especially interesting because such Lindbladians with long-range interactions are precisely those used for recent quantum Gibbs state preparation algorithms [CKBG25; DLL24; SA26]. Furthermore, for Lindbladians with power-law decay, as the decay strength γ increases, κ approaches 0, so we recover the scaling ttotal = O(log(n)/ε2 ) in the fast-decay limit. We illustrate the versatility of our framework further by applying it to not only learn different classes of Lindbladians, but also Hamiltonians. In this setting, we consider two different access models: the ability to evolve under the dynamics e−iHt and access to copies of the Gibbs state ρβ ≜ e−βH / tr(e−βH ), where H is the unknown Hamiltonian. In both cases, we obtain new, simple algorithms for structure learning Hamiltonians. Moreover, prior to our work, there were no results in the literature for structure learning Hamiltonians from any Gibbs state access model. Theorem 1.6 (Structure learning Hamiltonians from real-time evolution; Informal version of Theorem 5.13). Let ε > 0. Let H be a k-local Hamiltonian with bounded coefficients |λP | ≤ 1 and ∥λ∥B1 ≤ g. Then, given b − λ∥∞ ≤ ε with probability at b such that ∥λ access to e−iHt , there exists an algorithm which finds estimates λ least 0.99 using ttotal = O(g log(n)/ε2 ) and time resolution tmin = Θ(1/g). Theorem 1.7 (Structure learning Hamiltonians from high-temperature Gibbs states; Informal version of Theorem 5.18). Let H be a k-local Hamiltonian with bounded coefficients |λP | ≤ 1 and each qubit only interacts with O(1) nonzero terms. Let ε > 0, and let β < βc for some critical inverse temperature βc . b such Then, given access to copies of the Gibbs state ρβ , there exists an algorithm which finds estimates λ 2 b that ∥λ − λ∥∞ ≤ ε with probability at least 0.99 using O(log(n)/(βε) ) copies. The classical runtime of this algorithm is O(nk log(n)/(βε)2 ). Prior work on learning local Lindbladians either (1) only has guarantees for Lindbladians with single-qubit dissipative terms [SMDWB24; SMRW25] or (2) has a complexity dependent on the condition number of a large linear system [MECT25; IRGGHY26]. In the first case, [SMDWB24; SMRW25] also both have a time resolution scaling as tmin = O(1/ polylog(1/ε)). In the second case, a priori, this condition number may be exponentially large in n, and [MECT25; IRGGHY26] do not analyze it. In contrast, our work achieves structure learning of local Lindbladians, even for arbitrary dissipative terms and without condition number dependence. We discuss related work in more detail in Section 1.3. We also remark that the guarantees of our Theorem 1.1 are comparable to state-of-the-art Hamiltonian learning results [BLMT24c]. However, our total evolution time has a slightly worse dependence which is quadratic in the approximate degree versus [BLMT24c]’s linear dependence in the analogous sparsity parameter2 . We can also compare Theorem 1.6 to [BLMT24c]. Combining Remarks 3.2 and 5.3 of [BLMT24c] appears to yield the same total time evolution as our Theorem 1.6. However, our algorithm is significantly simpler and has an improved time resolution. We also note that the classical runtime of our algorithm is quasi-polynomial for arbitrary parameters g, d. This inefficiency stems from using the approach of [HKT24b] to compute a truncated series expansion of eLt . We remark that, if one is willing to pay a smaller time resolution and larger total evolution time, then one can improve the time complexity to polynomial. Remark 1.8. If our algorithm instead uses tmin = Θ(1/(g poly(d))) and ttotal = O(g poly(d) log(n)/ε2 ), then e k poly(d)/ε) via the time complexity analysis of [HKT24b]. the classical overhead can be reduced to O(n 2 The Heisenberg limit is not possible to attain for learning Lindbladians, as they define a quantum channel. Thus, the standard quantum limit of 1/ε2 is the optimal dependence on the error parameter.
4
1.2
Technical overview
We explain the ideas behind our algorithm for the well-studied special case of geometrically local Lindbladians. We also restrict to fully dissipative Lindbladians (i.e., αP = 0 for all P ) for simplicity, as the coherent case can be analyzed similarly. For the proof of our more general Theorem 1.1, we refer the reader to Section 5. At a high level, our paper extends the techniques of the prior works [HKT24a; BLMT24c], which were developed for learning Hamiltonians from their time evolutions, and adapts them to the problem of learning Lindbladians. Review of [HKT24b]. We begin by reviewing the approach P of [HKT24b], which performs parameter learning of a geometrically local Hamiltonian H = H(α) = P αP P . The algorithm of [HKT24b] follows bP (α) of expectation values of the form EP (α) ≜ a two-step procedure. First, it produces estimates E tr(O1,P e−iH(α)t O2,P eiH(α)t ) for particular choices of Pauli observables O1,P , O2,P for all P in the known structure of H. Then, it converts these estimates involving the time evolution e−iHt into estimates for the Hamiltonian H by finding zeros of the function FP (x) ≜
1b 1 EP (x) − E P (α), t t
(1.5)
using tools from convex optimization. In particular, it does so via an iterative method, where, at the j-th iteration, the updates to the estimated Hamiltonian parameters are 1 x(j+1) = x(j) − F(x(j) ). 2
(1.6)
Here, we display a simplified version of the Newton-Raphson iteration used in [HKT24b] which obtains the same guarantees3 . The point is that, for the specific choice of observables O1,P , O2,P , FP (x) = 2(xP − αP ) + O(t).
(1.7)
This can be seen by taking a linear approximation to the exponential. Thus, up to a linear order approximation, this is the correct update to apply. The key technical contributions of [HKT24b] are then bounding the contribution of the higher-order terms of F and showing how to efficiently approximate FP (x) via a truncated series expansion. Parameter learning Lindbladians. It is natural to attempt to generalize the algorithm of [HKT24b] to parameter learning for Lindbladians. Note that [HKT24b] does not perform structure learning, so we only focus on parameter learning for now, i.e., the algorithm knows a priori which Lindbladian coefficients are nonzero. Unfortunately, this approach quickly encounters obstacles. Namely, what observables O1 , O2 should we measure? In the Hamiltonian case, there exists a choice of O1,P , O2,P such that EP (α) = 2t · αP + O(t2 ), which is how we obtain Equation (1.7). In other words, for this choice of observables, the expectation values directly estimate the unknown Hamiltonian coefficients. However, even for a single-qubit, fully dissipative Lindbladian, no such choice of observables exists. Instead, expectation values yield a linear combination of Lindbladian coefficients, making it necessary to solve a linear system of equations to recover the coefficients using this approach. This is the source of unanalyzed condition numbers and restrictions to single-qubit dissipators in prior work [SMDWB24; SMRW25; IRGGHY26; MECT25]. One contribution of our work is to overcome this difficulty using tools from Fourier analysis. Inspired by Fourier inversion, we consider expectation values of the form EP1 ,P2 (D) ≜
1 E[tr(eLD t (R)P2 RP1 )], 2n R
(1.8)
3 This version is not present in the literature and was observed by the authors of the current manuscript. The key insight behind this simplification is that [HKT24b] does not leverage the fast quadratic convergence of the Newton-Raphson method. Instead, to obtain their result, they only require linear convergence, which can be attained by many other convex optimization methods, including Richardson iterations, upon which this simplified iteration is based. See Section 5.6 for an analysis of a variant of this iteration.
5
where LD denotes a (fully dissipative) Lindbladian with dissipative coefficients D. When R is uniformly random over n-qubit Paulis, then this quantity can be naturally interpreted as a Fourier coefficient of the channel eLD t . Moreover, these expectations satisfy similar properties as in the Hamiltonian case, where EP1 ,P2 (D) = t · DP1 ,P2 + O(t2 ). Now, one may hope that the guarantees of [HKT24b] apply when instead estimating the expectation values from Equation (1.8). However, when R is sampled uniformly over n-qubit Paulis, it is not possible to parallelize the measurements of these expectation values, resulting in a large total time evolution. Nevertheless, estimating the local Fourier coefficients of eLD , defined as in Equation (1.8) except where R is instead a uniformly random Pauli on supp(P1 ) ∪ supp(P2 ), turn out to be sufficient. Moreover, for a k-local Lindbladian, |supp(P1 ) ∪ supp(P2 )| ≤ k, so these local Fourier coefficients can be estimated efficiently. Luckily, using these local Fourier coefficients, it turns out that the iterative method of [HKT24b] can be extended to obtain a parameter learning algorithm for geometrically local Lindbladians. Structure learning Lindbladians. We now discuss pushing these ideas further to the problem of structure learning Lindbladians. While [BLMT24c] can be viewed as a way to extend [HKT24b] to structure learning for Hamiltonians, performing similar modifications to our parameter learning algorithm for Lindbladians is not straightforward. In particular, [BLMT24c] uses techniques which are specialized to the Hamiltonian b to the setting. Namely, it begins by running a base learning algorithm to produce a coarse approximation H b true Hamiltonian H. Then, to improve its estimate, it uses Trotterization to simulate access to e−i(H−H) b b and obtain an estimate of H − H, which can in turn be added back to H to produce a better estimate of H. However, Trotterization requires access to the inverse time evolution, and so this cannot be done for Lindbladians, because they are dissipative and their evolutions cannot be reversed. Thus, we attempt a different modification for structure learning. A simple approach one may take is to keep track of all k-local coefficients instead of only the coefficients in the known structure. In other words, an algorithm may update coefficient estimates using the errors FP1 ,P2 (x) ≜
1 1b EP1 ,P2 (x) − E P ,P (D) t t 1 2
(1.9)
for all k-local P1 , P2 . The problem is that the recovery of the Fourier coefficients of eLD t from the local Fourier coefficients then becomes more difficult, as many Lindbladian terms can interfere with each other. Moreover, the contribution of Lindbladian terms which are small but nevertheless still part of the structure can be obscured by the contribution of large terms. We quantify this “confusion” as follows. Because the expectation over R in Equation (1.8) is not taken over all n-qubit Paulis, the local Fourier coefficient EP1 ,P2 (D) does not precisely approximate the dissipative coefficient DP1 ,P2 . Instead, the local Fourier coefficients also include contributions from other Paulis (Q1 , Q2 ), which are “confused” with the correct term (P1 , P2 ). We write (Q1 , Q2 ) ⪰ (P1 , P2 ) to denote such Paulis, which are defined as (Q1 , Q2 ) that agree with (P1 , P2 ), respectively, on supp(P1 ) ∪ supp(P2 ) and that agree with each other outside of this support. In Section 2.2.1, we prove that X EP1 ,P2 (D) ≈ t DQ1 ,Q2 ≜ t(AD)P1 ,P2 , (1.10) (Q1 ,Q2 )⪰(P1 ,P2 )
up to a linear approximation of eLD t . The problem described above, that recovering the Lindbladian coefficients from the local Fourier coefficients becomes difficult, can be made precise in that the matrix A does not have a well-behaved inverse. In particular, ∥A−1 ∥∞→∞ can scale polynomially in n. Operationally, this means that if one attempts to perform an iteration of the form x(j+1) = x(j) − A−1 F(x(j) ),
(1.11)
which is a natural extension of Equation (1.6), the error in each iteration increases by a factor of poly(n). To counteract this, one would need to estimate the local Fourier coefficients to ε/ poly(n) error. However, this results in an abysmal total time evolution of O(poly(n)/ε2 ) for learning geometrically local Lindbladians, whereas one would typically expect an exponentially smaller total time evolution of O(log(n)/ε2 ). This is an obstacle unique to Lindbladian learning. In contrast, for Hamiltonian learning, there is no “confusion”: the local Fourier coefficients still approximate the corresponding Hamiltonian coefficient directly. 6
Thus, the matrix A in Equation (1.10) is simply the identity matrix, which has a bounded ∞ → ∞ norm. We refer the reader to Section 5.6 to see how this greatly simplifies the analysis. The critical problem here is that to estimate DP1 ,P2 , all possible pairs (Q1 , Q2 ) which can be confused with (P1 , P2 ), of which there can be roughly O(nk ), contribute some error, whereas we should only actually have O(d) large terms which contribute large error. To remedy this, we round small entries of both F(x(j) ) and our estimate after an update to zero, i.e., we consider the iteration x(j+1) = Roundτj (x(j) − A−1 Roundτj′ (F(x(j) ))),
(1.12)
for some carefully chosen thresholds τj , τj′ > 0. After rounding, the remaining nonzero coefficients correspond to the structure of L discovered so far. Rounding in this way ensures that our estimated Lindbladian in each iteration always has degree O(d), so we effectively only incur a total time evolution cost comparable to parameter learning a Lindbladian with degree O(d). Technically, analyzing this new rounded algorithm requires bounding ∥A−1 ∥B1 →B1 (instead of ∥A−1 ∥∞→∞ , see Section 3) and maintaining the error of our iterates in B1 -norm. Throughout this discussion, we have also been ignoring errors arising from linear approximation, and a significant portion of our analysis is dedicated to showing that this error is not too large (see Section 4).
1.3
Related work
Hamiltonian learning. For the simpler task of Hamiltonian learning, there is extensive literature for solving this task in a variety of access models, e.g., from copies of the Gibbs state [BAL19; AAKS21; HKT24b; BLMT24b; CAN25; QR19], access to the Hamiltonian’s real-time evolution [ZYLB21; HKT24b; HTFS23; BLMT24c; Zha24; MFPT24; Hu+25; ACGGMS25; DOS24; Car24; OKSM24; Gut24; ADG24; CW23; CLS25; ST25; BCGOR25; SLO26; MH24; MBGTWR25; LTGNY24; NLY24], and more restrictive settings [BZMT26; CCH26; CCH25]. The works most relevant to ours are those that consider Hamiltonian learning from access to dynamics, where the unknown Hamiltonian is promised to be (geometrically) local. We only detail the results of some of these works and refer to, e.g., [BLMT24c], for a more thorough review. In this setting, early works, e.g., [BAL19], designed an algorithm using ttotal = O(log(n)/ε3 ) and tmin = O(ε). This approach can be modified to achieve structure learning. For known structure, [HKT24b] improved this to ttotal = O(log(n)/ε2 ) and √ tmin = Ω(1). Moreover, [HTFS23] achieved the Heisenberg scaling with ttotal = O(log(n)/ε) and tmin = Ω( ε). Both [HKT24b; HTFS23] only work for parameter learning, not structure learning. [BLMT24c] later achieved Heisenberg scaling while maintaining a constant time resolution, i.e., ttotal = O(log(n)/ε) and tmin = Θ(1). Our work is most comparable to [HKT24b; BLMT24c], as we achieve the optimal scaling of ttotal = O(log(n)/ε2 ) with a constant time resolution tmin = Θ(1). At a high level, our algorithm resembles those of [HKT24b; BLMT24c], as ours is also an iterative procedure based on convex optimization algorithms. While [BLMT24c] can be seen as a way to extend [HKT24b] to learn the structure of Hamiltonians, performing a similar modification for Lindbladians is not straightforward. In particular, in each iteration, [BLMT24c] b where H b is the current updates the learned parameters based on expectations with respect to exp(−i(H − H)), b can be simulated via a new constant-time Trotterization estimate of the Hamiltonian, and access to H − H formula. However, a similar Trotterization formula is not expected to be possible for Lindbladians. Thus, one main conceptual contribution of our work is to overcome this barrier and find a new way to update our estimates of the Lindbladian parameters. The above discussion takes all parameters to be constant, i.e., k, g, d = O(1). When considering the scaling with g and d, our algorithm for Lindbladian learning achieves ttotal = O(gd2 log(n)/ε2 ). Meanwhile, [BLMT24c] has a better scaling of ttotal = O(d log(n)/ε)4 . While the Heisenberg scaling 1/ε is not possible for learning Lindbladians, the total time evolution of [BLMT24c] is still better by a factor of d and g. For our Hamiltonian learning result in Theorem 1.6, [BLMT24c] implicitly appears to obtain the same total time evolution (seen by combining their Remarks 3.2 and 5.3). However, our algorithm is significantly simpler and has an improved time resolution (our algorithm has tmin = Θ(1/g) compared to their Θ(1/d)). 4 [BLMT24b] considers a slightly different parameter than our approximate degree d, which they call the effective sparsity. These capture similar physical settings, so we state their complexity in terms of d for comparison.
7
Lindbladian learning. Early works studied the task of recovering a description of the Lindbladian from access to its dynamics or steady states but lacked rigorous guarantees [Buž98; BGPLA20]. Since then, interest in this problem has gained momentum with the development of many heuristic/numerical algorithms [LSKOC25; OKC23; OKKKZ25; POKZ22] and even some experimental demonstrations [Kra+25; Bir+26; LJFV26; BMWM25]. The works most relevant to the present manuscript are those with provable guarantees on the total time evolution required to learn the Lindbladian given access to its time evolution operator. However, no prior work achieves a rigorous guarantee for learning local Lindbladians with optimal performance in total evolution time and time resolution, even for the easier task of parameter learning. [SLP11] uses time derivative estimation, resulting in ttotal = O(ν log(n)/ε3 ) (when combined with randomized measurements as in [HKT24b]) and tmin = O(ε). Here, ν is a condition number factor, which is implicit in the complexity, as their algorithm requires inverting a linear system which is, a priori, not well-conditioned. In addition, recent works studying this problem fall into two main categories: 1. The work gives algorithms for structure learning, but with guarantees which only hold in very restricted settings. 2. The work gives algorithms for parameter learning. Moreover, the total time evolution ttotal depends on the condition number of a large linear system, which is not analyzed and could be exponentially large in n. In particular, [SMDWB24; SMRW25] both fall into the first category, where they achieve structure learning 2 e of local Lindbladians with ttotal = O(log(n)/ε ) and tmin = O(1/ polylog(1/ε)), but their guarantees only apply to Lindbladians with single-qubit dissipative terms. Their algorithm’s dependence on d is not clear, as they take d = O(1) throughout. [SMDWB24] also considers Lindbladians with single-qubit dissipators and algebraic decay on a p-dimensional lattice. In this case, their algorithm achieves a similar complexity as ours. However, for a decay rate of γ > 0, their guarantee only holds for γ ≥ 5p. [MECT25; IRGGHY26] belong to the second category5 . [MECT25; IRGGHY26] achieve ttotal = e 2 ν log(n)/ε2 ) and tmin = Θ(1/d) albeit only for parameter learning. Moreover, their complexities O(d hide the cost of solving an a priori ill-conditioned linear system, which is quantified here via the condition number factor ν. It is also worth noting that for general k-local Lindbladians, [IRGGHY26] presents a 4k 4 e structure learning algorithm with ttotal = O(νn /ε ), but this complexity still hides a condition number factor. In contrast, our work achieves structure learning of local Lindbladians with ttotal = O(gd2 log(n)/ε2 ) and tmin = Θ(1/g), even for arbitrary dissipative terms, and without condition number dependence. Moreover, for Lindbladians with algebraic decay, our guarantee holds for decay rates γ > p, beyond which there is a natural barrier [DDMPRT23]. Our algorithm also applies to general k-local Lindbladians, where we achieve ttotal = O(gn2k−2 log(n)/ε2 ). Concurrent work. While preparing this manuscript, we became aware of several independent and concurrent works [Rom+26; ACGRY26; Sin26; MBFR26] which study the problem of learning Lindbladians. Two of these works [Rom+26; Sin26] operate in the sparse setting and are not comparable to our work. On the one hand, they obtain structure learning guarantees for a broader class of Lindbladians without sparsity assumptions, but, on the other hand, their algorithms utilize a total evolution time which scales polynomially in the system size and require at least n ancillary qubits. In contrast, our algorithm’s total evolution time scales logarithmically in system size, and we use zero ancillas. [ACGRY26] addresses the same setting as our work, but their algorithm requires worse parameter 2 dependencies. In particular, for k = O(1), [ACGRY26] uses a total evolution time of ttotal = Õ(Λd2k dis log(n)/ε ) and time resolution tmin = Θ(1/Λ) for Lindbladians with “dissipative site degree” ddis and “local dynamical strength” Λ. Their parameters of ddis and Λ are comparable to our parameters d and g, respectively. Moreover, [ACGRY26] does not consider learning Lindbladians with decaying long-range interactions. In contrast, our algorithm achieves ttotal = O(gd2 log(n)/ε2 ), which is fixed-parameter tractable, and a similar time resolution 5We only quote the result from [IRGGHY26] relevant to our paper, which is parameter learning for local Lindbladians. In the sparse setting, they obtain a structure learning algorithm.
8
of tmin = Θ(1/g). We also obtain guarantees for learning Lindbladians with exponentially decaying and power-law decaying interactions. [MBFR26] also considers the local setting, and their algorithm has similar scalings as [ACGRY26]. Namely, for k, g = O(1), their result uses ttotal = Õ(d2k log(n)/ε2 ). It is not clear how their algorithm scales with g. Also, for Lindbladians with algebraic decaying interactions, [MBFR26] does not perform structure learning: the set of terms with large interaction strengths are provided as input to the algorithm. In comparison, our result is fixed-parameter tractable, achieving a significantly better d dependence. Even with Remark 1.8 for improving our classical runtime, our dependence on d is still poly(d), rather than dk . Both of our results also hold for learning general k-local Lindbladians and achieve similar complexities. For long-range interactions, our results are strictly stronger than [MBFR26], as we perform structure learning and are not given the set of large terms. The above comparisons are made with respect to the quantum resources required, i.e., the total time evolution and time resolution. However, for classical time complexity, our algorithm is only quasi-polynomial e 2 nk d2k ), which for arbitrary parameters g, d. Meanwhile, [ACGRY26] has a classical time complexity of O(Λ dis k is polynomial for k = O(1). Also, [MBFR26] has a classical time complexity of O(n ).
1.4
Discussion
In this work, we develop a general framework for structure learning an unknown Lindbladian given access to its dynamics. We show that our result can be instantiated in several physically motivated settings, including geometrically local, general k-local, quasi-local, and power-law Lindbladians. We can also specialize our proof to apply to two problems in Hamiltonian learning. Namely, we give a new algorithm for structure learning Hamiltonians from real-time evolution, where we obtain a surprising total time evolution which is independent of the approximate degree or effective sparsity. In addition, we design the first algorithm for structure learning Hamiltonians from high-temperature Gibbs states. There are still several interesting open questions to explore. 1. What is the optimal total time evolution scaling for learning Lindbladians? Even for parameter learning Lindbladians, our algorithm achieves a total time evolution of ttotal = O(gd2 log(n)/ε2 ). Meanwhile, the state-of-the-art Hamiltonian learning results [BLMT24c] are able to attain O(r log(n)/ε2 ) for an effective sparsity parameter r, which is analogous to our approximate degree d. Is the d2 dependence for learning Lindbladians fundamental? 2. What is the optimal complexity for structure learning Hamiltonians? By adapting our framework, we show a new guarantee of ttotal = O(g log(n)/ε2 ) and tmin = Θ(1/g) for Hamiltonian learning. Is the dependence on effective sparsity r in [BLMT24c] required to obtain the Heisenberg scaling? Could this be related to the lack of inverse access in this model [TW25]? 3. Our techniques for structure learning Hamiltonians from high-temperature Gibbs states do not extend to the low-temperature regime because the cluster expansion diverges at low temperatures. Can one efficiently learn the structure of Hamiltonians from Gibbs states at any temperature? [BLMT24b] achieves the analogous result in the parameter learning case.
2
Preliminaries
e ) = O(f polylog(f )). We use the Iverson bracket: JGK = 1 if the Throughout, [n] = {1, . . . , n}, and O(f proposition G is true and 0 otherwise. We denote the complement of a set S by S ⊆ [n], where S ≜ [n] \ S. For a matrix M , we use M † to denote its conjugate transpose. Given τ ≥ 0, we define ( ( 0 |x| ≤ τ x |x| ≤ τ Roundτ (x) = , and Truncτ (x) = . (2.1) x otherwise 0 otherwise Throughout, we also let n denote the number of qubits, N ≜ 2n , and tr ≜ tr /N . We let Pn ≜ {I, X, Y, Z}⊗n denote the set of 4n n-qubit Pauli matrices. For a Pauli P ∈ Pn , we denote SP ≜ supp(P ), where supp(P ) ≜ {i ∈ [n] : Pi ̸= I}. For a pair of Paulis P1 , P2 ∈ Pn , we sometimes write P ≜ (P1 , P2 ) for 9
brevity. Moreover, similarly to the single Pauli case, we write SP ≜ supp(P1 ) ∪ supp(P2 ). We also write sP ≜ |SP | and sP ≜ |SP |. Finally, consider a matrix M = a · P , where a ∈ {±1, ±i} and P ∈ Pn ; such matrices arise from products of Pauli matrices. Then we will write c(M ) = a and P(M ) = P , so that c(M ) · P(M ) = M.
(2.2)
For a vector x, we use ∥x∥∞ ≜ maxi |xi | to denote the ∞-norm. We use ∥x∥22 ≜ i |xi |2 to denote x∥2 the 2-norm. For a matrix M , we use ∥M ∥op ≜ maxx̸=0 ∥M ∥x∥2 to denote the spectral norm. We also write √ ∥M ∥tr ≜ tr( M † M ) to denote the trace norm. For a matrix M ∈ CN ×N , we define the operator norm corresponding to norms ∥·∥a and ∥·∥b on CN as P
∥M ∥a→b ≜ sup ∥M x∥b .
(2.3)
∥x∥a ≤1
Note that ∥M x∥b ≤ ∥M ∥a→b ∥x∥a , and ∥LM ∥a→c ≤ ∥L∥b→c ∥M ∥a→b . Moreover, we define a superoperator as a linear map from CN ×N to CN ×N .
2.1
Lindbladians
First, we define a Lindbladian, which defines the Markovian dynamics of an open quantum system via the Lindblad master equation. Definition 2.1 (Lindbladian). A Lindbladian is a linear map L : CN ×N → CN ×N that, applied to a quantum state ρ ∈ CN ×N , can be written as X 1X 1 L(ρ) = αP [P, ρ] + DP1 ,P2 P1 ρP2 − {P2 P1 , ρ} , (2.4) 2 2 P ̸=I
P1 ,P2 ̸=I
where D is a positive semidefinite matrix and α is purely imaginary. This Lindbladian is k-local if every P satisfies |supp(P )| ≤ k and every P1 , P2 satisfies |supp(P1 ) ∪ supp(P2 )| ≤ k. Throughout, we assume that we are working with nonzero Lindbladians. In particular, we always assume that k and g are nonzero. We often find it convenient in our analysis to consider combining the coherent and dissipative coefficients into a single vector, i.e., writing X 1 L(ρ) = λP1 ,P2 P1 ρP2 − {P2 P1 , ρ} , (2.5) 2 P1 ,P2 P1 ̸=I
where now P2 is allowed to be the identity. This simply allows us to index into the vector λ with a pair of Paulis instead of one, making our notation simpler throughout. These are clearly equivalent definitions,6 as one can set λP,I = αP and λP1 ,P2 = DP1 ,P2 . To make the dependence on the coefficients explicit, we sometimes write Lα,D or Lλ . When the subscript is a single vector, one should consider the representation in Equation (2.5). Definition 2.2 (Parameterized Lindbladian). For a Lindbladian and a coefficient vector (α, D) ∈ Cm , we let Lα,D refer to the corresponding Lindbladian in Equation (2.4). Moreover, we denote by Lλ the superoperator of the form Equation (2.5), which we call a parameterized Lindbladian. For notational simplicity, we will sometimes index into λ, not by an explicit Lindbladian term P1 , P2 (where P2 is arbitrary and P1 = ̸ I), but by an index a ∈ [m]. The term associated to a will then be denoted P a . Though we sometimes refer to Lλ as a Lindbladian, for general λ, this may not define a valid physical Lindbladian and may only refer to a superoperator mapping CN ×N → CN ×N . We also often use the adjoint of the Lindbladian. 6 The reason we include the factor of 1/2 in the definition in Equation (2.4) is so that λ
10
P,I = αP , instead of 2αP .
Definition 2.3 (Adjoint of a Lindbladian). The adjoint L† of a Lindbladian L is defined such that tr(XL(Y )) = tr(L† (X)Y ) for any operators X, Y . Explicitly, using the Pauli expansion above, one can write the adjoint as X 1X 1 L† (O) = − αP [P, O] + DP1 ,P2 P2 OP1 − {P2 P1 , O} , (2.6) 2 2 P ̸=I
P1 ,P2 ̸=I
where O ∈ CN ×N is an observable. The analogous definition for Lλ is clear. Similarly to [BLMT24c], our results depend on a “local norm,” which is defined as follows. Definition 2.4 (Local norm of a Lindbladian). Let Lλ be a superoperator. Then, we define the B1 -norm of λ as X ∥λ∥B1 ≜ max |λP |. (2.7) i∈[n]
P :SP ∋i
We sometimes write ∥Lλ ∥B1 = ∥λ∥B1 . Moreover, note that ∥Lλ ∥B1 = ∥Lα,D ∥B1 . Definition 2.5 (Degree and approximate degree of a Lindbladian). The degree deg(Lλ ) of a superoperator Lλ is the maximum number of terms supported on a site, i.e., deg(Lλ ) ≜ max |{P : i ∈ SP , λP ̸= 0}|
(2.8)
i∈[n]
The approximate degree degε (Lλ ) of a superoperator Lλ is the minimum deg(Lbig ) over ways to split Lλ = Lbig + Lsmall such that ∥Lsmall ∥B1 < ε, i.e., degε (Lλ ) ≜
min
Lbig :∥Lλ −Lbig ∥B1 <ε
deg(Lbig ).
(2.9)
Here, Lsmall = L − Lbig . Sometimes, we also write deg(λ) = deg(Lλ ) and degε (λ) = degε (Lλ ). Consider the optimal splitting Lλ = Lbig + Lsmall , and let λbig and λsmall denote the coefficients of Lbig and Lsmall , respectively. We observe that, without loss of generality, the supports of λbig and λsmall are disjoint. This is because, if a coefficient is nonzero in both, then one could remove it from λsmall and keep it in λbig while maintaining the same approximate degree. Remark 2.6. We note that deg(Roundτ (λ)) ≤ ∥λ∥B1 /τ .
2.2
Fourier analysis of quantum channels
We begin by recalling standard facts about the Fourier analysis of quantum channels. For a more thorough introduction to this topic, we refer the reader to [BY23]. Consider a superoperator Φ which acts on n-qubit states. We can expand Φ in terms of Pauli matrices by defining the functions χP1 ,P2 (ρ) = P1 ρP2 , (2.10) where P1 , P2 ∈ {I, X, Y, Z}⊗n . Such χP1 ,P2 form a basis for the set of superoperators [BY23, Proposition 6], so Φ can be written as X b 1 , P2 ) · χP ,P , Φ= Φ(P (2.11) 1 2 P1 ,P2
b 1 , P2 )’s are referred to as the Fourier coefficients of the superoperator Φ. The next lemma where the Φ(P shows that these functions are actually orthonormal to each other. Lemma 2.7 (The Fourier basis is orthonormal). Let P1 , P2 , Q1 , Q2 ∈ Pn . Then ( h i 1 if P1 = Q1 and P2 = Q2 E tr χP1 ,P2 (R)† · χQ1 ,Q2 (R) = , R∼Pn 0 otherwise. 11
(2.12)
Proof. If P1 = Q1 and P2 = Q2 , h i E tr χP1 ,P2 (R)† · χQ1 ,Q2 (R) = E tr(P2 RP1 · Q1 RQ2 ) = E[tr(RR)] = 1. R
R
R
(2.13)
Otherwise, let us assume without loss of generality that P1 ̸= Q1 . Then, P1 Q1 is some nonidentity Pauli matrix, and a random Pauli R will commute with P1 Q1 with probability 1/2 and anticommute with probability 1/2. Thus, half of the time, we have tr χP1 ,P2 (R)† · χQ1 ,Q2 (R) = tr(P2 RP1 Q1 RQ2 ) = tr(P2 P1 Q1 RRQ2 ) = tr(P2 P1 Q1 Q2 ), (2.14) and the other half of the time, we have tr χP1 ,P2 (R)† · χQ1 ,Q2 (R) = tr(P2 RP1 Q1 RQ2 ) = − tr(P2 P1 Q1 RRQ2 ) = − tr(P2 P1 Q1 Q2 ).
(2.15)
These two average out to zero, which completes the proof. 2.2.1
Local Fourier coefficients
By linearity, Lemma 2.7 gives us the following Fourier inversion formula for the Fourier coefficients of Φ: h i b 1 , P2 ) = E tr χP ,P (R)† · Φ(R) . Φ(P (2.16) 1 2 R
In principle, this suggests a natural experiment that we could carry out to learn the Fourier coefficient b 1 , P2 ) of Φ in the case that Φ is a quantum channel. However, for our application we will be interested in Φ(P efficiently estimating many expectations of this form, and in this case it will only be possible to do so if we restrict the random Paulis R in these expectations to have local support. Motivated by this, let us first show an analogue of Lemma 2.7 for Paulis with local support. Lemma 2.8 (Local orthonormality relations). Let P1 , P2 , Q1 , Q2 ∈ Pn . Let S ⊆ [n] contain SP . Then ( i h 1 if QS1 = P1 , QS2 = P2 , and QS1 = QS2 , † E tr χP1 ,P2 (R) · χQ1 ,Q2 (R) = (2.17) R∼PS 0 otherwise. Proof. By definition, h i E tr χP1 ,P2 (R)† · χQ1 ,Q2 (R) =
tr(P2 RP1 Q1 RQ2 ) = E tr(P2S RP1S QS1 RQS2 ) · tr(QS1 QS2 ) R∼PS h i = E tr χP1S ,P2S (R)† · χQS1 ,QS2 (R) · tr(QS1 QS2 ).
R∼PS
E
R∼PS
R∼PS
(2.18) (2.19) (2.20)
By Lemma 2.7, the expectation is 1 if QS1 = P1S and QS2 = P2S and 0 otherwise. In addition, the second term is 1 if QS1 = QS2 and 0 otherwise. This completes the proof. It will turn out that we only ever care about the case in which S = SP . In this case, Lemma 2.8 says that the functions χP1 ,P2 are no longer orthonormal over Paulis R restricted to SP ; instead, χP1 ,P2 can be “confused” for certain other functions χQ1 ,Q2 which agree with it on SP . We will use Q ⪰ P to mean that Q can be confused with P , i.e. Q⪰P
⇐⇒
S
S
S
S
Q1 P = P1 , Q2 P = P2 , and Q1 P = Q2 P
(2.21)
We define the “local Fourier coefficients” of Φ via b loc (P1 , P2 ) ≜ Φ
E
R∼PS
h i tr χP1 ,P2 (R)† · Φ(R) .
P
12
(2.22)
Note that, unlike the actual Fourier coefficients of Φ, these do not have an interpretation as the coefficients of Φ when expanded in a particular basis of superoperators. Instead, we only have that X b loc (P1 , P2 ) = b 1 , Q2 ). Φ Φ(Q (2.23) Q⪰P
However, it turns out that the actual Fourier coefficients of Φ can still be recovered from these “local Fourier coefficients” by a simple linear transformation. We discuss this in Section 3 below. 2.2.2
Estimating the local Fourier coefficients
We conclude by designing an algorithm to estimate the local Fourier coefficients via random Pauli measurements. First, let us establish some notation which will help to specify the algorithm. Given a single qubit Pauli P ∈ {X, Y, Z} and a bit b ∈ {±1}, we write |P, b⟩ for the eigenstate of P with eigenvalue b. We also set |I, +1⟩ ≜ |0⟩ and |I, −1⟩ ≜ |1⟩, and note that both of these are +1 eigenstates of I. For an n-qubit Pauli Q and a vector v ∈ {±1}n , we define the vector |Q, v⟩ by extending the single-qubit definition via the tensor product. We note that |Q, v⟩ is an eigenvector of Q with eigenvalue χSQ (v), where Y χSQ (v) = vi (2.24) i∈SQ
is the standard Boolean Fourier character. (The fact that we take the product of the vi ’s only within the support of Q accounts for the fact that Q may have coordinates which are equal to I.) As a result, we have the eigendecomposition X Q= χSQ (v) · |Q, v⟩⟨Q, v|. (2.25) v∈{±1}n
We are now ready to state our algorithm. Proposition 2.9 (Single sample unbiased estimators). Let P1 , P2 ∈ Pn with sP ≤ k, and let 1 ≤ i ≤ M . (i) b loc (P1 , P2 ). Then the quantity EP1 ,P2 from Algorithm 1 is an unbiased estimator for Φ Proof. We drop the (i) superscript from our Pauli matrices for notational convenience. We also write S for SP and s for sP . (i) (i) First, we consider the expectation of EP1 ,P2 conditioned on a fixed A, B. If Z ̸= P(P2 RP1 ), then EP1 ,P2 is set to 0, so the expectation is 0. Otherwise, we receive the measurement outcome |B, w⟩ with probability tr |B, w⟩⟨B, w| · Φ(|A, v⟩⟨A, v|) . (2.26) Thus, we can write the expectation as X (i) E[EP1 ,P2 | A, B] = E n 4s · c(P2 RP1 ) · χSR (v) · χSZ (w) · tr |B, w⟩⟨B, w| · Φ(|A, v⟩⟨A, v|) v∈{±1}
w∈{±1}n
= 2−n 4s · c(P2 RP1 ) · tr
X
χSZ (w) · |B, w⟩⟨B, w| · Φ(
w
X
(2.27) χSR (v) · |A, v⟩⟨A, v| (2.28)
v
= 2−n 4s · c(P2 RP1 ) · tr(Z · Φ(R))
(2.29)
s
= 4 tr(P2 RP1 · Φ(R)),
(2.30)
where we used Equations (2.2) and (2.25) in the last two steps. Now, note that as A varies uniformly over Pn , R varies uniformly over PS . Furthermore, conditioned on the value of A, we have Z = P(P2 RP1 ) with probability exactly 1/4s . Hence, (i) tr(χP1 ,P2 (R)† · Φ(R)) . (2.31) E[EP1 ,P2 ] = E R∼PS
P
This concludes the proof. 13
Algorithm 1: Local Fourier coefficient estimation Input: Black-box ability to apply an unknown n-qubit quantum channel Φ; locality parameter k; error parameters ε, δ > 0. bP ,P for all P1 , P2 ∈ Pn with s ≤ k such that Output: Estimates E 1 2 P bP ,P − Φ b loc (P1 , P2 )| ≤ ε. |E 1 2 Set M = Θ(C k log(n/δ)/ε2 ) for some absolute constant C > 0. 2 for i = 1, . . . , M do 3 Initialize a uniformly random Pauli eigenstate |A(i) , v (i) ⟩. 4 Apply Φ to this state. 5 Sample a uniformly random B (i) ∈ Pn and measure in this basis. 6 Record the outcome w(i) ∈ {±1}n . 1
7 8
for P1 , P2 with sP ≤ k do For each i = 1, . . . , M , define the Pauli matrices R(i) = (A(i) )SP ⊗ I SP ,
and
Z (i) = (B (i) )SP ⊗ I SP .
If Z (i) = P(P2 R(i) P1 ), set (i)
EP1 ,P2 = 4sP · c(P2 R(i) P1 ) · χSR(i) (v (i) ) · χSZ (i) (w(i) ). 9
Otherwise, set it to 0. Compute the estimate M X (i) bP ,P ≜ 1 · E . E 1 2 M i=1 P1 ,P2
14
Lemma 2.10. Let k > 0 be a locality parameter. Let M = Θ(C k log(n/δ)/ε2 ), where C > 0 is an absolute constant. Then, the outputs of Algorithm 1 satisfy bP ,P − Φ b loc (P1 , P2 )| ≤ ε |E 1 2
(2.32)
for all P1 , P2 ∈ Pn with sP ≤ k with probability at least 1 − δ, using M queries to the channel Φ. Moreover, Algorithm 1 runs in time O(M nk ). (i) b loc (P1 , P2 ) which is bounded in Proof. Consider a fixed P1 and P2 . Each EP1 ,P2 is an unbiased estimator for Φ k magnitude by 4 . As a result, Hoeffding’s inequality, applied to the real and imaginary parts of the estimator, implies that bP ,P − Φ b loc (P1 , P2 )| ≥ ε] ≤ 4e−4M ε2 /16k . Pr[|E (2.33) 1 2 Now, the number of P1 , P2 ∈ Pn with sP ≤ k is at most nk · 16k ≤ (16n)k . Hence, by the union bound, the probability that there exists an EP1 ,P2 with error more than ε is at most 2
k
(16n)k · 2e−2M ε /16 ≤ δ,
(2.34)
by our choice of M . This completes the proof.
3
Local Fourier coefficients of the time evolution operator
An important ingredient of this work is the local Fourier coefficients of the time evolution operator corresponding to our Lindbladian eLx t (·). We begin by introducing some notation we will use to represent these local Fourier coefficients. Definition 3.1 (Vector of expectation values). For a (parameterized) Lindbladian Lx and a time t ∈ R, we define the vector E : Cm → Cm as follows: (E(x))P1 ,P2 =
[tr(χP1 ,P2 (R)† · eLx t (R))].
E
R∼PS
(3.1)
P
By a Taylor series expansion, we can write eLx t (R) =
∞ ℓ X t
Lℓx (R),
(3.2)
E [tr(χP1 ,P2 (R)† · Lℓx (R))] ℓ! R∼PSP
(3.3)
t E [tr(χP1 ,P2 (R)† · Lℓx (R))]. ℓ! R∼PSP
(3.4)
ℓ=0
ℓ!
where R ∈ Pn . Hence, we have (E(x))P1 ,P2 = =
∞ ℓ X t ℓ=0 ∞ ℓ X ℓ=1
In the second step, we used the fact that the ℓ = 0 expectation is E
R∼PS
[tr(χP1 ,P2 (R)† · L0x (R))] =
P
E
R∼PS
[tr(χP1 ,P2 (R)† · R)] = 0,
(3.5)
P
due to Lemma 2.8, and the fact that at least one of P1 , P2 is non-identity. Precisely understanding the infinite sum in Equation (3.4) is challenging; however, we show in Section 4 that it is well-approximated by its linear term (the ℓ = 1 term). Motivated by this, we dedicate this section to understanding this linear term. The ℓ = 1 term in Equation (3.4) is given by t·
E
R∼PS
[tr(χP1 ,P2 (R)† · Lx (R))]
P
15
(3.6)
Recalling the definition of a (parameterized) Lindbladian, we have X 1 xQ1 ,Q2 Q1 RQ2 − {Q2 Q1 , R} Lx (R) = 2 Q1 ̸=I,Q2 X 1 1 xQ1 ,Q2 χQ1 ,Q2 (R) − c(Q2 Q1 ) · χP(Q2 Q1 ),I (R) − c(Q2 Q1 ) · χI,P(Q2 Q1 ) (R) . = 2 2
(3.7) (3.8)
Q1 ̸=I,Q2
Now, we can use Lemma 2.8 to calculate the expectation. The simplest case is when P1 , P2 ̸= I. In this case, X E [tr(χP1 ,P2 (R)† · Lx (R))] = xQ . (3.9) R∼PS
P
Q⪰P
On the other hand, let P = (P, I), with P ̸= I (we do not use the P = (I, P ) case). Of the three terms in Equation (3.8), the first can be “confused” with (P, I) exactly when Q ⪰ (P, I); the second can be “confused” with (P, I) when P(Q2 Q1 ) = P ; and the third can never be “confused” with (P, I). In the second of these cases, note that Q2 = P(P Q1 ) and c(Q2 Q1 ) = c(P Q1 ); this is because ⇒
P Q1 = c(P Q1 )P(P Q1 ) = c(P Q1 )Q2
Q2 Q1 = c(P Q1 )P = c(Q2 Q1 )P.
(3.10)
1 2
(3.11)
As a result, R∼PS
X
[tr(χP,I (R)† · Lx (R))] =
E
P
xQ −
Q⪰(P,I)
X
=
xQ −
X
c(Q2 Q1 ) · xQ1 ,Q2
Q1 ̸=I, Q2 , P(Q2 Q1 )=P
1X c(P Q) · xQ,P(P Q) . 2
(3.12)
Q̸=I
Q⪰(P,I)
To help us analyze these expressions, we introduce the following notation. Notation 3.2. For 0 ≤ k ≤ n, we write Ak for the square matrix whose rows and columns are indexed by pairs P1 , P2 with P1 ̸= I which acts as follows: (Ak x)P ≜
[tr(χP (R)† · Lx (R))].
E
R∼PS
(3.13)
P
From Equations (3.9) and (3.11), we have ( P x PQ⪰P Q P (Ak x)P = 1 Q̸=I c(P Q) · xQ,P(P Q) Q⪰(P,I) xQ − 2
if P1 , P2 ̸= I, if P = (P, I).
(3.14)
Applying this notation to Equation (3.4), we have that EP (x) = t · (Ak x)P +
∞ ℓ X t ℓ=2
E [tr(χP1 ,P2 (R)† · Lℓx (R))]. ℓ! R∼PSP
(3.15)
The linear ℓ = 1 term of our expansion has a nice interpretation as the local Fourier coefficients of the Lindbladian Lx . From Equation (3.14), we see that these expressions are a linear combination of the true Lindbladian parameters x. We are interested in the inverse of A, i.e., how to recover the true Lindbladian parameters if we know either the local Fourier coefficients or approximations of them. To understand this, we first introduce the following matrix. Notation 3.3. For 0 ≤ k ≤ n, we write Vk for the m × m matrix which acts as follows. For P1 , P2 ̸= I, X (Vk y)P ≜ (−1)sQ −sP yQ . (3.16) Q⪰P
Otherwise, if P = (P, I), (Vk y)P,I ≜ 2yP,I +
X Q,SQ ⊆SP Q̸=I,P
c(P Q) · yQ,P(P Q) +
X
(−1)sQ −sP yQ −
Q⪰(P,I) Q̸=(P,I)
16
X
(−1)sQ −sP yQ .
Q⪰(I,P ) Q̸=(I,P )
(3.17)
Next, we show that the Ak matrix is invertible and that its inverse is equal to Vk . This implies that it is possible to recover the true Lindbladian parameters if we know the local Fourier coefficients exactly. To begin, we need the following helper lemma. Lemma 3.4 (Helper lemma). Suppose that R ⪰ P . Then X 1 sQ −sP (−1) = 0
if R = P , otherwise.
(3.18)
Q:Q⪰P R⪰Q
Proof. The pairs Q which satisfy R ⪰ Q and Q ⪰ P are exactly those which (i) agree with P on SP , (ii) agree with R on some subset T of SR \ SP , and (iii) are identity on the remaining qubits in SR . As R is non-identity on every qubit in T , we have sQ − sP = |T |. Thus, X
sQ −sP
(−1)
=
X
|T |
(−1)
=
T ⊆SR \SP
Q:Q⪰P R⪰Q
1 if R = P , 0 otherwise.
(3.19)
This completes the proof. To prove that Vk is the inverse of Ak , we first show that it successfully recovers any Lindbladian parameter xP1 ,P2 with P1 , P2 ̸= I. Lemma 3.5. For any P1 , P2 ̸= I, we have (Vk Ak x)P = xP . Proof. To see this, (Vk Ak x)P =
X
(−1)sQ −sP (Ak x)Q
(3.20)
Q⪰P
=
X
(−1)sQ −sP
Q⪰P
=
X R⪰P
X
xR
(3.21)
(−1)sQ −sP = xP ,
(3.22)
R⪰Q
xR
X Q:Q⪰P R⪰Q
where the last step used Lemma 3.4. This completes the proof. Next, we show that Vk also recovers any Lindbladian parameter xP with P = (P, I). Lemma 3.6. For any P ̸= I, we have (Vk Ak x)P,I = xP,I . Proof. To see this, note that (Vk Ak x)P,I is equal to X X 2(Ak x)P,I + c(P Q) · (Ak x)Q,P(P Q) + (−1)sQ −sP (Ak x)Q − Q,SQ ⊆SP Q̸=I,P
Q⪰(P,I) Q̸=(P,I)
X
(−1)sQ −sP (Ak x)Q . (3.23)
Q⪰(I,P ) Q̸=(I,P )
Note that only the first term involves indexing (Ak x) by a pair of Paulis, one of which is I; the other terms always index by two non-identity Paulis. Hence, the first two terms are equal to X X X X 2 xQ − c(P Q) · xQ,P(P Q) + c(P Q) · xR . (3.24) Q⪰(P,I)
Q,SQ ⊆SP Q̸=I,P
Q̸=I
R⪰(Q,P(P Q))
Note that in the second summation, the pair (Q, P(P Q)) has support equal to SP and is identity outside of it. This means that R is equal to (Q, P(P Q)) within SP , and R1 and R2 agree outside of SP . This means that
17
(i) R = (R1 , P(P R1 )), (ii) c(P R1 ) = c(P Q), and (iii) R ranges over all possible pairs of this form, subject to R1 |SP not being I or P . Hence, X X X xQ − c(P Q) · xQ,P(P Q) + c(P R) · xR,P(P R) (3.25) (3.24) = 2 X
=2
X
X
c(P Q) · xQ,P(P Q)
(3.26)
xQ,P Q
(3.27)
Q̸=I, Q|SP =I or P
X
xQ −
Q̸=I, Q|SP =I or P
Q⪰(P,I)
=2
X
xQ −
Q⪰(P,I)
=2
R:R|SP ̸=I,P
Q̸=I
Q⪰(P,I)
X
xQ − xP,I −
xQ,P Q ,
(3.28)
Q̸=I,P Q|SP =I or P
Q⪰(P,I)
where in the second-to-last line we use the fact that if QSP = I or P , then P Q is a Pauli, and so P(P Q) = P Q and c(P Q) = 1. Now, the third term in Equation (3.23) is equal to X X X xR (3.29) (−1)sQ −sP (Ak x)Q = (−1)sQ −sP Q⪰(P,I) Q̸=(P,I)
Q⪰(P,I) Q̸=(P,I)
X
=
R⪰Q
X
(−1)sQ −sP
Q⪰(P,I)
xR −
R⪰Q
X
xR .
(3.30)
R⪰(P,I)
The first of these terms is equal to X
X
(−1)sQ −sP = xP,I ,
(3.31)
by Lemma 3.4. Hence, the third term in Equation (3.23) is equal to X X xP,I − xR = − xR,P R .
(3.32)
R⪰(P,I)
xR
Q:Q⪰(P,I) R⪰Q
R̸=I,P, R|SP =P
R⪰(P,I)
Similarly, the fourth term in Equation (3.23) is equal to X xR,P R .
(3.33)
R̸=I,P, R|SP =I
Plugging everything back into Equation (3.23), we get that X X (Vk Ak x)P,I = 2 xQ − xP,I − xQ,P Q − Q̸=I,P, Q|SP =I or P
Q⪰(P,I)
=2
X
xQ − xP,I − 2
Q⪰(P,I)
X
xQ,P Q
X R̸=I,P, R|SP =P
xR,P R +
X
xR,P R
(3.34)
R̸=I,P, R|SP =I
(3.35)
Q̸=I,P, Q|SP =P
= 2xP,I − xP,I
(3.36)
= xP,I ,
(3.37)
where in the third step we used that if Q = (Q, P Q) and Q|SP = P , then Q ⪰ (P, I). This completes the proof. 18
Combining these two lemmas, we have the following corollary. Corollary 3.7 (Inverse of Ak ). Ak is invertible, and its inverse is A−1 k = Vk . Not only do we want A to be invertible, we also want both it and its inverse to be well-behaved. We show that they are well-behaved in a precise technical sense in the following lemma. Lemma 3.8 (Ak is well-behaved). Let A = Ak for 0 ≤ k ≤ n. Then ∥A∥B1 →∞ ≤ ∥A∥B1 →B1 ≤ 4k , −1
∥A
−1
∥B1 →∞ ≤ ∥A
(3.38) k
∥B1 →B1 ≤ 4 .
(3.39)
As the rows of A and A−1 in which P1 , P2 = ̸ I behave much differently than the rows in which P2 = I, we will handle these two cases separately. First, we consider the case of P1 , P2 ̸= I. Lemma 3.9. Let A = Ak . Let Π be the projector onto the rows (P1 , P2 ) with P1 , P2 ̸= I. Then ∥ΠA∥B1 →∞ ≤ ∥ΠA∥B1 →B1 ≤ 2k , −1
∥ΠA
(3.40)
−1
k
∥B1 →B1 ≤ 2 .
(3.41)
X
|xQ1 ,Q2 |.
(3.42)
∥B1 →∞ ≤ ∥ΠA
Proof. For P1 , P2 ̸= I, we consider bounding the entry |(Ax)P1 ,P2 | ≤
(Q1 ,Q2 )⪰(P1 ,P2 )
Note that |(A−1 x)P1 ,P2 | is also bounded by the same quantity, so the entirety of the following argument will work for it as well. We will now show an equivalent way to write this expression which allows us to derive our desired bound. For each set S ⊆ [k], we will construct a vector xS as follows: 1. Initialize xSP1 ,P2 = 0 for all pairs of Paulis with sP ≤ k. 2. For each Q = (Q1 , Q2 ) with sQ ≤ k, consider the ℓ ≤ k qubits in which the Paulis are identical, and number them from 1 to ℓ. 3. For each i ∈ S, set the i-th identical pair in Q1 and Q2 to the identity I. Call the resulting Paulis R1 and R2 . If there is some i ∈ S for which there is no corresponding identical pair (meaning that i > ℓ), do not update xS . 4. Otherwise, update xSR1 ,R2 ← xSR1 ,R2 + |xQ1 ,Q2 |. Note that for each (Q1 , Q2 ) ⪰ (P1 , P2 ), there is exactly one choice of S so that (R1 , R2 ) = (P1 , P2 ). Thus, X X |(Ax)P1 ,P2 | ≤ |xQ1 ,Q2 | = (xS )P1 ,P2 . (3.43) (Q1 ,Q2 )⪰(P1 ,P2 )
S⊆[k]
Moreover, xS has the same locality properties as x. In particular, ∥xS ∥B1 ≤ ∥x∥B1 . This gives the desired bounds. Lemma 3.10. Let A = Ak . Let Π be the projector onto the rows (P1 , P2 ) with P1 , P2 ̸= I. Then ∥ΠA∥B1 →∞ ≤ ∥ΠA∥B1 →B1 ≤ 2, ∥ΠA
−1
∥B1 →∞ ≤ ∥ΠA
−1
(3.44)
∥B1 →B1 ≤ 2.
(3.45)
Proof. Let x be a vector such that ∥x∥B1 ≤ 1. Let i ∈ [n]. Then, we can bound the B1 norm of (ΠA−1 )x associated with site i ∈ [n] using Equation (3.17) as X X X X X |(A−1 x)P,I | ≤ 2xP,I + c(P Q)·xQ,P(P Q) + (−1)sQ −sP xQ − (−1)sQ −sP xQ P :SP ∋i
P :SP ∋i
Q,SQ ⊆SP Q̸=I,P
Q⪰(P,I) Q̸=(P,I)
Q⪰(I,P ) Q̸=(I,P )
(3.46) 19
To understand this expression, note that the second term ranges over all (Q, P Q) with SQ ⊆ SP and Q = ̸ I, P ; the third term ranges over all (Q, P Q) where Q agrees with P on SP (and arbitrary outside) and Q ̸= P ; and the fourth term ranges over all (Q, P Q) where Q agrees with I on SP (and arbitrary outside) and Q = ̸ I. Together, every term of the form (Q, P Q) appears at most once, and Q = I and P both appear zero times. Hence, X X X |(A−1 x)P,I | ≤ 2|xP,I | + |xQ,P(P Q) | (3.47) P :SP ∋i
P :SP ∋i
X
≤
Q̸=I,P
2|xP,I | +
P :SP ∋i
X
|xP |
(3.48)
P :i∈SP , P1 ,P2 ̸=I
≤ 2∥x∥B1
(3.49)
≤ 2.
(3.50)
As for (ΠA)x, the B1 norm associated with site i ∈ [n] is X
|(Ax)P,I | ≤
P :SP ∋i
≤
X X
|xQ | +
P :SP ∋i
Q⪰(P,I)
X
|xP1 ,P2 | +
P :i∈SP
1X |xQ,P(P Q) | 2
(3.51)
Q̸=I
1 2
X
|xP |
(3.52)
P :i∈SP , P1 ,P2 ̸=I
≤ 2∥x∥B1
(3.53)
≤ 2.
(3.54)
Both of these held for all i ∈ [n], so this completes the proof. Combining the previous two lemmas with the triangle inequality and the fact that 2k + 2 ≤ 4k for k ≥ 1 yields Lemma 3.8 as a consequence.
4
Series expansions
We will continue our study of the local Fourier coefficients of the time evolution operator, defined in Section 3 as (E(x))P1 ,P2 = E [tr(χP1 ,P2 (R)† · eLx t (R))]. R∼PS
P
We showed in Equation (3.15) that these coefficients, when Taylor expanded as a function of t, can be expressed as ∞ ℓ X t EP (x) = t · (Ak x)P + E [tr(χP1 ,P2 (R)† · Lℓx (R))]. (4.1) ℓ! R∼PSP ℓ=2
In this section, we will show that this Taylor series converges and concentrates around its first-order term t · (Ak x)P for sufficiently small t. Specifically, we will show that for Lindbladians with ∥L∥B1 ≤ g, t only needs to be smaller than roughly 1/g for this series to converge. Our main goal is to show the following lemma. Lemma 4.1 (Operator norm bound on higher-order terms). Suppose that ∥x∥B1 ≤ g. Suppose t > 0 satisfies t < 1/(4ekg). Let A be the matrix defined in Notation 3.2. Then, the Jacobian J(x) of E(x)/t satisfies ∥J(x) − A∥B1 →B1 ≤ (23k)!gt. Note that E(x) = tAx + B(x), where the Jacobian of B(x)/t is J(x) − A.
20
(4.2)
We would like a bound on the Jacobian, because the analysis of our algorithm will ultimately compare the expectations of an estimate E(x) to the true expectations, E(λ). This Jacobian bound tells us that, if x is close to λ, then the difference in the corresponding expectations can be explained by a linear term, along with a (smaller) higher-order term. To prove this statement, we will bound Equation (4.1) in a term-by-term manner. In particular, let us define the expression (ℓ)
EP (x) ≜
E
R∼PS
[tr(χP1 ,P2 (R)† · Lℓx (R))]
(4.3)
[tr(P2 RP1 · Lℓx (R))]
(4.4)
[tr(L†ℓ x (P2 RP1 ) · R)].
(4.5)
P
=
E
R∼PS
P
=
E
R∼PS
P
We will prove a bound on the Jacobian of this expression, which will then extend to a bound on the Jacobian of E(x) itself. To do so, we will carefully control the complexity of the superchannel L†ℓ x (·) as a function of the growing parameter ℓ. Intuitively, if Lx (and therefore L†x ) is local, then L†ℓ x (·) should remain reasonably local provided that ℓ is reasonably small; this will require us to show an expression for L†ℓ x (·) known as a cluster expansion, which is an expansion of L†ℓ (·) in terms of local components known as clusters. For this, it x will be crucial that we study the adjoint L†x rather than the Lindbladian Lx itself. For brevity, in this section, we often index terms by a instead of P , as described in Definition 2.2.
4.1
Cluster expansions
Definition 4.2 (Multisets and clusters). We refer to unordered multisets with the notation a = {a1 , . . . , aℓ }, where the elements need not Q be distinct. We denote its cardinality as |a| = ℓ, and we denote its support as Sa . We also denote xa ≜ a∈a xa . We call a a cluster if it is connected in the dual interaction graph (there is an edge between a and b if SP a ∩ SP b ̸= ∅). We call a a cluster from S if a ∪ S is a cluster in the modified dual interaction graph where there is an additional term for S. Lemma 4.3 (Cluster expansion of Lindbladians). Let O be an operator whose support is contained in (though not necessarily equal to) a set S ⊆ [n]. Then L†ℓ x (O) is a degree-ℓ matrix-valued polynomial which we can write in the following way: X ℓ L†ℓ xa GS,a (O)Ja ∪ S is a clusterK, x (O) = 2 ℓ! a:|a|=ℓ
P d where GS,a is a superoperator with bounded Fourier weight, Q1 ,Q2 |Gd S,a (Q1 , Q2 )| ≤ 1, and GS,a (P1 , P2 ) is only nonzero provided SP ⊆ Sa . Note that GS,a depends on the subset S but not O. Proof. We prove the lemma by induction on ℓ. For the base case of ℓ = 0, we have L†ℓ (O) = O. Then, consider GS,a (O) = O, i.e., Ga is the identity superoperator, for a = ∅. This satisfies the required conditions on its Fourier coefficients, because the only nonzero Fourier coefficient is Gd S,a (I, I) = 1. In addition, SI,I = ∅ = Sa . Moreover, a ∪ S is vacuously a cluster. For the inductive step, suppose the result holds for ℓ. Then, L†x ℓ+1 (O) = L†x (L†ℓ x (O)) X 1 †ℓ †ℓ = xP1 ,P2 P2 L (O)P1 − {P2 P1 , L (O)} 2
(4.6) (4.7)
P1 ,P2 P1 ̸=I ℓ
= 2 ℓ!
X
X
P1 ,P2 a:|a|=ℓ P1 ̸=I
xP1 ,P2 · x
a
1 P2 GS,a (O)P1 − {P2 P1 , GS,a (O)} Ja ∪ S is a clusterK. 2
In the second line, we use the definition of L† . In the last line, we use the inductive hypothesis. 21
(4.8)
Now, suppose a satisfies that a ∪ S is a cluster, and let us consider the corresponding term 1 2ℓ ℓ! · xP1 ,P2 · xa P2 GS,a (O)P1 − {P2 P1 , GS,a (O)} Ja ∪ S is a clusterK 2 1 = 2ℓ ℓ! · xP1 ,P2 · xa [χP2 ,P1 − (χP2 P1 ,I + χI,P2 P1 )](GS,a (O)) Ja ∪ S is a clusterK. 2
(4.9) (4.10)
This term corresponds to the new cluster b = a ∪ {SP } inside L†x ℓ+1 (O). Indeed, note that xP1 ,P2 · xa = xb . The induction hypothesis tells us that GS,a is a linear combination of terms of the form Q1 OQ2 with Q ⊆ Sa . Hence, Equation (4.10) is a linear combination of terms of the form 1 ℓ a 2 ℓ! · xP1 ,P2 · x [χP2 ,P1 − (χP2 P1 ,I + χI,P2 P1 )](Q1 OQ2 ) Ja ∪ S is a clusterK. (4.11) 2 Note that this is in turn a linear combination of terms of the form Q′1 OQ′2 with SQ′ ⊆ SQ ∪SP ⊆ Sa ∪SP = Sb . Furthermore, note that since the total support of Q1 OQ2 is contained in Sa ∪S, Equation (4.11) is nonzero only if SP overlaps with Sa ∪ S. Since we know that a ∪ S is a cluster, this is equivalent to Sa ∪ {SP } ∪ S = Sb ∪ S being a cluster. By the triangle inequality, the expression in Equation (4.10) has Fourier weight at most 2 · 2ℓ ℓ! = 2ℓ+1 ℓ!. Now, we collect all the terms in the sum associated to the monomial xb . There are at most ℓ + 1 of them, corresponding to the clusters formed by removing one element from a along with the element removed. Hence, their total Fourier weight is at most (ℓ + 1) · 2ℓ+1 ℓ! = 2ℓ+1 (ℓ + 1)!. This gives the desired bound by collecting the corresponding (matrix) coefficient and labeling it Gb (O). This sum is bounded because of the following statement bounding the number of clusters in a boundeddegree (weighted) graph. Lemma 4.4 (Cluster count). Let ℓ ≥ 0. Let x satisfy ∥x∥B1 ≤ g, and let i ∈ [n]. Define X Zi (x) ≜ xa Ja is a cluster from iK.
(4.12)
a:|a|=ℓ
Then |Zi (x)| ≤ (egk)ℓ .
(4.13)
Proof. We first show Equation (4.13). Let r be an integer satisfying 1 ≤ r ≤ k. We begin with the following standard fact (Lemma 4 of [MM24]): Let G = (V, E) be a multihypergraph with maximum degree at most g and rank at most k; then the number of connected subgraphs (sets of edges) of size r containing a vertex v ∈ V in its support is at most (eg(k − 1))r . The analogous statement also holds for weighted graphs: P let we be the nonnegative weight associated to hyperedge e, and let g be now the weighted degree, maxi∈V e∋i |we |. Then X Y JS is a connected subgraph of size r containing vK we ≤ (eg(k − 1))r . (4.14) e∈S
S⊆E
To prove this, let us note that it suffices to prove this when the we ’s are rational, by a continuity argument. But for rational we ’s, we can multiply each we by a scalar such that the weights become integral, and then apply the unweighted statement to the analogous hypergraph. To apply this to our setting, let Gx be the multihypergraph with the vertex set V = [n] and, for each Lindbladian term a ∈ [m], a hyperedge Sa with weight |xa |. Then Equation (4.14) implies that for each vertex i ∈ [n], X Y JS is a connected subgraph in Gx of size r containing iK |xa | ≤ (eg(k − 1))r . (4.15) a∈S
S⊆[m]
ℓ−1 From there, we now consider clusters. For every subgraph S = {a1 , . . . , ar } ⊆ [m] of size r, there are r−1 ways to assign positive integer weights to the r elements of S which sum up to ℓ. Each of these corresponds 22
to a unique cluster a of cardinality |a| = ℓ consisting of r distinct elements; we write a ∼ S for a cluster formed in this manner. Moreover, since ∥x∥B1 ≤ g, the weight xa of the cluster a is at most the weight of the subgraph times g ℓ−r . Therefore, we can bound the number of clusters using the number of subgraphs, giving |Zi (x)| =
ℓ X X
(4.16)
|xa |
(4.17)
r=1 S⊆[m]
a∼S
ℓ X X
XY
(4.18)
r=1 S⊆[m]
≤
≤
=
ℓ X X
JS is a connected subgraph in Gx of size r containing iK · JS is a connected subgraph in Gx of size r containing iK ·
a∼S
X
r=1 S⊆[m]
a∼S a∈S
ℓ X X
Y
JS is a connected subgraph in Gx of size r containing iK ·
r=1 S⊆[m]
=
X
xa
JS is a connected subgraph in Gx of size r containing iK ·
|xa | · g ℓ−r
|xa | · g
a∈S
ℓ−r
ℓ−1 r−1
ℓ X ℓ−1 (eg(k − 1))r · g ℓ−r , r−1 r=1
ℓ X
(4.19)
(4.20)
where we used Equation (4.15) in the last step. Now, using the fact that for a nonnegative number x, we have (4.20) = g ℓ ·
(e(k − 1))r ·
r=1
ℓ−1 r−1
Pℓ
r ℓ−1 r=1 x · r−1
= x(x+1)ℓ−1 ≤ (x+1)ℓ
≤ g ℓ (e(k − 1) + 1)ℓ ≤ g ℓ (ek)ℓ .
(4.21)
This completes the proof. We will also need the following consequence of Lemma 4.4. Lemma 4.5. Let ℓ ≥ 2. Let x satisfy ∥x∥B1 ≤ g, and let v satisfy ∥v∥B1 ≤ 1. Let i ∈ [n]. Then X
va ∂a Zi (x) ≤ e2 kℓ(egk)ℓ−1
(4.22)
X
(4.23)
a
Proof. To begin, let us compute ∂a Zi (x) =
a:a∈a, |a|=ℓ
xa\{a} Ja is a cluster from iK,
where here the notation “a \ {a}” refers to removing a single instance of a from a. Thus, if v is a vector satisfying ∥v∥B1 ≤ 1, we have X X X va ∂a Zi (x) = va · xa\{a} Ja is a cluster from iK. (4.24) a
a
a:a∈a, |a|=ℓ
Note that the absolute value of this quantity is largest when x and v are nonnegative, and so we will henceforth x 1 make this assumption without loss of generality. Now, take u = ℓ−1 ℓ g + ℓ v, and notice that ∥u∥B1 ≤ 1 by construction. Then X (4.25) Zi (u) = ua Ja is a cluster from iK a:|a|=ℓ
=
X ℓ − 1 x 1 a · + · v Ja is a cluster from iK. ℓ g ℓ
a:|a|=ℓ
23
(4.26)
Note that if a = {a1 , . . . , aℓ }, then we can expand ℓ − 1 x 1 a ℓ − 1 ℓ ℓ − 1 ℓ−1 1 X xa + · · + ·v = xa\{a} va + · · · ℓ g ℓ ℓ·g ℓ·g ℓ a∈a ℓ − 1 ℓ−1 1 X · xa\{a} va . ≥ ℓ·g ℓ a∈a
(4.27) (4.28)
1 x a Here, in the first step, we use the binomial formula to expand ( ℓ−1 ℓ · g + ℓ · v) . In the second step, we use the fact that all the terms in this expansion are nonnegative (which follows from the fact that x is nonnegative), to lower bound the expression by only those terms which use a single coordinate of v. Plugging this in to Equation (4.26), we have that X ℓ − 1 ℓ−1 1 X Zi (u) ≥ · xa\{a} va · Ja is a cluster from iK (4.29) ℓ·g ℓ a∈a a:|a|=ℓ
ℓ − 1 ℓ−1 1 X X ≥ · · va xa\{a} Ja is a cluster from iK ℓ·g ℓ a a:a∈a,
(4.30)
|a|=ℓ
ℓ − 1 ℓ−1 1 X · · va ∂a Zi (x). = ℓ·g ℓ a
(4.31)
Rearranging, we have ℓ ℓ−1 ℓ ℓ−1 X va ∂a Zi (x) ≤ ℓ g ℓ−1 · Zi (u) ≤ ℓ g ℓ−1 · (ek)ℓ ≤ ℓ(ek)ℓ eg ℓ−1 . ℓ − 1 ℓ − 1 a
(4.32)
In the second step, we used Lemma 4.4 and the fact that ∥u∥B1 ≤ 1. This concludes the proof. The following corollary follows directly from combining Lemmas 4.3 and 4.4. Corollary 4.6 (Operator norm bound). Let P ∈ Pk , and let ℓ ≥ 2. Suppose Lλ is a k-local superoperator with ∥Lλ ∥B1 ≤ g. Then, ℓ ∥L†ℓ (4.33) λ (P )∥op ≤ ℓ!(2ekg) .
4.2
Bounds on derivatives
We can derive the following as a corollary of Lemma 4.3. Lemma 4.7 (Cluster expansion of Fourier expectations). We can write the function X ℓ E [tr(L†ℓ γa xa Ja is a clusterKJSa ⊇ SP K. x (P2 RP1 )R)] = 2 ℓ! R∼PS
P
(4.34)
a:|a|=ℓ
where γa are some coefficients satisfying |γa | ≤ 1. Proof. We use Lemma 4.3 to expand out L†ℓ x (P2 RP1 ) into a polynomial for every R ∈ PSP : h i X ℓ a E [tr(L†ℓ (P RP )R)] = 2 ℓ! E tr R x G (P RP )Ja ∪ S is a clusterK 2 1 S ,a 2 1 x P P R∼PS
R∼PS
P
= 2ℓ ℓ!
P
X a:|a|=ℓ
(4.35)
a:|a|=ℓ
xa Ja ∪ SP is a clusterK
E
R∼PS
h i tr RGSP ,a (P2 RP1 ) ,
(4.36)
P
|
{z
≜γa
}
where in the first equality, we used the fact that the support of P2 RP1 is contained in SP to apply Lemma 4.3. The coefficients of this expansion can be bounded: h i X X \ \ E G (4.37) |γa | ≤ G (Q , Q ) tr Rχ (P RP ) ≤ 1 2 Q1 ,Q2 2 1 SP ,a (Q1 , Q2 ) ≤ 1, SP ,a Q1 ,Q2
R∼PS
P
Q1 ,Q2
24
where in the final inequality we use the bound on the Fourier weight of GSP ,a . Moreover, because GSP ,a only acts on sites contained in Sa , γa is only nonzero provided that Sa does not merely overlap SP , but contains it. This concludes the proof. Lemma 4.8. Let ℓ ≥ 2. Let J (ℓ) (x) be the Jacobian of the vector-valued function E (ℓ) (x)/t defined in Equation (4.5). Then ∥J (ℓ) (x)∥B1 →B1 ≤ 1t 2ℓ kek+2 ℓk 16k (ℓ + 1)!(egk)ℓ−1 . Proof. Consider some v such that ∥v∥B1 ≤ 1. Further consider some site i ∈ [n]. Then, by Lemma 4.7, X a:Sa ∋i
1 X X vb ∂b Ea(ℓ) (x) t a:Sa ∋i b X X X ℓ 1 = 2 ℓ! vb ∂b γa,a xa Ja is a clusterKJSa ⊇ Sa K t
|(J (ℓ) (x)v)a | =
a:Sa ∋i
= 2ℓ ℓ!
b
1 X ≤ 2ℓ ℓ! t
(4.39)
a:|a|=ℓ
X 1 X X vb ∂b γa,a xa Ja is a cluster from iKJSa ⊇ Sa K t a:Sa ∋i
(4.38)
b
(4.40)
a:|a|=ℓ
X X
a:Sa ∋i a:|a|=ℓ
b
|vb ∂b xa |Ja is a cluster from iKJSa ⊇ Sa K
X 1 X X |vb ∂b xa |Ja is a cluster from iK JSa ⊇ Sa K t a:Sa ∋i a:|a|=ℓ b X X ℓ 1 ℓ(k − 1) 16k ℓ! |vb ∂b xa |Ja is a cluster from iK ≤2 t k−1 a:|a|=ℓ b 1 ℓ(k − 1) ≤ 2ℓ 16k e2 kℓ!ℓ(egk)ℓ−1 t k−1 1 ≤ 2ℓ (16eℓ)k e2 kℓ · ℓ!(egk)ℓ−1 t = 2ℓ ℓ!
(4.41) (4.42) (4.43) (4.44) (4.45)
In the second line, we use Lemma 4.7. In the third line, we use the fact that i ∈ Sa ⊆ Sa , and so a ∪ {{i}} is a cluster. In the fourth line, we use triangle inequality. In the fifth line, we move the order of the sums. The sixth line uses that a cluster with ℓ elements has support size at most ℓ(k − 1) + 1, so the number of possible subsets Sa with k elements (but still containing i) is at most ℓ(k−1) Paulis P k−1 , and the number of possible supported on Sa is at most 16k . The seventh line uses Lemma 4.5. The last line uses that cb ≤ (be/c)c . Since this holds for all i, we have the desired bound. Proof of Lemma 4.1. Let P and Q be k-local. Recall from Equation (4.1) that E(x) can be written as EP (x) = t(Ax)P +
∞ ℓ X t ℓ=2
E [tr(Lℓx (R)P2 RP1 )], ℓ! R∼PSP
(4.46)
where A = Ak is the matrix defined in Notation 3.2. Thus, J(x) − A is given by JP ,Q (x) − AP ,Q =
∞
∞
ℓ=2
ℓ=2
X tℓ (ℓ) 1 X tℓ (x). E [∂xQ tr(Lℓx (R)P2 RP1 )] = J t ℓ! R∼PSP ℓ! P ,Q
(4.47)
Thus, using Lemma 4.8, we have ∥J(x) − A∥B1 →B1 ≤ ≤
∞ ℓ X t
ℓ!
∥J (ℓ) (x)∥B1 →B1
ℓ=2 ∞ ℓ X
1 t
ℓ=2
t ℓ k+2 k k 2 ke ℓ 16 (ℓ + 1)!(egk)ℓ−1 ℓ!
25
(4.48) (4.49)
∞
X 1 = · 4ek+3 16k k 2 gt2 (ℓ + 1)ℓk (2egkt)ℓ−2 t ≤
1 · 4ek+3 16k k 2 gt2 t
ℓ=2 ∞ X
(ℓ + 1)ℓk 2−ℓ+2
ℓ=1 ∞ k+1 X 1 ℓ ≤ · 32ek+3 16k k 2 gt2 t 2ℓ ℓ=1
1 · 32ek+3 16k k 2 gt2 · 2k+1 (k + 1)! t 1 ≤ · 27k+12 gt2 · (k + 3)! t ≤ gt · (23k)!. ≤
(4.50) (4.51) (4.52) (4.53) (4.54) (4.55)
In the second line, we use Lemma 4.8. In the fourth line, we use t < 1/(4egk). In the fifth line, we use that ℓ + 1 ≤ 2ℓ and enlarge the sum to include ℓ = 1.
5
Algorithm
The goal of this section is to prove Theorem 1.1. The detailed version of this theorem is given in Theorem 5.1. Throughout the following section, we will assume that (1) k > 1, (2) g > 0, and (3) ε < g. If (1) fails, then the Lindbladian is easy to learn, since it decomposes into a product of Lindbladians on each qubit, which can be learned separately and in parallel. If either (2) or (3) fails, then outputting the zero Lindbladian suffices.
5.1
Overview of the algorithm
We begin by giving an overview of our algorithm for learning local Lindbladians. Let α, D denote the true parameters that we want to learn. Let L ≜ Lα,D be a Lindbladian with bounded local one-norm ∥L∥B1 ≤ g and approximate degree d ≜ degε/(100·16k ) (L). We assume access to the time evolution operator eLt for a time t satisfying 1 t < tmax ≜ . (5.1) 200 · 4k · (23k)!g Given a vector of coefficients x, we define tr(χP1 ,P2 (R)† eLx t (R)) = EP1 ,P2 (x) ≜ E R∼PS
P
tr(P2 RP1 eLx t (R))
E
R∼PS
(5.2)
P
to be the local Fourier coefficients (as in Sections 2.2.1 and 3) of the corresponding time evolution. Note that EP1 ,P2 (λ) are the local Fourier coefficients corresponding to the true Lindbladian’s time evolution eLt . Estimating the local Fourier coefficients. Our algorithm begins by running Algorithm 1 to produce bP ,P for all P1 , P2 ∈ Pn with s ≤ k such that estimates E 1 2 P bP ,P − EP ,P (λ)| ≤ tη. |E 1 2 1 2
(5.3)
Here, η is an error parameter which we set to η≜
ε . 24000 · 256k d
(5.4)
To accomplish this, we set the “ε” parameter of Algorithm 1 to tη and the “δ” parameter to 0.01, so that Algorithm 1 performs Ck Ck g 2 d2 log(n) Θ log(n) = Θ (5.5) (tη)2 ε2
26
applications of the time evolution eLt , where Ck is a constant that depends only on k. Since each application costs time t, this leads to a total time evolution of Ck gd2 log(n) . (5.6) ttotal = Θ ε2 Estimating the Lindbladian coefficients. The main challenge the algorithm faces is to convert these estimates of the local Fourier coefficients into estimates of the actual Lindbladian parameters. To do so, it maintains a vector x of its estimates for the Lindbladian coefficients and evaluates the quality of its estimates by comparing the local Fourier coefficients of its guessed (parameterized) Lindbladian eLx t with those of the true Lindbladian eLt . Formally, it considers the errors FP1 ,P2 (x) ≜
1 1 EP ,P (x) − EP1 ,P2 (λ). t 1 2 t
(5.7)
We denote the vector of these values as F(x) ≜ (FP1 ,P2 (x))P1 ,P2 .
(5.8)
If all of these errors are small, then x should be close to the true Lindbladian parameters, but if one of these errors is large, then the algorithm updates x in the direction needed to reduce the error. In this way, the algorithm starts with a poor estimate of the true Lindbladian parameters and iteratively improves it until the result is a good estimate. Our algorithm is inspired by the Newton-Raphson root-finding algorithm, as a perfect solution x = λ will cause Equation (5.8) to be equal to 0 and is therefore a root of F(x). Our algorithm can also discover the structure of L. Different Lindbladian terms can interact in ways which are complicated and hard to understand, e.g., the “confusion” Paulis in the sense of Section 2.2.2. Moreover, the presence of large Lindbladian terms can overshadow the contribution of Lindbladian terms which are small but nevertheless still part of the structure. However, we have no trouble extracting information about such large Lindbladian terms, unobscured by the noise of other Lindbladian terms. Inspired by this, our algorithm proceeds in rounds: In the j-th round, the algorithm maintains O(εj )-accurate estimates for every Lindbladian term with magnitude εj = 2−j g or larger. If an estimate is smaller in magnitude than εj /(4d), the algorithm rounds it down to 0. The remaining nonzero coordinates of the current estimate then reflect the structure of L discovered by this iteration. By iteratively decreasing the error threshold, in a given round, we already have good enough estimates of the larger Lindbladian coefficients so that we can effectively filter out their contribution and only detect the smaller terms. One (minor) technical wrinkle is that the algorithm is not able to access the errors in Equation (5.8) exactly. This is for two reasons. First, given access to eLt , we can only approximate the local Fourier coefficients EP1 ,P2 (λ), not compute them exactly. Second, although the algorithm has access to its own estimates x, it still cannot compute EP1 ,P2 (x) exactly, as this involves taking a matrix exponential of the Lindbladian Lx . Instead, the algorithm Taylor expands EP1 ,P2 (x) and truncates at a sufficiently high degree. b As a result, the algorithm works with an approximation F(x) to the error rather than the true error F(x). We describe how the algorithm obtains such an approximation in more detail in Section 5.4. For the purposes b of this section, it suffices to know that we can obtain an approximation F(x) such that the error is bounded as b η(x) ≜ F(x) − F(x), ∥η(x)∥∞ ≤ η always. (5.9)
5.2
The algorithm and guarantee
We now state our algorithm for learning local Lindbladians. The full algorithm is detailed in Algorithm 2. Notably, our algorithm only uses simple experiments of the form: prepare a Pauli eigenstate, apply the unknown evolution eLt , and measure in a Pauli eigenstate. A schematic diagram of these simple circuits is presented in Figure 1. Our algorithm has the following guarantee. We do not attempt to optimize the performance of our algorithm with respect to the locality k.
27
Algorithm 2: Structure learning Lindbladians 1 ε Input: Accuracy ε > 0; time t satisfying t < 200·4k ·(23k)!g ; expectation accuracy η = 24000·256 kd . b = (b b B ≤ ε. b such that ∥λ − λ∥ Output: Estimates λ α, D) 1 (0) 1 Initialize x = 0 ∈ CM and T = ⌈log2 (g/ε)⌉. m l bP ,P (x) defined in Equation (5.81) for degree Γ = log(4/(tη)) − 1 2 Compute the coefficients of E 1
log(1/(2kgt))
2
via [HKT24b]. b b 3 Use Algorithm 1 to obtain estimates E(λ) of E(λ) such that ∥E(λ) − E(λ)∥∞ ≤ ηt/2. 4 for j = 0, . . . , T − 1 do 5 Set εj ≜ 2−j g. 6 Set τj ≜ εj /(800 · 16k d). 7 Update b (j) ) , x(j+1) = Roundε /(4d) x(j) − A−1 Roundτ F(x j
j
(5.10)
where A is defined in Notation 3.2. (T ) b 8 Set λ ≜ x . bP,I and D bP ,P . b P ,P ≜ λ 9 return α bP,I ≜ λ 1
2
1
2
|0⟩ .. .
U
eLt
V
.. .
|0⟩
Figure 1: Quantum experiments in our learning algorithm. All quantum circuits used in our Lindbladian learning algorithm take this form. Here, L is the unknown Lindbladian, and U, V are layers of single-qubit Clifford gates.
28
Theorem 5.1. Let ε, δ > 0. Let α, D be the true parameters, and let L = Lα,D be a k-local Lindbladian with bounded local one-norm ∥L∥B1 ≤ g. Let d = degε/(100·16k ) (L) be the approximate degree of L. Let t, η > 0 be such that 1 ε t< , η< . (5.11) 200 · 4k · (23k)!g 24000 · 256k d b = (b b of the Lindbladian coefficients λ = (α, D) such that Then, Algorithm 2 outputs estimates λ α, D) b − λ∥B ≤ ε with probability at least 1 − δ using a total time evolution of ttotal = O(Ck gd2 log(n/δ)/ε2 ) and ∥λ 1 classical runtime O(nk d log d + (4d)Ck log(dg/ε) + g 2 d2 nk log(n)/ε2 ), where Ck is a constant that depends only b − D∥∞ ≤ ε. on k. We take k = O(1) in the classical runtime. This also implies that ∥b α − α∥∞ ≤ ε and ∥D
5.3
Proof of correctness
b First, to prove the correctness of our algorithm, we assume that we are given estimates F(x) of F(x) such that b η(x) = F(x) − F(x), ∥η∥∞ ≤ η. (5.12) We describe how one may obtain such estimates in Section 5.4. With these estimates, our algorithm simplifies to the form in Algorithm 3. We analyze the algorithm in Theorem 5.2. Algorithm 3: Structure learning algorithm for simplified case 1 b Input: Accuracy ε > 0; time t satisfying t < 200·4k ·(23k)!g ; estimates F(x) satisfying Equation (5.9) ε for given inputs x; η < 24000·256k d . b = (b b B ≤ ε. b such that ∥λ − λ∥ Output: Estimates λ α, D) 1 (0) M 1 Initialize x = 0 ∈ C and T = ⌈log2 (g/ε)⌉. 2 for j = 0, . . . , T − 1 do 3 Set εj ≜ 2−j g. 4 Set τj ≜ εj /(800 · 16k d). 5 Update b (j) ) , (5.13) x(j+1) = Roundε /(4d) x(j) − A−1 Roundτ F(x j
j
where A is defined in Notation 3.2. b ≜ x(T ) . Set λ bP,I and D bP ,P . b P ,P ≜ λ 7 return α bP ≜ λ 1 2 1 2 6
Theorem 5.2. Let ε > 0. Consider Lindbladian parameters α, D, and let L = Lα,D be a k-local Lindbladian with bounded local one-norm ∥L∥B1 ≤ g. Let d = degε/(100·16k ) (L) be the approximate degree of L. Let t, η > 0 be such that 1 ε t< , η< . (5.14) 200 · 4k · (23k)!g 24000 · 256k d b Suppose we can compute estimates F(x) for given inputs x such that b ∥F(x) − F(x)∥ ∞ ≤ η.
(5.15)
b = (b b − λ∥B ≤ ε. This also implies that ∥b b such that ∥λ Then, Algorithm 3 finds estimates λ α, D) α − α∥∞ ≤ ε 1 b and ∥D − D∥∞ ≤ ε. For the sake of analysis, we consider writing the true unknown Lindbladian Lα,D as the parameterized Lindbladian Lλ , where λP,I = αP and λP = DP . Clearly, these representations are equivalent, but Lλ allows us to index into the parameter vector more simply. Before we prove our main theorem, we need to prove some properties of F(x). In particular, one can show that the first order term of A−1 F(x) corresponds precisely to the ideal update we want to perform in the algorithm. To see this, consider expanding each of the entries of F(x) in a Taylor series. By Lemma 4.1, 29
we can write the vector of expectation values as E(x) = tAx + B(x) for some higher order terms B(x). Thus, we have 1 1 1 F(x) = E(x) − E(λ) = A(x − λ) + (B(x) − B(λ)). (5.16) t t t Here, we see that, if we could ignore the higher order terms denoted by B, we would be done. Namely, the update x ← x − A−1 F(x) would directly reveal the unknown parameters λ. Of course, we cannot simply throw away the higher order terms, so one key technical step is to bound the contribution of the higher order terms in the expansion of F(x). We can do so by using the bounds on the higher order terms of the Jacobian J(x) ≜ ∂xQ FP (x) of F(x), which we developed in Section 4. In particular, we can use Lemma 4.1 to bound the higher order terms of F via the Fundamental Theorem of Calculus. Corollary 5.3. Let ∥x∥B1 ≤ g, and let ∆ ≜ x − λ. Let t < 1/(8ekg). Let A = Ak be the matrix defined in Notation 3.2. Then ∥F(x) − (A∆)∥B1 ≤ ct∥∆∥B1 for c ≜ (23k)!2g. (5.17) Proof. Consider fP : [0, 1] → C defined by fP (s) ≜ EP (λ + s∆)/t. Then, by the Fundamental Theorem of Calculus, Z 1
fP (1) − fP (0) = Expanding both sides and using ∂s = FP (x) =
∂s fP (s) ds.
(5.18)
P
Q ∆Q ∂Q , we see that
Z 1X 0
0
Z 1 ∆Q JP ,Q (λ + s∆) ds =
Q
0
(J(λ + s∆)∆)P ds.
(5.19)
Subtracting (A∆) from both sides, we have Z 1 FP (x) − (A∆)P =
0
((J(λ + s∆) − A)∆)P ds.
Bounding the absolute value of this, we have Z 1 ∥F(x) − (A∆)∥B1 ≤ ∥J(λ + s∆) − A∥B1 →B1 ∥∆∥B1 ds ≤ (23k)!2gt∥∆∥B1 ,
(5.20)
(5.21)
0
where in the last inequality, we used Lemma 4.1 applied to λ + s∆, which has ∥λ + s∆∥B1 ≤ 2g. Now, we are ready to prove Theorem 5.2. Proof of Theorem 5.2. Let j ∈ {0, . . . , T − 1}. We prove this via induction on j, where at each iteration, we maintain the invariants ∥x(j) − λ∥B1 ≤ εj , (5.22) deg(x(j) ) ≤ 3d. For the base case of j = 0, recall that x(0) = 0 and ε0 = g. Thus, we have ∥x(0) − λ∥B1 = ∥λ∥B1 ≤ g. Moreover, it is vacuously true that deg(x(0) ) ≤ 3d. For the inductive step, suppose that the inductive hypotheses hold at iteration j. We will prove that they still hold at iteration j + 1. To simplify notation, we drop the iteration index. Let x ≜ x(j) denote the current iterate, x+ ≜ x(j+1) the next iterate, ∆ ≜ x − λ the error vector of the current iterate, ∆+ ≜ x+ − λ the error vector of the next iterate, ε ≜ εj the current error, ε+ ≜ εj+1 = ε/2 the desired error of the next iterate, and τ ≜ τj the current threshold. (Note that setting ε ≜ εj creates a notational conflict with the “ε” used as input to this algorithm. However, in this proof we will only ever use the form ε and never the latter input “ε”.) We use y to denote the next iterate before rounding: b x+ = Roundε/(4d) (y) where y ≜ x − A−1 Roundτ F(x) . (5.23)
30
By the inductive hypothesis, x satisfies Equation (5.22). We will show that x+ satisfies the inductive hypotheses with error parameter ε+ = ε/2. It will suffice to analyze the unrounded vector y and show that ∥y − λ∥B1 ≤
ε . 10
(5.24)
To see why, we will first show that this implies all of the inductive hypotheses for x+ for error parameter ε+ . Recall that the definition of approximate degree splits the superoperator Lλ into two parts Lλ = big small small Lbig , where ∥Lsmall ∥B1 < ε/(100 · 16k ), and minimizes deg(Lbig be λ λ + Lλ λ ). From here on, let Lλ , Lλ the parts attained in this minimization, i.e., Lbig λ =
deg(Lbig ),
argmin Lbig :∥Lλ −Lbig ∥B1 <ε/(100·16k )
Lsmall = Lλ − Lbig λ λ .
(5.25)
As discussed in Section 2.1, without loss of generality, splitting the Lindbladian in this way simply selects small a subset of the coefficients to include in either Lbig . In other words, if λbig are the coefficients of λ or Lλ big small small Lλ and λ are the coefficients of Lλ , then we may assume that the supports of λbig and λsmall are disjoint. Let W big denote the set of pairs of Paulis which are nonzero in λbig , and let W small denote the set of remaining pairs of Paulis. Let Πλbig denote the coordinate projection onto W big , i.e., the indices of small coefficients included in Lbig . Note that Πλsmall = I − Πλbig . λ , and let Πλsmall denote the projection onto W For the first hypothesis in Equation (5.22), we can bound the contributions after projecting onto Πλbig and Πλsmall separately: ∥Πλbig (x+ − λ)∥B1 ≤ ∥Πλbig (x+ − y)∥B1 + ∥Πλbig (y − λ)∥B1 +
≤ ∥Πλbig ∥∞→B1 ∥x − y∥∞ + ∥y − λ∥B1 ε ε ≤d + 4d 10 7ε = . 20
(5.26) (5.27) (5.28) (5.29)
In the second line, we use that, since Πλbig is a coordinate projection, then ∥Πλbig ∥B1 →B1 ≤ 1. In the + third line, we use that d = deg(Lbig λ ) so that ∥Πλbig ∥∞→B1 ≤ d. We also use that x = Roundε/(4d) (y) and Equation (5.24). For the Πλsmall part, let Uy ≜ {P : |yP | ≤ ε/(4d)}, and define ΠUy to be the coordinate projection onto this set. In this way, then x+ = (I − ΠUy )y so that Πλsmall (x+ − λ) = Πλsmall (I − ΠUy )y − Πλsmall λ = Πλsmall (I − ΠUy )(y − λ) − Πλsmall ΠUy λ
(5.30)
= Πλsmall (I − ΠUy )(y − λ) − ΠUy Πλsmall λ,
(5.31)
where in the last step we used the fact that ΠUy and Πλsmall are coordinate projections and hence commute. Bounding the B1 -norm of this, we have ∥Πλsmall (x+ − λ)∥B1 ≤ ∥Πλsmall (I − ΠUy )(y − λ)∥B1 + ∥ΠUy Πλsmall λ∥B1 ≤ ∥Πλsmall ∥B1 →B1 ∥I − ΠUy ∥B1 →B1 ∥y − λ∥B1 + ∥ΠUy ∥B1 →B1 ∥Πλsmall λ∥B1 ε ε + ≤ 10 100 · 16k 41ε ≤ . 400
(5.32) (5.33) (5.34) (5.35)
Here, we use that ∥Πλsmall ∥B1 →B1 , ∥ΠUy ∥B1 →B1 , ∥I − ΠUy ∥B1 →B1 ≤ 1 (since they’re coordinate projections), ∥B1 ≤ ε/(100 · 16k ). Combining these, we have Equation (5.24), and ∥Lsmall λ ∥x+ − λ∥B1 ≤ ∥Πλbig (x+ − λ)∥B1 + ∥Πλsmall (x+ − λ)∥B1 ≤ ε/2.
(5.36)
For the second hypothesis in Equation (5.22), note that deg(x+ ) ≤ deg(Πλbig x+ ) + deg(Πλsmall x+ ) ≤ d + deg(Πλsmall x+ ), 31
(5.37)
where the first inequality uses that I = Πλbig + Πλsmall . The second inequality uses that d = deg(Lbig ). To bound the second term, consider a fixed qubit i ∈ [n]. Notice that ∥Πλsmall x+ ∥B1 ≤ ∥Πλsmall (x+ − λ)∥B1 + ∥Πλsmall λ∥B1 ≤
ε ε 2ε ≤ , + 10 100 · 16k 2
(5.38)
where the last inequality uses Equation (5.34) and ∥Lsmall ∥B1 ≤ ε/(100 · 16k ). Then, because the overall λ + B1 -norm of Πλsmall x is less than ε/2, the number of entries of Πλsmall x+ of magnitude at least ε/(4d) (which, due to the rounding in x+ , are the only nonzero entries in Πλsmall x+ ) whose support contains a fixed qubit i is at most 2d. The degree of Πλsmall x+ is bounded by the number of elements of x+ . Thus, together with Equation (5.37), deg(x+ ) ≤ 3d, so the second hypothesis in Equation (5.22) is satisfied. Now, it remains to prove Equation (5.24). We define two sets of coordinates. Define n o Uτ ≜ (Q1 , Q2 ) : FbQ1 ,Q2 (x) ≥ τ , (5.39) V ≜ (Q1 , Q2 ) : (Ax)Q1 ,Q2 ̸= 0 or (Aλbig )Q1 ,Q2 ̸= 0 , (5.40) b and let ΠUτ and ΠV be the coordinate projections onto Uτ and V , respectively. In particular, ΠUτ (F(x)) = b Roundτ (F(x)). Claim 5.4. Let η > 0 be such that η ≤ τ /(120 · 16k ). Then, ∥ΠUτ ∥∞→B1 ≤
4k+1 ε ε ≤ , τ (30 · 4k η)
∥ΠV ∥∞→B1 ≤ 4k (4d).
(5.41) (5.42)
Proof. Since ∥ΠUτ ∥∞→B1 = maxi∈[n] |{P ∈ Uτ | i ∈ SP }|, we aim to bound the number of elements of Uτ whose supports contain a fixed qubit i ∈ [n]. Call this set Ui,τ , i.e., Ui,τ ≜ {P ∈ Uτ : i ∈ SP }.
(5.43)
b F(x) = F(x) + η(x),
(5.44)
n o Ui,τ ⊆ Q : FQ (x) ≥ τ − η ≥ τ /2 ,
(5.45)
By Equation (5.9), we know that where ∥η(x)∥∞ ≤ η, so
where we use that η ≤ τ /2. By Corollary 5.3, for c = (23k)!2g, ∥F(x)∥B1 ≤ ∥A∆∥B1 + ct∥∆∥B1 ≤ (4k + ct)ε ≤ 2 · 4k ε,
(5.46)
where in the second to last inequality, we use Lemma 3.8 and the inductive hypothesis that ∥∆∥B1 ≤ ε. In the last inequality, we use our choice of t. Thus, because the B1 -norm of F(x) is bounded by 2 · 4k ε, the number of elements of F(x) that are larger than τ /2 in magnitude and which contain qubit i in their support must be at most 4k+1 (ε/τ ). In particular, the size of Ui,τ is bounded by 4k+1 (ε/τ ). This can also be seen via Remark 2.6. Now, we prove the bound on ∥ΠV ∥∞→B1 . By the inductive hypothesis, deg(x) ≤ 3d, and we know deg(λbig ) ≤ d by definition. Then, by Lemma 3.8, deg(Ax) ≤ 4k 3d and deg(Aλbig ) ≤ 4k d. The norm of ΠV is bounded by the sum of these two degree bounds. We can write y − λ in terms of these projectors: b y − λ = x − A−1 ΠUτ F(x) −λ b = ∆ − A−1 ΠU F(x)
(5.47) (5.48)
τ
= ∆ − A−1 ΠUτ (F(x)) −A−1 ΠUτ (η(x)) | {z } ≜err1
32
(5.49)
= ∆ − A−1 ΠUτ (A∆) + A−1 ΠUτ (A∆ − F (x)) + err1 | {z }
(5.50)
≜err2
=∆−A
−1
−1
A∆ + A |
(I − ΠUτ )A∆ + err2 + err1 {z }
(5.51)
≜err3
= err3 + err2 + err1 .
(5.52)
Thus, there are three forms of error we need to bound. The first error err1 comes from only having approximate access to Fb (as in Equation (5.9)). The second error err2 comes from the impact of the higher-order terms of b F(x), which we bounded in Corollary 5.3. The third error err3 comes from the rounding of each F(x). In the following, we bound each error separately. Now, we can bound each of the error terms in Equation (5.52). First, consider err1 , the error from the b approximation to F(x), which we can bound as follows. ∥err1 ∥B1 = A−1 ΠUτ η(x) B1 ≤ ∥A−1 ∥B1 →B1 ∥ΠUτ ∥∞→B1 ∥η(x)∥∞ ≤ 4k
ε ε η= , 30 · 4k η 30
(5.53)
where we use Lemma 3.8, Claim 5.4, and Equation (5.9). In the last inequality, we also use our choice of η. Next, consider err2 , the error from the higher-order terms of F(x). ∥err2 ∥B1 = ∥A−1 ΠUτ (A∆ − F (x))∥B1 −1
≤ ∥A
(5.54)
∥B1 →B1 ∥ΠUτ ∥B1 →B1 ∥(A∆ − F (x))∥B1
(5.55)
k
(5.56)
k
(5.57)
≤ 4 ct∥∆∥B1 ≤ 4 ctε ε . ≤ 30
(5.58)
where c = (23k)!2g. In the third line, we use Lemma 3.8, ∥ΠUτ ∥B1 →B1 ≤ 1 since ΠUτ is a coordinate projection, and Corollary 5.3. In the second to last line, we also use the inductive hypothesis that ∥∆∥B1 ≤ ε. In the last line, we use our choice of t. Finally, consider err3 , the error from rounding. We further break this error up into two parts, corresponding to whether the terms are in V or not. ∥err3 ∥B1 = ∥A−1 (I − ΠUτ )A∆∥B1
(5.59)
−1
≤ ∥A ∥B1 →B1 ∥(I − ΠUτ )A∆∥B1 ≤ 4k ∥(I − ΠUτ )(I − ΠV )A∆∥B1 + ∥(I − ΠUτ )ΠV A∆∥B1 = 4k ∥(I − ΠUτ )(I − ΠV )A∆∥B1 + ∥ΠV (I − ΠUτ )A∆∥B1 , | {z } | {z } ≜err4
(5.60) (5.61) (5.62)
≜err5
where in the third line, we used Lemma 3.8, and in the fourth line, we used the fact that ΠV and (I − ΠUτ ) are coordinate projections and hence commute. To bound err4 , note that (I − ΠV )A∆ = (I − ΠV )A(x − λbig ) − (I − ΠV )Aλsmall = −(I − ΠV )Aλsmall
(5.63)
by definition of V . Then, we can bound ∥err4 ∥B1 ≤ ∥(I − ΠV )Aλsmall ∥B1 ≤ 4k ·
ε ε = , k 100 · 16 100 · 4k
(5.64)
where in the first inequality, we use ∥I −ΠUτ ∥B1 →B1 ≤ 1. In the second inequality, we use ∥I −ΠV ∥B1 →B1 ≤ 1, Lemma 3.8, and ∥Lsmall ∥B1 ≤ ε/(100 · 16k ). To bound err5 , recall that for coordinates not in Uτ , λ |FQ1 ,Q2 (x)| ≤ |FbQ1 ,Q2 (x)| + η ≤ τ + η ≤ 2τ, 33
(5.65)
where we use that η ≤ τ . Consequently, we can conclude that ∥err5 ∥B1 = ∥ΠV (I − ΠUτ )A∆∥B1
(5.66)
≤ ∥ΠV (I − ΠUτ )F(x)∥B1 + ctε
(5.67)
≤ ∥ΠV ∥∞→B1 ∥(I − ΠUτ )F(x)∥∞ + ctε
(5.68)
k
≤ 4 (4d)2τ + ctε ε , ≤ 50 · 4k
(5.69) (5.70)
where c = (23k)!2g. In the second line, we use Corollary 5.3. In the fourth line, we use Claim 5.4 and Equation (5.45). In the last line, we use our choice of t and τ . Overall, plugging the bounds on the three errors back into Equation (5.52), we have ∥y − λ∥B1 ≤ ∥err1 ∥B1 + ∥err2 ∥B1 + ∥err3 ∥B1 ≤ ∥err1 ∥B1 + ∥err2 ∥B1 + 4k (∥err4 ∥B1 + ∥err5 ∥B1 ) ε ε ε ε ≤ + + + 30 30 100 50 ε , ≤ 10
(5.71) (5.72) (5.73) (5.74)
as required. As discussed after Equation (5.24), this completes the proof.
5.4
Sample and time complexity analysis
We analyze the total time evolution, time resolution, and classical runtime of our algorithm from Algorithm 3. b In the previous section, we proved Theorem 5.2, which states that, as long as we can produce estimates F(x) of F(x) to error η in ∞-norm, then we can learn the Lindbladian parameters well. Thus, it suffices to analyze the resources required to obtain such an approximation of F(x). Recall that, for P = (P1 , P2 ), F (x) is defined as FP (x) =
1 1 1 1 E (x) − EP (λ) = E [tr(eLx t (R)P2 RP1 )] − E [tr(eLt (R)P2 RP1 )], t P t t R∼PSP t R∼PSP
(5.75)
where λ are the true parameters of the unknown Lindbladian. We can estimate EP (λ) using our access to eLt for the unknown Lindbladian L, as in Section 2.2.2. Moreover, we can approximate EP (x) by approximating eLx t via a truncated series expansion and computing the terms in this series, similarly to [HKT24b]. The full algorithm is given in Algorithm 2. Using this approach, we have the following guarantee. Theorem 5.5. Let λ = (α, D), and let L = Lλ be a k-local Lindbladian with bounded local one-norm ∥L∥B1 ≤ g. Let d = degε/(100·16k ) (L). Let t, η > 0, where t < 1/(4ekg). Then, for iterates x in Algorithm 3, b there exists an algorithm for computing estimates F(x) of F(x) such that b ∥F(x) − F(x)∥ ∞ ≤η
(5.76)
with probability at least 1 − δ which uses a total evolution time of ttotal = O(C k log(n/δ)/(tη 2 )) for some absolute constant C > 0 andm classical runtime O(knk d log d + dΓ+2 (4Γ + k) poly(Γ) + (Cn)k log(n)/(ηt)2 ), l where Γ =
log(4/(tη)) log(1/(2ekgt)) − 1
.
With our choices of parameters from Theorem 5.2, this gives us Theorem 5.1. Proof of Theorem 5.1. This follows by instantiating Theorem 5.5 with our choice of t, η, Γ. We have that Γ≤
log(4/(tη)) = O(Ck log(dg/ε)), log(1/(2ekgt))
(5.77)
where Ck is some constant that depends only on k. Then, taking k = O(1), the time complexity simplifies to the claimed quantity. 34
Now, we prove Theorem 5.5. Proof of Theorem 5.5. First, we can estimate the expectation values EP (λ) using our query access to eLt . In bP ,P (λ) such that particular, by Lemma 2.10, we can obtain estimates E 1 2 bP ,P (λ) − EP ,P (λ)| ≤ |E 1 2 1 2
tη 2
(5.78)
with probability at least 1 − δ using Θ(C k log(n/δ)/(tη)2 ) queries to the evolution operator eLt , for some absolute constant C > 0. Because each query evolves for time t, this corresponds to the stated total evolution bP ,P (x) such that time. It remains to show that we can obtain estimates E 1 2 bP ,P (x) − EP ,P (x)| ≤ |E 1 2 1 2
tη . 2
To do so, consider Taylor expanding out eLx t in the expression for EP1 ,P2 (x) up to degree log(4/(tη)) Γ≜ −1 . log(1/(2ekgt))
(5.79)
(5.80)
In other words, we approximate bP ,P (x) ≜ E 1 2
Γ X tℓ
E [tr(RL†ℓ x (P2 RP1 ))]. ℓ! R∼PSP
ℓ=1
(5.81)
To show that this indeed satisfies Equation (5.79), by triangle inequality, we have bP ,P (x) − EP ,P (x)| ≤ |E 1 2 1 2
∞ X tℓ E [| tr(RL†ℓ x (P2 RP1 ))|]. ℓ! R∼PSP
(5.82)
ℓ=Γ+1
We can bound each of these terms using Corollary 4.6: | tr(RL†ℓ (P2 RP1 ))| ≤
1 ℓ ∥R∥tr ∥L†ℓ x (P2 RP1 )∥op ≤ ℓ!(2ekg) . N
(5.83)
Plugging back into the above, we get bP ,P (x) − EP ,P (x)| ≤ |E 1 2 1 2
∞ X
(2ekgt)ℓ =
ℓ=Γ+1
tη (2ekgt)Γ+1 ≤ 2(2ekgt)Γ+1 ≤ , 1 − 2ekgt 2
(5.84)
where in the equality, we use the sum of a geometric series when 2ekgt ≤ 1. In the second to last inequality, we use 2ekgt < 1/2. In the last inequality, we use our choice of Γ. We want to bound the time complexity of computing this truncated series expansion. The time complexity of such operations is analyzed in detail in [HKT24b]. The same analysis applies here because the Lindbladian for which we compute the series expansion has degree 3d. This is because, in the proof of Theorem 5.2, we maintain the invariant that x has degree 3d. Thus, the analysis in [HKT24b] for computing the series yields a time complexity of O(kmd log d + ed(1 + e(d − 1))Γ (4Γ + k) poly(Γ)), (5.85) where m is the number of terms in the Lindbladian. In addition, we have the time complexity from Algorithm 1, which by Lemma 2.10, costs O((Cn)k log(n/δ)/(tη)2 ) for an absolute constant C.
5.5
Applications to specific settings
Our main result in Theorem 5.1 is stated in terms of general parameters, namely the approximate degree d and local one-norm bound g. To exemplify the applicability of our result, we instantiate our theorem for various well-studied cases of Lindbladians. In particular, we consider the cases of geometrically local, general 35
k-local, quasi-local, and power-law Lindbladians. We show explicit choices of g, d for each of these settings, resulting in a performance guarantee for our algorithm. Moreover, in Section 5.6, we show that a simplification of our algorithm can be applied to structure learning Hamiltonians with bounded local one-norm. Surprisingly, our result demonstrates that dependence on the effective sparsity parameter of [BLMT24c] is not necessary for structure learning Hamiltonians. The analysis of the algorithm is greatly simplified in this case compared with Lindbladian learning. First, we consider the setting of geometrically local Lindbladians, which is arguably the most well-studied setting in both Hamiltonian and Lindbladian learning. Corollary 5.6 (Learning of strictly local Lindbladians). Let Lλ be a Lindbladian whose coefficients are bounded, |λP | ≤ 1. Further suppose that L is strictly local with respect to some geometry: every term has support size k = O(1), and every qubit interacts with at most d nonzero terms. Then, Algorithm 2 uses a total b − λ∥B < ε with probability at b such that ∥λ evolution time of ttotal = O(d3 log(n)/ε2 ) to learn an estimate λ 1 least 0.99. Moreover, the time resolution is tmin = Θ(1/d). Proof. In this case, ∥λ∥B1 = max
X
i∈[n]
|λP | ≤ max |{P : SP ∋ i}| ≤ d,
(5.86)
i∈[n]
P :SP ∋i
where we use that every qubit interacts with at most d nonzero terms. Moreover, note that degε (λ) ≤ d. We small small can decompose Lλ = Lbig , where Lbig = 0. In this way, ∥Lsmall ∥B1 = 0 < ε, and λ λ + Lλ λ = Lλ and Lλ degε (Lλ ) ≤ deg(Lλ ) = d. Taking g = d in Theorem 1.1 gives the result. For the following two applications, we need the following fact. Fact 5.7 (Volume bounds on lattices). Consider a p-dimensional lattice for p ≥ 2, and let ℓ ≥ 1. The number p of j such that dist(i, j) = ℓ is ≤ 2p p+ℓ−1 ≤ 2 (e(ℓ + 1))p−1 . Consequently, there are at most (2eℓ)pk sets S p−1 of size ≤ k where i ∈ S and dist(i, j) < ℓ for all j ∈ S. Next, we consider learning quasi-local Lindbladians. Such Lindbladians are especially interesting due to the recent surge of interest in quantum Gibbs samplers, which are constructed using quasi-local Lindbladians [CKBG25; CKG23; DLL24; SA26; RSA26; BLMT26; BLMT24a; BCL24; BC26]. Corollary 5.8 (Learning of quasi-local Lindbladians). Let p ≥ 2. Consider a system of n qubits on a p-dimensional lattice, and let Lλ be a k-local Lindbladian. Suppose Lλ is also quasi-local with respect to the lattice, i.e., there is a parameter γ > 0 such that, for all i, j ∈ [n], X |λP | ≤ e− dist(i,j)/γ . (5.87) P {i,j}⊆SP
b such that Then, Algorithm 2 uses a total evolution time of ttotal = O(d2 log(n)/ε2 ) to learn an estimate λ b ∥λ − λ∥B1 < ε with probability at least 0.99, where d = Ck (Cγ(p log(1 + γ) + log(1/ε)))pk ,
(5.88)
where Ck is a constant that depends only on k, and C is an absolute constant. Moreover, the time resolution is tmin = Θ(1). Proof. First, we can easily bound the local one-norm. Let i ∈ [n]. Then, X |λP | ≤ 1,
(5.89)
P :SP ∋i
where in the inequality, we use Equation (5.87) for j = i. This holds for all i ∈ [n], so we see that ∥λ∥B1 ≤ 1. small To bound the approximate degree, consider decomposing Lλ = Lbig , where Lsmall is the part of λ + Lλ Lλ with coefficients indexed by P such that diam(SP ) > r, where r = 4γ(p log(2(1 + 2γ)) + log(100 · 16k /ε)), i.e., λsmall ≜ λP Jdiam(SP ) > rK, λbig ≜ λ − λsmall . (5.90) P 36
Then, if ∥λsmall ∥B1 ≤ ε/(100 · 16k ), we have degε/(100·16k ) (Lλ ) ≤ deg(Lbig ). First, we show that ∥λsmall ∥B1 ≤ ε/(100 · 16k ). Let i ∈ [n]. Then, X X small |λP |= (5.91) |λP | P :SP ∋i
P :SP ∋i diam(SP )>r
X
≤
X
|λP |
(5.92)
e− dist(i,j)/γ
(5.93)
j∈[n] P dist(i,j)≥r/2 {i,j}⊆SP
X
≤
j∈[n] dist(i,j)≥r/2
=
∞ X ℓ=r/2
X
e−ℓ/γ
(5.94)
j∈[n] dist(i,j)=ℓ
∞ X p + ℓ − 1 −ℓ/γ ≤2 e p−1 ℓ=r/2 ∞ X p + ℓ − 1 −ℓ/(2γ) p −r/(4γ) ≤2 e e p−1 p
(5.95)
(5.96)
ℓ=0
= 2p e−r/(4γ) (1 − e−1/(2γ) )−p −r/(4γ)
≤e
(2(1 + 2γ)) ε . ≤ 100 · 16k
p
(5.97) (5.98) (5.99)
In the first line, we use the definition of λsmall . In the second line, we use that if i ∈ SP and diam(SP ) > r, there must exist two points j, a ∈ SP such that dist(j, a) > r, by definition of diameter. Then, by triangle inequality, one of the two points j or a, say j, must satisfy dist(i, j) ≥ r/2. In the third line, we use Equation (5.87). In the fifth line, we use Fact 5.7. In the sixth line, we use that ℓ ≥ r/2 and expand the number of terms in the sum. In the seventh line, we use the identity ∞ X p+ℓ−1 ℓ x = (1 − x)−p (5.100) p−1 ℓ=0
for x = e−1/(2γ) . In the eighth line, we use (1 − e−1/(2γ) )−1 ≤ 1 + 2γ and that p > 0. In the last line, we use our choice of r. Since this holds for any i ∈ [n], this completes the proof that ∥λsmall ∥B1 ≤ ε/(100 · 16k ). Now, we need to bound degε/(100·16k ) (Lλ ) ≤ deg(Lbig ). By definition of λbig , for a given site i ∈ [n], if i ∈ SP and λbig = ̸ 0, then diam(SP ) ≤ r. By Fact 5.7, there are at most (2er)pk possible supports that P satisfy this. This corresponds to at most 16k (2er)pk possible Pauli terms P . Thus, d = degε/(100·16k ) (Lλ ) ≤ 16k (2er)pk ≤ Ck (Cγp log(1 + γ) + k + log(1/ε))pk ,
(5.101)
where C > 0 is an absolute constant, and Ck is a constant that depends only on k. Plugging in these parameters to Theorem 1.1 gives the result. In addition to Lindbladians with exponentially decaying correlations, we can also consider Lindbladians with weaker long-range correlations, where the decay is inverse polynomial rather than inverse exponential. Corollary 5.9 (Learning Lindbladians satisfying power-law decay). Let p ≥ 2. Consider a system of n qubits on a p-dimensional lattice, and let Lλ be a k-local Lindbladian. Suppose Lλ obeys power-law decay with respect to the lattice, i.e., there is a parameter γ > 0 such that, for all i, j ∈ [n], X 1 |λP | ≤ . (5.102) max(1, dist(i, j))γ P {i,j}⊆SP
37
Suppose that p, k = O(1) and γ > p with γ − p = Ω(1). Let κ≜
2pk . γ−p
(5.103)
b such that ∥λ−λ∥ b Then, Algorithm 2 uses a total time evolution of ttotal = O(2γκ log(n)/ε2+κ ) to learn a λ B1 < ε with probability at least 0.99. Moreover, the time resolution is tmin = Θ(1). Proof. As in the quasi-local case, ∥λ∥B1 ≤ 1. To bound the approximate degree, consider decomposing small Lλ = Lbig , where Lsmall is the part of Lλ with coefficients indexed by P such that diam(SP ) > r, λ + Lλ −1/(γ−p) ε(γ−p) where r = 100·16 , i.e., k ·2γ (2e)p λbig ≜ λ − λsmall .
λsmall ≜ λP Jdiam(SP ) > rK, P
(5.104)
We first show that ∥λsmall ∥B1 < ε/(100 · 16k ). The calculation proceeds similarly to the quasi-local case, so we combine a few steps. Let i ∈ [n]. Then, X X X |≤ |λsmall |λP | (5.105) P j∈[n] P dist(i,j)≥r/2 {i,j}⊆SP
P :SP ∋i
≤
X j∈[n] dist(i,j)≥r/2
=
∞ X ℓ=r/2
1 dist(i, j)γ
(5.106)
1 ℓγ
(5.107)
X j∈[n] dist(i,j)=ℓ
≤ 2p ep−1
∞ X (ℓ + 1)p−1
ℓγ
ℓ=r/2
≤ (4e)p ≤ (4e)p
∞ X ℓ=r/2 Z ∞
(5.108)
ℓp−1−γ
(5.109)
xp−1−γ dx
(5.110)
r/2−1 p−γ
r ≤ (4e)p p−γ 2 (γ − p) ε ≤ . 100 · 16k
(5.111) (5.112)
In the second line, we use Equation (5.102). In the fourth line, we use Fact 5.7. In the fifth line, we use that ℓ + 1 ≤ 2ℓ for ℓ ≥ 1. In the last line, we use our choice of r. Now, we need to bound degε/(100·16k ) (Lλ ) ≤ deg(Lbig ). By definition of λbig , for a given site i ∈ [n], if i ∈ SP and λbig = ̸ 0, then diam(SP ) ≤ r. By Fact 5.7, there are at most (2er)pk possible supports that P satisfy this. This corresponds to at most 16k (2er)pk possible Pauli terms P . Thus, for γ − p = Ω(1) and k, p = O(1), we have γ pk/(γ−p) 2 d = degε/(100·16k ) (Lλ ) ≤ 16k (2er)pk = O . (5.113) ε Plugging into Theorem 1.1 gives the result. Finally, we can also instantiate our theorem for general k-local Lindbladians with no geometric constraints. Corollary 5.10 (Learning local Lindbladians). Let Lλ be a k-local Lindbladian on n qubits such that ∥λ∥B1 ≤ g. Then, Algorithm 2 uses a total time evolution of ttotal = O(gn2k−2 log(n)/ε2 ) to learn an estimate b such that ∥λ b − λ∥B < ε with probability at least 0.99. Moreover, the time resolution is tmin = Θ(1/g). λ 1 38
Proof. In this case, we can naively bound d as d ≤ deg(Lλ ) ≤ 16k
n−1 k−1
≤ 16k nk−1 .
(5.114)
Plugging this into Theorem 1.1 gives the result.
5.6
Structure learning Hamiltonians
In this section, we show that our framework can be applied to the problem of structure learning Hamiltonians from both real-time evolution and high-temperature Gibbs states. While we cannot simply instantiate our theorem for new choices of g and d in this case, we find that analyzing our algorithm for Hamiltonians is significantly simpler than that for Lindbladians. The reason for this simplification is the absence of “confusion” (in the sense of Section 2.2.1). Thus, the matrix A defined in Notation 3.2 is the identity matrix, and we can keep track of the P error of our iterates in ∞ → ∞ norm instead. Let H = i P λP P be a k-local Hamiltonian. The factor of i ensures that the coefficients λP are purely imaginary, which is consistent with our definition of the coherent part of a Lindbladian in Section 2.1. We consider two different access models to this Hamiltonian: the ability to evolve under e−iHt for a chosen time t > 0 or access to many copies of high-temperature Gibbs states ρβ = e−βH / tr(e−βH ) for a small inverse temperature β > 0. Similarly to the Lindbladian case, we define a vector of expectation values with entries EP (x) ≜
E
R∼PSP
[tr(e−iH(x)t ReiH(x)t RP )],
or
EP (x) ≜ tr(P ρβ (x)),
(5.115)
for the real-time evolution and Gibbs state cases, respectively7 . Moreover, as before, we define FP (x) ≜
1 1 EP (x) − EP (λ), t t
1 1 FP (x) ≜ − EP (x) + EP (λ), β β
or
(5.116)
respectively8 . Define the Jacobian J(x) of F(x) entrywise as JP,Q (x) ≜ ∂xQ FP (x). As in Section 5, we b consider an approximation F(x) of F(x) such that b η(x) = F(x) − F(x),
∥η(x)∥∞ ≤ η.
(5.117)
One can obtain such an approximation via the same approach as in Section 5.4. Our algorithm is detailed in Algorithm 4. Algorithm 4: Structure learning algorithm for Hamiltonians Input: Accuracy ε > 0; time t satisfying t < tmax or inverse temperature β satisfying β < βmax ; b estimates F(x) satisfying Equation (5.117) for given inputs x. b such that ∥λ − λ∥ b ∞ ≤ ε. Output: Estimates λ (0) M 1 Initialize x = 0 ∈ R and T = ⌈log2 (1/ε)⌉. 2 for j = 0, . . . , T − 1 do 3 Set εj ≜ 2−j . 4 Update b (j) ) . x(j+1) = Roundε /4 x(j) − F(x (5.118) j
5
b ≜ x(T ) . return λ
b First, we prove a general theorem stating that, if we have an approximation for F(x) and a bound on the Jacobian of F(x), then the algorithm in Algorithm 4 obtains good estimates of the Hamiltonian parameters. In the following sections, we prove both of these hypotheses. In particular, we can obtain an approximation b for F(x) similarly to Section 5.4. 7 For the real-time evolution case, we could also use the canonical choice of observables typically used in the Hamiltonian learning literature. However, we use this choice for a closer analogy with the Lindbladian case. 8 Note that the sign is switched for the Gibbs state version to make the first order terms match. Namely, if the signs are not flipped, the first order term of FP for the dynamics version is xP while for the Gibbs state version it is −xP .
39
P Theorem 5.11 (Structure learning of Hamiltonians). Let ε > 0. Let H = i P λP P be a k-local Hamiltonian with ∥λ∥B1 ≤ g and ∥λ∥∞ ≤ 1. Let SH be the support of the Hamiltonian, which is unknown to the algorithm. b Let cg > 0 be a constant depending on g. Let 0 < t, β < 1/(20c2g ). Suppose we can compute estimates F(x) for given inputs x such that ε b ∥F(x) − F(x)∥ . (5.119) ∞ ≤ 20 Also suppose that, for any x such that ∥x∥B1 ≤ g and for any P ∈ Pk , X X |JP,Q (x) − δP,Q | ≤ cg β, (5.120) |JP,Q (x) − δP,Q | ≤ cg t, or Q∈SH
Q∈SH
b such that for the real-time and Gibbs state settings, respectively. Then, Algorithm 4 finds estimates λ b ∥λ − λ∥∞ ≤ ε. First, under the hypothesis of Theorem 5.11, we can also bound the higher order terms of F via the Fundamental Theorem of Calculus. This is the analogue of Corollary 5.3. Lemma 5.12. Let ∥x∥B1 ≤ g, and let ∆ = x − λ. Suppose that ∆P = 0 unless P ∈ SH , and suppose for any P ∈ Pk that X |JP,Q (x) − δP,Q | ≤ cg t. (5.121) Q∈SH
Then, ∥F(x) − ∆∥∞ ≤ c2g t∥∆∥∞ ,
or
∥F(x) − ∆∥∞ ≤ c2g β∥∆∥∞ ,
(5.122)
for the real-time and Gibbs state cases, respectively. Proof. As in the proof of Corollary 5.3, consider fP : [0, 1] → C defined by fP (s) ≜ EP (λ + s∆). Then, by the Fundamental Theorem of Calculus, Z 1 fP (1) − fP (0) = ∂s fP (s) ds. (5.123) 0
Expanding both sides and using ∂s =
P
Q ∆Q ∂Q , we see that
FP (x) =
Z 1X 0
∆Q JP,Q (λ + s∆) ds.
Subtracting ∆ from both sides, we have Z 1X Z 1 X FP (x) − ∆P = ∆Q (JP,Q (λ + s∆) − δP,Q ) ds = ∆Q (JP,Q (λ + s∆) − δP,Q ) ds, 0
(5.124)
Q
(5.125)
0 Q∈S H
Q
where in the last equality, we use that ∆Q = 0 unless Q ∈ SH . Taking the absolute value of both sides, we have Z 1 X |FP (x) − ∆P | ≤ |JP,Q (λ + s∆) − δP,Q ||∆Q | ds (5.126) 0 Q∈S H
≤
Z 1 X
|JP,Q (λ + s∆) − δP,Q | ds · ∥∆∥∞
(5.127)
0 Q∈S H
≤ c2g t∥∆∥∞ .
(5.128)
In the last line, we use the bound on the Jacobian for λ + s∆, where ∥λ + s∆∥B1 ≤ 2g. This holds for all P ∈ Pk , so this gives the claim. The proof is the same for the Gibbs state case. With this, we can prove Theorem 5.11. 40
Proof of Theorem 5.11. In this case, the proof is short. Let j ∈ {0, . . . , T − 1}. We prove this via induction on j, where at each iteration, we maintain the invariants ∥x(j) − λ∥∞ ≤ εj ,
(5.129)
(j)
|xP | ≤ 2|λP |
For the base case of j = 0, recall that x(0) = 0 and ε0 = 1. Thus, we have ∥x(0) − λ∥∞ = ∥λ∥∞ ≤ 1, as (0) required. Moreover, |xP | = 0 ≤ 2|λP | is trivially satisfied. (j) For the inductive step, suppose that ∥x(j) − λ∥∞ ≤ εj and |xP | ≤ 2|λP |. We prove that this holds for iteration j + 1. For brevity, we drop the iteration index. Let x ≜ x(j) denote the current iterate, x+ ≜ x(j+1) the next iterate, ∆ ≜ x − λ the error vector of the current iterate, ∆+ ≜ x+ − λ the error vector of the next iterate, ε ≜ εj the current error, and ε+ ≜ εj+1 = ε/2 the desired error of the next iterate. We use y to denote the next iterate before rounding: x+ ≜ Roundε/4 (y), Note that it suffices to show that ∥y − λ∥∞ ≤
b y ≜ x − F(x).
(5.130)
ε . 10
(5.131)
This implies Equation (5.129) with error parameter ε+ : ∥x+ − λ∥∞ ≤ ∥x+ − y∥∞ + ∥y − λ∥∞ ≤
ε ε ε + ≤ = ε+ . 4 10 2
(5.132)
For the second hypothesis in Equation (5.129), consider a parameter indexed by a Pauli P . Then, either the rounding kicks in so that x+ P = 0, in which case the bound is immediate, or |yP | > ε/4. In the latter case, by the reverse triangle inequality, then ε 3ε |λP | ≥ |yP | − > , (5.133) 10 20 so ε |x+ ≤ 2|λP |, (5.134) P | = |yP | ≤ |λP | + 10 so that Equation (5.129) holds. Thus, it suffices to show that ∥y − λ∥∞ ≤ ε/10. We show this as follows: b y − λ = x − F(x) − λ = ∆ − F(x) −η(x) = ∆ − ∆ + (∆ − F(x)) +err1 = err2 + err1 . | {z } {z } | ≜err1
(5.135)
≜err2
Then, we can bound each of the errors as follows: ∥err1 ∥∞ ≤
ε 20
(5.136)
by the hypothesis of the theorem. Also, ∥err2 ∥∞ = ∥∆ − F (x)∥∞ ≤ c2g t∥∆∥∞ ≤
ε , 20
(5.137)
where we use the inductive hypothesis and Lemma 5.12. Note that the hypothesis of Lemma 5.12, that ∆P = 0 unless P ∈ SH , holds because we maintain the invariant |xP | ≤ 2|λP |. Putting everything together, ∥y − λ∥∞ ≤ ε/10, as required. The proof is the same in the Gibbs state case via Lemma 5.12. 5.6.1
Real-time evolution
We can instantiate Theorem 5.11 when we are given access to real-time evolution under the unknown Hamiltonian H. We do so by proving the hypotheses of Theorem 5.11 hold in this setting. Theorem 5.13 (StructurePlearning of Hamiltonians from real-time evolution). Let ε, δ > 0, and let 0 < t < 1/(162kg). Let H = i P λP P be a k-local Hamiltonian with ∥λ∥B1 ≤ g and ∥λ∥∞ ≤ 1. Then, there b such that ∥λ b − λ∥∞ ≤ ε with probability at least 1 − δ using exists an algorithm that finds estimates λ 2 ttotal = O(g log(n/δ)/ε ). 41
In order to prove this theorem, we require the series expansion properties we proved in Section 4. Define the Jacobian of F as JP,Q (x) ≜ ∂xQ FP (x). We want to prove a bound on the higher order terms of this Jacobian. The reason the bounds in Section 4 do not apply is because we will want to keep track of the error of our estimates in Algorithm 4 in ∞-norm rather than B1 -norm. Thus, the bounds in Section 4 do not apply, as there we only bound the higher order terms of the Jacobian in B1 → B1 norm (Lemma 4.1), not ∞ → ∞ norm. In fact, in the Lindbladian case, the ∞ → ∞ norm of the Jacobian can scale with the system size. Luckily, in the Hamiltonian case, the ∞ → ∞ norm of the Jacobian is bounded. P Lemma 5.14. Let H(λ) = i P λP P be a k-local Hamiltonian with bounded local one-norm ∥λ∥B1 ≤ g. Suppose t > 0 satisfies t < 1/(162kg). Then, for any x such that ∥x∥B1 ≤ g, ∥J(x) − I∥∞→∞ ≤ cg t, where cg = 81kg/20. First, we show that, with this lemma, we can obtain Theorem 5.13. Proof of Theorem 5.13. We proceed as in the proof of Theorem 5.5. First, by Lemma 2.10, we can obtain bP (λ) such that estimates E bP (λ) − EP (λ)| ≤ tε |E (5.138) 40 for all k-local Paulis P with probability at least 1 − δ using Θ(C k log(n/δ)/(tε)2 ) queries to the time evolution operator e−iHt , for some absolute constant C > 0. This corresponds to a total time evolution of Θ(C k log(n/δ)/(tε2 )) = Θ(Ck g log(n/δ)/ε2 ) for some constant Ck > 0 that depends only on the locality k. bP (x) such that As in Theorem 5.5, we can obtain estimates E bP (x) − EP (x)| ≤ |E
tε 40
(5.139)
also by Taylor expanding e−iHt up to degree
log(160/(tε)) Γ≜ −1 , log(1/(2ekgt))
(5.140)
by the same analysis as before. This gives us an approximation of F(x) up to ε/20 error. Moreover, the Jacobian bound needed in Theorem 5.11 is clearly implied by Lemma 5.14: X X |JP,Q (x) − δP,Q | ≤ |JP,Q (x) − δP,Q | = ∥J(x) − I∥∞→∞ ≤ cg t. (5.141) Q∈SH
Q
It remains to prove Lemma 5.14, and we spend the rest of this section proving it. First, we require an inductive bound on the operator norm of a nested commutator. This is an analogue of Corollary 4.6. P Lemma 5.15. Let P Q ∈ Pk . Suppose H = i P λP P is a k-local Hamiltonian with |λP | ≤ 1 and ∥λ∥B1 ≤ g. Then, [H, Q]ℓ = T dT,ℓ T , where supp(T ) ≤ k(ℓ + 1) and X ∥[H, Q]ℓ ∥P,1 ≜ |dT,ℓ | ≤ ℓ!(2kg)ℓ . (5.142) T
Proof.PWe prove this by induction on ℓ. For the base case of ℓ = 0, the claim is clear because [H, Q]0 = Q so that T |dT | = 1 and supp(Q) ≤ k. For the inductive step, suppose the result holds for ℓ. Then, X X X [H, Q]ℓ+1 = [H, [H, Q]ℓ ] = i λP [P, [H, Q]ℓ ] = i λP dT,ℓ [P, T ], (5.143) P
P
T
where in the last equality we use the inductive hypothesis. By the inductive hypothesis, supp(T ) ≤ k(ℓ + 1). Since |SP | ≤ k, |supp([P, T ])| ≤ k(ℓ + 2). We can bound the 1-norm of the Pauli coefficients. Note that λP dT,ℓ is only included when [P, T ] ̸= 0. X XX |dT,ℓ+1 | ≤ 2 |λP ||dT,ℓ |JT overlaps with P K (5.144) T
P
T
42
Let supp(T ) = {i1 , . . . , ik(ℓ+1) } (the argument also works for fewer qubits). Then, we can break up the sum over P according to whether P is supported on some qubit ij . X X X X |dT,ℓ+1 | ≤ 2 |dT,ℓ | |λP | + · · · + |λP | (5.145) T
P :SP ∋i1
T
≤ 2gk(ℓ + 1)
X
P :SP ∋ik(ℓ+1)
|dT,ℓ |
(5.146)
≤ 2gk(ℓ + 1)ℓ!(2kg)ℓ
(5.147)
T
= (ℓ + 1)!(2kg)
ℓ+1
.
(5.148)
In the second line, we use that ∥λ∥B1 ≤ g. In the third line, we use the inductive hypothesis. We also find the following two claims useful. The first claim provides an explicit expression for the derivative of a nested commutator. Claim 5.16. Let Q, T ∈ Pn . For all ℓ ≥ 1 and x ∈ [−1, 1]m , ∂xQ [H(x), T ]ℓ = i
ℓ−1 X
[H(x), [Q, [H(x), T ]ℓ−1−j ]]j .
(5.149)
j=0
Proof. This is essentially the Leibniz rule, but we prove it for completeness. We proceed via induction. For the base case of ℓ = 1, ∂xQ [H(x), T ] = [∂xQ H(x), T ] = [Q, T ], which is clearly the same as the right-hand side (only the j = 0 term survives and ℓ − 1 − j = 0). For the inductive case, suppose the claim holds for ℓ. ∂xQ [H(x), T ]ℓ+1 = ∂xQ [H(x), [H(x), T ]ℓ ]
(5.150)
= [∂xQ H(x), [H(x), T ]ℓ ] + [H(x), ∂xQ [H(x), T ]ℓ ] = i[Q, [H(x), T ]ℓ ] + i
ℓ−1 X
(5.151)
[H(x), [H(x), [Q, [H(x), T ]ℓ−1−j ]]j ]
(5.152)
[H(x), [Q, [H(x), T ]ℓ−j ]]j
(5.153)
j=0
= i[Q, [H(x), T ]ℓ ] + i
ℓ X j=1
ℓ X =i [H(x), [Q, [H(x), T ]ℓ−j ]]j .
(5.154)
j=0
Here, the second line follows from the Leibniz rule. The third line follows by the inductive hypothesis. The fourth line follows by shifting the index j → j + 1. Claim 5.17. Let E1 , . . . , EM ∈ Pn be pairwise distinct Paulis. Then, for any observables A, B, M X
tr(A[Eb , B]) ≤ 2∥A∥P,1 ∥B∥P,1 ,
(5.155)
b=1
where, for A =
P
Proof. Write A =
R∈Pn αR R, ∥A∥P,1 =
P
P
R∈Pn αR R and B = M X
R |αR |.
P
S∈Pn βS S. Then, we have
tr(A[Eb , B]) ≤
X
|αR ||βS |
R,S
b=1
M X
tr(R[Eb , S]) .
(5.156)
b=1
Thus, it suffices to show that, for any Paulis R, S, M X
tr(R[Eb , S]) ≤ 2.
b=1
43
(5.157)
Notice that
( [Eb , S] =
0 2ωT
if [Eb , S] = 0 , if [Eb , S] ̸= 0
where T ∈ Pn such that ωT = Eb S, where ω ∈ {±1, ±i}. Then, ( 0 if [Eb , S] = 0 tr(R[Eb , S]) = . 2ω1{T = R} if [Eb , S] ̸= 0
(5.158)
(5.159)
We claim that there can be only one b ∈ [M ] such that [Eb , S] ̸= 0 and T = R. This means that only one b contributes to the sum over b, so the result follows. This is true because for [Eb , S] ̸= 0 and T = R to hold, Eb = ωRS (since ωT = Eb S when [Eb , S] ̸= 0). Thus, given R and S, Eb is uniquely determined amongst the set of pairwise distinct (phaseless) Paulis. With all of these lemmas, we can now prove Lemma 5.14. Proof of Lemma 5.14. Let x be such that ∥x∥B1 ≤ g. Recall by the definition of FP , 1 1 EP (x) − EP (λ) t t 1 1 = E [tr(e−iHt ReiHt RP )] − EP (λ) t R∼PSP t 1 1 = E [tr(ReiHt RP e−iHt )] − EP (λ) t R∼PSP t ∞ X 1 1 (it)ℓ = − E [tr(R[H(x), RP ])] + E [tr(R[H(x), RP ]ℓ )] − EP (λ) R∼PSP R∼P t ℓ! t SP
FP (x) =
(5.160) (5.161) (5.162) (5.163)
ℓ=2
∞
1 X (it)ℓ 1 E [tr(R[H(x), RP ]ℓ )] − EP (λ) (5.164) R∼PSP t ℓ! R∼PSP t Q ℓ=2 ∞ X 1 X (it)ℓ 1 =− E [tr(R[H(x), RP ]ℓ )] − EP (λ) xQ E [tr(RQRP )] − E [tr(P Q)] + R∼PSP R∼P R∼PSP t ℓ! t SP =−
X
xQ
E
[tr(R[Q, RP ])] +
Q
ℓ=2
(5.165) = xP +
∞ 1 X (it)ℓ
t
ℓ=2
ℓ!
E
R∼PSP
1 [tr(R[H(x), RP ]ℓ )] − EP (λ), t
(5.166)
where in the last line, we use Lemma 2.8. We highlight that, unlike in the Lindbladian setting, there is no confusion (in the sense of Sections 2.2.1 and 3). This is what makes the analysis of the algorithm much simpler because the A matrix from Notation 3.2 is now just the identity matrix. Then, the Jacobian is given by ∞ 1 X (it)ℓ JP,Q (x) = ∂xQ FP (x) = δP,Q + E [∂x tr(R[H(x), RP ]ℓ )]. (5.167) t ℓ! R∼PSP Q ℓ=2
Thus, the leading order term of the Jacobian is I. It suffices to show that the higher order terms are small, i.e., ∥J(x) − I∥∞→∞ ≤ cg t. Recall that ∥M ∥∞→∞ = sup∥x∥∞ ≤1 ∥M x∥∞ . Consider u = (u1 , . . . , uM ) such that ∥u∥∞ ≤ 1. We will bound ∥(J(x) − I)u∥∞ . |((J(x) − I)u)P | =
X
(J(x) − I)P,Q uQ
(5.168)
Q
! ∞ X (it)ℓ 1 X uQ E [∂x tr(R[H(x), RP ]ℓ )] = t ℓ! R∼PSP Q Q
ℓ=2
44
(5.169)
∞ 1 X tℓ X ≤ E [ ∂xQ tr(R[H(x), RP ]ℓ ) ], R∼PSP t ℓ!
(5.170)
Q
ℓ=2
where in the second line we use Equation (5.167), and in the third line, we use |uQ | ≤ 1 and the triangle inequality. Now, using Claim 5.16, we can expand the derivative of the nested commutator: ∞
|((J(x) − I)u)P | ≤
ℓ−1
1 X tℓ X X E tr(R[H(x), [Q, [H(x), RP ]ℓ−1−j ]]j ) R∼PSP t ℓ! j=0
(5.171)
Q
ℓ=2
To get this into the correct form to apply Claim 5.17, we can repeatedly use tr(A[H, B]) = tr([A, H]B). Thus, we have ∞ ℓ−1 X 1 X tℓ X |((J(x) − I)u)P | ≤ tr([R, H(x)]j [Q, [H(x), RP ]ℓ−1−j ]) (5.172) E t ℓ! j=0 R∼PSP Q
ℓ=2
≤
1 t
∞ ℓ X ℓ−1 X ℓ=2 ∞
≤
≤
t E [2∥[R, H(x)]j ∥P,1 ∥[H(x), RP ]ℓ−j−1 ∥P,1 ] ℓ! j=0 R∼PSP
(5.173)
ℓ−1
2 X tℓ X j!(2kg)j (ℓ − 1 − j)!(2kg)ℓ−1−j t ℓ! j=0
(5.174)
2 t
tℓ (2kg)ℓ−1
(5.175)
(2kgt)ℓ
(5.176)
=2
ℓ=2 ∞ X
ℓ=2 ∞ X ℓ=1
2kgt (5.177) 1 − 2kgt 81kg ≤t· . (5.178) 20 In the second line, we use Claim 5.17. In the third line, we use Lemma 5.15, which can be applied because ∥x∥B1 ≤ g. In the fourth line, we use that =2·
ℓ−1 X
j!(ℓ − 1 − j)! =
j=0
ℓ−1 X
ℓ−1
X 1 (ℓ − 1)! ℓ−1 ≤ (ℓ − 1)! 1 = ℓ!. j
j=0
(5.179)
j=0
In the fifth line, we shift the index ℓ → ℓ − 1. In the sixth line, we use the sum of a geometric series when 2kgt < 1. Finally, in the last line, we use t < 1/(162kg). Thus, this implies that ∥J(x) − I∥∞→∞ ≤ cg t. 5.6.2
High-temperature Gibbs states
We can also apply Theorem 5.11 to the problem of structure learning a bounded degree Hamiltonian from copies of its Gibbs state. This is the first algorithm for structure learning Hamiltonians from any Gibbs state access model. Again, we prove this by showing that the hypotheses of Theorem 5.11 hold in this setting. Theorem 5.18 (Structure P learning of Hamiltonians from high-temperature Gibbs states). Let ε, δ > 0, and let k = O(1). Let H = i P λP P be a k-local Hamiltonian with ∥λ∥∞ ≤ 1 and bounded-degree interactions, i.e., every qubit interacts with at most d nonzero terms. Let β > 0 satisfy β≤
1 . 1000e6 (2kd + 1)8
(5.180)
b such that ∥λ b − λ∥∞ ≤ ε with probability at least Then, there exists an algorithm that finds estimates λ 2 −βH 1 − δ using O(log(n/δ)/(βε) ) copies of the Gibbs state ρβ = e / tr(e−βH ). The classical runtime of this k 2 algorithm is O(n poly(d) log(n/δ)/(βε) ). 45
Before proving this theorem, we recall a result from [HKT24b], which bounds the ∞ → ∞ norm of the Jacobian for high-temperature Gibbs states. Lemma 5.19 (Lemma 4.3 in [HKT24b]). Let JP,Q (x) = ∂xQ FP (x) be the Jacobian, where P, Q range over k-local Paulis such that the graph with these Paulis as vertices, and edges between P, Q when SP ∩ SQ ̸= ∅, has degree at most D. Then, ∥J(x) − I∥∞→∞ ≤ cD β, (5.181) where cD = 50e6 (D + 1)8 . We use this to instantiate the hypothesis on the Jacobian bound in Theorem 5.11. Thus, with this, we can prove Theorem 5.18. bP (λ) such that Proof of Theorem 5.18. As with Theorem 5.13, we first need to obtain estimates E bP (λ) − EP (λ)| ≤ |E
βε 40
(5.182)
for all k-local Paulis P . We can do so with probability at least 1 − δ using O(4k log(n/δ)/(βε)2 ) copies of the Gibbs state ρβ via classical shadows [HKP20]. Moreover, for k = O(1), this has a classical runtime of O(nk log(n/δ)/(βε)2 ). For our choice of β, this corresponds to O(d16 log(n/δ)/ε2 ) copies. bP (x) such that Also, we can obtain estimates E bP (x) − EP (x)| ≤ |E
βε 40
(5.183)
by computing a truncated cluster expansion. The proof of Theorem 4.6 in [HKT24b] shows that this takes classical runtime O(nk poly(d, log(1/(βε))/ε). This applies to our setting because we maintain the invariant in (j) our algorithm that |xP | ≤ 2|λP |. This means that deg(x) ≤ deg(λ) ≤ d, which fits in the setting of [HKT24b]. This gives us an approximation of F(x) up to ε/20 error. It remains to argue that the Jacobian bound holds. Consider the statement we want to prove: X |JP,Q (x) − δP,Q | ≤ cd β, (5.184) Q∈SH
for deg(x) ≤ d. This closely resembles Lemma 5.19. Notice that Q ∈ SH , so since H has bounded degree, the graph with Q ∈ SH as its vertices has degree at most d. However, the key difference is that P can be any k-local Pauli and does not necessarily lie in this graph. Nevertheless, the graph with vertices consisting of SH ∪ {P } still has degree at most D = kd. Thus, we can apply Lemma 5.19 with D = kd to obtain our claim with cd = 50e6 (kd + 1)8 . This completes the proof.
Acknowledgments L.L. thanks Aram Harrow, Jordi Montana-Lopez, Quynh Nguyen, Umesh Vazirani, and Thomas Vidick for helpful discussions. L.L. is supported by the U.S. Department of Energy, Office of Science, Office of Advanced Scientific Computing Research, Department of Energy Computational Science Graduate Fellowship under Award Number DE-SC0026073. E.T. is supported by the Miller Institute for Basic Research in Science, University of California Berkeley. J.W. is supported by the NSF CAREER award CCF-233971 and a Sloan Fellowship. This report was prepared as an account of work sponsored by an agency of the United States Government. Neither the United States Government nor any agency thereof, nor any of their employees, makes any warranty, express or implied, or assumes any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that its use would not infringe privately owned rights. Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise does not necessarily constitute or imply its endorsement, recommendation, or favoring by the United States Government or any agency thereof. The views and opinions of authors expressed herein do not necessarily state or reflect those of the United States Government or any agency thereof. 46
References [AAKS21]
Anurag Anshu, Srinivasan Arunachalam, Tomotaka Kuwahara, and Mehdi Soleimanifar. “Sample-efficient learning of interacting quantum systems”. In: Nature Physics (2021) (page 7).
[Aba+25]
Dmitry Abanin, Rajeev Acharya, Laleh Aghababaie-Beni, Georg Aigeldinger, Ashok Ajoy, Ross Alcaraz, Igor Aleiner, Trond Anderson, Markus Ansmann, Frank Arute, et al. “Observation of constructive interference at the edge of quantum ergodicity”. In: Nature (2025) (page 2).
[ACGGMS25]
Amira Abbas, Nunzia Cerrato, Francisco Escudero Gutiérrez, Dmitry Grinko, Francesco Anna Mele, and Pulkit Sinha. “Nearly optimal algorithms to learn sparse quantum Hamiltonians in physically motivated distances”. In: arXiv preprint arXiv:2509.09813 (2025) (page 7).
[ACGRY26]
Itai Arad, Zhili Chen, Naixu Guo, Patrick Rebentrost, and Zhan Yu. “Near-Optimal Learning of Local Lindbladians”. In: arXiv preprint arXiv:2606.20535 (2026) (pages 8, 9).
[ADG24]
Srinivasan Arunachalam, Arkopal Dutt, and Francisco Escudero Gutiérrez. “Testing and learning structured quantum Hamiltonians”. In: arXiv preprint arXiv:2411.00082 (2024) (page 7).
[AKL16]
Itai Arad, Tomotaka Kuwahara, and Zeph Landau. “Connecting global and local energy distributions in quantum spin models on a lattice”. In: Journal of Statistical Mechanics: Theory and Experiment (2016). doi: 10.1088/1742-5468/2016/03/033301. arXiv: 1406.3898 [quant-ph] (page 3).
[Alh23]
Álvaro M Alhambra. “Quantum many-body systems in thermal equilibrium”. In: PRX Quantum (2023) (page 3).
[Aru+19]
Frank Arute, Kunal Arya, Ryan Babbush, Dave Bacon, Joseph C Bardin, Rami Barends, Rupak Biswas, Sergio Boixo, Fernando GSL Brandao, David A Buell, et al. “Quantum supremacy using a programmable superconducting processor”. In: Nature (2019) (page 2).
[BAL19]
Eyal Bairey, Itai Arad, and Netanel H Lindner. “Learning a local Hamiltonian from local measurements”. In: Physical review letters (2019) (page 7).
[BC26]
Thiago Bergamaschi and Chi-Fang Chen. “Fast Mixing of Quantum Spin Chains at All Temperatures”. In: Proceedings of the 58th Annual ACM Symposium on Theory of Computing. 2026 (pages 2, 36).
[BCGOR25]
Andreas Bluhm, Matthias C Caro, Francisco Escudero Gutiérrez, Aadil Oufkir, and Cambyse Rouzé. “Certifying and learning quantum Ising Hamiltonians”. In: arXiv preprint arXiv:2509.10239 (2025) (page 7).
[BCL24]
Thiago Bergamaschi, Chi-Fang Chen, and Yunchao Liu. “Quantum computational advantage with constant-temperature Gibbs sampling”. In: 2024 IEEE 65th Annual Symposium on Foundations of Computer Science (FOCS). IEEE. 2024 (pages 2, 36).
[BGPLA20]
Eyal Bairey, Chu Guo, Dario Poletti, Netanel H Lindner, and Itai Arad. “Learning the dynamics of open quantum systems from their steady states”. In: New Journal of Physics (2020). doi: 10.1088/1367-2630/ab73cd (page 8).
[Bir+26]
Rune Thinggaard Birke, Johann Bock Severin, Malthe A Marciniak, Emil Hogedal, Andreas Nylander, Irshad Ahmad, Amr Osman, Janka Biznarova, Marcus Rommel, Anita Fadavi Roudsari, et al. “Demonstrating and benchmarking classical shadows for lindblad tomography”. In: arXiv preprint arXiv:2602.14694 (2026) (page 8).
[BLMT24a]
Ainesh Bakshi, Allen Liu, Ankur Moitra, and Ewin Tang. “High-temperature Gibbs states are unentangled and efficiently preparable”. In: 2024 IEEE 65th Annual Symposium on Foundations of Computer Science (FOCS). IEEE. 2024 (pages 2, 36).
[BLMT24b]
Ainesh Bakshi, Allen Liu, Ankur Moitra, and Ewin Tang. “Learning quantum Hamiltonians at any temperature in polynomial time”. In: Proceedings of the 56th Annual ACM Symposium on Theory of Computing. 2024 (pages 7, 9). 47
[BLMT24c]
Ainesh Bakshi, Allen Liu, Ankur Moitra, and Ewin Tang. “Structure learning of Hamiltonians from real-time evolution”. In: 2024 IEEE 65th Annual Symposium on Foundations of Computer Science (FOCS). IEEE. 2024 (pages 3–7, 9, 11, 36).
[BLMT26]
Ainesh Bakshi, Allen Liu, Ankur Moitra, and Ewin Tang. “A Dobrushin condition for quantum Markov chains: Rapid mixing and conditional mutual information at high temperature”. In: Proceedings of the 58th Annual ACM Symposium on Theory of Computing. 2026 (pages 2, 36).
[BMWM25]
Ewout van den Berg, Brad Mitchell, Ken Xuan Wei, and Moein Malekakhlagh. “Largescale Lindblad learning from time-series data”. In: arXiv preprint arXiv:2512.08165 (2025) (page 8).
[Buž98]
Vladimír Bužek. “Reconstruction of Liouvillian superoperators”. In: Physical Review A (1998) (page 8).
[BY23]
Zongbo Bao and Penghui Yao. “On testing and learning quantum junta channels”. In: The Thirty Sixth Annual Conference on Learning Theory. PMLR. 2023 (page 11).
[BZMT26]
Shrigyan Brahmachari, Shuchen Zhu, Iman Marvian, and Yu Tong. “Learning Hamiltonians in the Heisenberg limit with static single-qubit fields”. In: arXiv preprint arXiv:2601.10380 (2026) (page 7).
[CAN25]
Chi-Fang Chen, Anurag Anshu, and Quynh T. Nguyen. “Learning quantum Gibbs states locally and efficiently”. In: arXiv preprint arXiv:2504.02706 (2025) (page 7).
[Car24]
Matthias C Caro. “Learning quantum processes and Hamiltonians via the Pauli transfer matrix”. In: ACM Transactions on Quantum Computing (2024) (page 7).
[CCH25]
Sitan Chen, Jordan Cotler, and Hsin-Yuan Huang. “Quantum probe tomography”. In: arXiv preprint arXiv:2510.08499 (2025) (page 7).
[CCH26]
Constantin Cedillo Vayson de Pradenne, Jordan Cotler, and Hsin-Yuan Huang. “Learning Hamiltonians at long times”. In: arXiv preprint arXiv:2606.05690 (2026) (page 7).
[CKBG25]
Chi-Fang Chen, Michael Kastoryano, Fernando GSL Brandão, and András Gilyén. “Efficient quantum thermal simulation”. In: Nature (2025) (pages 2, 4, 36).
[CKG23]
Chi-Fang Chen, Michael J Kastoryano, and András Gilyén. “An efficient and exact noncommutative quantum Gibbs sampler”. In: arXiv preprint arXiv:2311.09207 (2023) (pages 2, 36).
[CLS25]
Ziyun Chen, Jerry Li, and Joseph Slote. “Lower Bounds for Learning Hamiltonians from Time Evolution”. In: arXiv preprint arXiv:2509.20665 (2025) (page 7).
[CW23]
Juan Castaneda and Nathan Wiebe. “Hamiltonian learning via shadow tomography of pseudo-choi states”. In: arXiv preprint arXiv:2308.13020 (2023) (page 7).
[DDMPRT23]
Nicolo Defenu, Tobias Donner, Tommaso Macri, Guido Pagano, Stefano Ruffo, and Andrea Trombettoni. “Long-range interacting quantum systems”. In: Reviews of Modern Physics (2023) (page 8).
[DLL24]
Zhiyan Ding, Bowen Li, and Lin Lin. “Efficient quantum Gibbs samplers with Kubo–Martin– Schwinger detailed balance condition”. In: arXiv preprint arXiv:2404.05998 (2024) (pages 2, 4, 36).
[DOS24]
Alicja Dutkiewicz, Thomas E O’Brien, and Thomas Schuster. “The advantage of quantum control in many-body Hamiltonian learning”. In: Quantum (2024) (page 7).
[GKS76]
Vittorio Gorini, Andrzej Kossakowski, and Ennackal Chandy George Sudarshan. “Completely positive dynamical semigroups of N-level systems”. In: Journal of Mathematical Physics (1976) (page 2).
[Gut24]
Francisco Escudero Gutiérrez. “Simple algorithms to test and learn local Hamiltonians”. In: arXiv preprint arXiv:2404.06282 (2024) (page 7).
48
[HKP20]
Hsin-Yuan Huang, Richard Kueng, and John Preskill. “Predicting many properties of a quantum system from very few measurements”. In: Nature Physics (2020) (page 46).
[HKT24a]
Jeongwan Haah, Robin Kothari, and Ewin Tang. “Learning quantum Hamiltonians from high-temperature Gibbs states and real-time evolutions”. In: Nature Physics (2024). doi: 10.1038/s41567-023-02376-x. arXiv: 2108.04842 [quant-ph] (page 5).
[HKT24b]
Jeongwan Haah, Robin Kothari, and Ewin Tang. “Learning quantum Hamiltonians from high-temperature Gibbs states and real-time evolutions”. In: Nature Physics (2024) (pages 4– 8, 28, 34, 35, 46).
[HTFS23]
Hsin-Yuan Huang, Yu Tong, Di Fang, and Yuan Su. “Learning many-body Hamiltonians with Heisenberg-limited scaling”. In: Physical Review Letters (2023) (page 7).
[Hu+25]
Hong-Ye Hu, Muzhou Ma, Weiyuan Gong, Qi Ye, Yu Tong, Steven T Flammia, and Susanne F Yelin. “Ansatz-free Hamiltonian learning with Heisenberg-limited scaling”. In: arXiv preprint arXiv:2502.11900 (2025) (page 7).
[IRGGHY26]
Petr Ivashkov, Nikita Romanov, Weiyuan Gong, Andi Gu, Hong-Ye Hu, and Susanne F Yelin. “Ansatz-Free Learning of Lindbladian Dynamics In Situ”. In: arXiv preprint arXiv:2603.05492 (2026) (pages 4, 5, 8).
[Kra+25]
Tristan Kraft, Manoj K Joshi, William Lam, Tobias Olsacher, Florian Kranzl, Johannes Franke, Lata Kh Joshi, Rainer Blatt, Augusto Smerzi, Daniel Stilck França, et al. “BoundedError Quantum Simulation via Hamiltonian and Lindbladian Learning”. In: arXiv preprint arXiv:2511.23392 (2025) (page 8).
[Lin76]
Goran Lindblad. “On the generators of quantum dynamical semigroups”. In: Communications in mathematical physics (1976) (page 2).
[LJFV26]
William T Lam, Manoj K Joshi, Daniel Stilck França, and Benoit Vermersch. “Pairwise Liouvillian learning from randomized measurements: practical aspects and guidelines for operating the protocol in large-scale experiments”. In: arXiv preprint arXiv:2605.26953 (2026) (page 8).
[LSKOC25]
Yinchen Liu, James R Seddon, Tamara Kohler, Emilio Onorati, and Toby S Cubitt. “Robust Lindbladian Estimation for Quantum Dynamics”. In: arXiv preprint arXiv:2507.07912 (2025) (page 8).
[LTGNY24]
Haoya Li, Yu Tong, Tuvia Gefen, Hongkang Ni, and Lexing Ying. “Heisenberg-limited Hamiltonian learning for interacting bosons”. In: npj Quantum Information (2024) (page 7).
[MBFR26]
Tim Mobus, Thiago Bergamaschi, Daniel Stilck França, and Cambyse Rouze. “Robust Structure Learning of k-local Lindbladians”. In: arXiv preprint arXiv:2606.20706 (2026) (pages 8, 9).
[MBGTWR25]
Tim Möbus, Andreas Bluhm, Tuvia Gefen, Yu Tong, Albert H. Werner, and Cambyse Rouzé. “Heisenberg-limited Hamiltonian learning continuous variable systems via engineered dissipation”. In: arXiv preprint arXIv:2506.00606 (2025) (page 7).
[MECT25]
Jordi A. Montana-Lopez, Andreas Elben, Joonhee Choi, and Rahul Trivedi. “Efficiently learning non-Markovian noise in many-body quantum simulators”. In: arXiv preprint arXiv:2511.16772 (2025) (pages 4, 5, 8).
[MFPT24]
Muzhou Ma, Steven T Flammia, John Preskill, and Yu Tong. “Learning k-body Hamiltonians via compressed sensing”. In: arXiv preprint arXiv:2410.18928 (2024) (page 7).
[MH24]
Arjun Mirani and Patrick Hayden. “Learning interacting fermionic Hamiltonians at the Heisenberg limit”. In: Physical Review A (2024) (page 7).
[MM24]
Ryan L Mann and Romy M Minko. “Algorithmic cluster expansions for quantum problems”. In: PRX Quantum (2024) (page 22).
[NLY24]
Hongkang Ni, Haoya Li, and Lexing Ying. “Quantum hamiltonian learning for the fermihubbard model”. In: Acta Applicandae Mathematicae (2024) (page 7).
49
[OKC23]
Emilio Onorati, Tamara Kohler, and Toby S Cubitt. “Fitting quantum noise models to tomography data”. In: Quantum (2023) (page 8).
[OKKKZ25]
Tobias Olsacher, Tristan Kraft, Christian Kokail, Barbara Kraus, and Peter Zoller. “Hamiltonian and Liouvillian learning in weakly-dissipative quantum many-body systems”. In: Quantum Science and Technology (2025) (page 8).
[OKSM24]
Tatsuki Odake, Hlér Kristjánsson, Akihito Soeda, and Mio Murao. “Higher-order quantum transformations of Hamiltonian dynamics”. In: Physical Review Research (2024) (page 7).
[POKZ22]
Lorenzo Pastori, Tobias Olsacher, Christian Kokail, and Peter Zoller. “Characterization and verification of Trotterized digital quantum simulation via Hamiltonian and Liouvillian learning”. In: PRX Quantum (2022) (page 8).
[QR19]
Xiao-Liang Qi and Daniel Ranard. “Determining a local Hamiltonian from a single eigenstate”. In: Quantum (2019) (page 7).
[Rom+26]
Nikita Romanov, Petr Ivashkov, Weiyuan Gong, Ishaan Kannan, Andi Gu, Hong-Ye Hu, and Susanne F Yelin. “Learning Arbitrary Lindbladians with Quantum Error Correction”. In: arXiv preprint arXiv:2606.18188 (2026) (page 8).
[RSA26]
Cambyse Rouzé, Daniel Stilck França, and Álvaro M Alhambra. “Optimal quantum algorithm for Gibbs state preparation”. In: Physical Review Letters (2026) (pages 2, 36).
[SA26]
Matteo Scandi and Álvaro M Alhambra. “Thermalization in open many-body systems and KMS detailed balance”. In: Physical Review X (2026) (pages 2, 4, 36).
[Sin26]
Savar Sinha. “Efficient and SPAM-Robust Ansatz-Free Lindbladian Learning”. In: arXiv preprint arXiv:2606.20706 (2026) (page 8).
[SLO26]
Myeongjin Shin, Junseo Lee, and Changhun Oh. “Heisenberg-limited Hamiltonian learning without short-time control”. In: arXiv preprint arXiv:2604.27838 (2026) (page 7).
[SLP11]
Marcus P da Silva, Olivier Landon-Cardinal, and David Poulin. “Practical characterization of quantum devices without tomography”. In: Physical Review Letters (2011) (page 8).
[SMDWB24]
Daniel Stilck França, Liubov A Markovich, VV Dobrovitski, Albert H Werner, and Johannes Borregaard. “Efficient and robust estimation of many-qubit Hamiltonians”. In: Nature Communications (2024) (pages 4, 5, 8).
[SMRW25]
Daniel Stilck Franca, Tim Mobus, Cambyse Rouze, and Albert H. Werner. “Learning and certification of local time-dependent quantum dynamics and noise”. In: arXiv preprint arXiv:2510.08500 (2025) (pages 4, 5, 8).
[ST25]
Savar D Sinha and Yu Tong. “Improved Hamiltonian learning and sparsity testing through Bell sampling”. In: arXiv preprint arXiv:2509.07937 (2025) (page 7).
[TW25]
Ewin Tang and John Wright. “Amplitude amplification and estimation require inverses”. Technical report, arXiv:2507.23787. 2025 (page 9).
[Wu+21]
Yulin Wu, Wan-Su Bao, Sirui Cao, Fusheng Chen, Ming-Cheng Chen, Xiawei Chen, TungHsun Chung, Hui Deng, Yajie Du, Daojin Fan, et al. “Strong quantum computational advantage using a superconducting quantum processor”. In: Physical review letters (2021) (page 2).
[Zha24]
Andrew Zhao. “Learning the structure of any Hamiltonian from minimal assumptions”. In: arXiv preprint arXiv:2410.21635 (2024) (page 7).
[Zho+20]
Han-Sen Zhong, Hui Wang, Yu-Hao Deng, Ming-Cheng Chen, Li-Chao Peng, Yi-Han Luo, Jian Qin, Dian Wu, Xing Ding, Yi Hu, et al. “Quantum computational advantage using photons”. In: Science (2020) (page 2).
[Zhu+23]
Daiwei Zhu, Gregory D Kahanamoku-Meyer, Laura Lewis, Crystal Noel, Or Katz, Bahaa Harraz, Qingfeng Wang, Andrew Risinger, Lei Feng, Debopriyo Biswas, et al. “Interactive cryptographic proofs of quantumness using mid-circuit measurements”. In: Nature Physics (2023) (page 2).
50
[ZYLB21]
Assaf Zubida, Elad Yitzhaki, Netanel H Lindner, and Eyal Bairey. “Optimal short-time measurements for Hamiltonian learning”. In: arXiv preprint arXiv:2108.08824 (2021) (page 7).
51