With the goal to study dissipative Landau–Zener (LZ) sweeps in realistic solid-state qubits, we utilize novel methods from non-Markovian open quantum system dynamics that enable reliable long-time simulations for sub-Ohmic environments. In particular, we combine a novel representation of the dynamical propagator, the uniform time evolving matrix product operator method, with a stochastic realization of finite temperature fluctuations. The latter greatly reduces the computational cost for the matrix product operator approach, enabling convergence in the experimentally relevant deeply sub-Ohmic regime. Our method allows the exact simulation of dynamical protocols with long operation times, such as the LZ sweep, in challenging parameter regimes that are realized in current experimental platforms.
I. INTRODUCTION
Quantum dynamics are profoundly affected by environmental influences, leading to damping and coherence loss, which pose significant challenges for the design of scalable quantum hardware.1 Despite advancements, the limitations of current noisy intermediate-scale quantum devices have led to a shift in research focus from universal quantum computing toward more attainable goals, such as quantum simulation and annealing.2–4 These applications still require a precise control over individual quantum systems, necessitating a comprehensive theoretical understanding of the underlying noise processes. For instance, in quantum annealers, the dissipation poses a limit on the adiabaticity constraint, leading to compromises in the final state fidelity.2,5,6 However, in many currently available solid-state quantum devices, a thorough characterization of dissipation remains difficult, and the microscopic origins of the noise largely remain unclear.7
Previous exact numerical simulations of the dissipative Landau–Zener problem relied on the quasi-adiabatic path integral (QUAPI) approach and were restricted to environments with relatively short memory times. However, in recent experimental realizations of the Landau–Zener transition with flux qubits, strong low frequency noise contributions7,18,19 lead to highly non-Markovian dynamics for which the QUAPI method becomes unfeasible and master equation approaches become inaccurate. Here, we demonstrate that exact numerical simulations of such hard open system problems can be achieved by combining novel simulation techniques based on matrix product operator (MPO) methods20–31 with a stochastic realization of thermal fluctuations.32–34 This powerful numerical tool can be used in future studies to realize precise agreements of qubit dynamics in real world quantum devices with simulations and to deliver important insights into the characteristics of dissipation.
This paper is structured as follows: In Sec. II, we provide a brief review of open system dynamics and introduce the numerical tools that we use for our simulations. In Sec. III, we present our main results, exact simulations of the Landau–Zener transition in dissipative environments with very long memory times. In Sec. IV, we provide detailed introductions to the numerical methods and their practical implementation. Finally, we close with our conclusions and a review on future perspectives.
II. OPEN SYSTEM DYNAMICS
The crucial numerical bottleneck for the uniTEMPO computation is the final bond dimension of the MPO influence functional, as the computation time scales with the third power of this bond dimension. A long memory of the bath together with a strong system–bath coupling will require large bond dimensions, thus limiting the applicability of the approach. As demonstrated below, this problem arises in particular for strongly sub-Ohmic environments at high temperatures. Environments with such strong low-frequency contributions occur generically in solid state quantum devices. Thus, exact simulation techniques that can faithfully describe the open system evolution for large times are highly desired.
We demonstrate in the following that such difficult environments can still be handled numerically by combining the uniTEMPO method for zero temperature with an exact realization of finite temperature fluctuations in terms of stochastic unitary evolution. TEMPO approaches are particularly well suited for this stochastic approach because the main computational step, obtaining a MPO form for the influence functional, has to be performed only once. As has been observed in Ref. 21, when the MPO is available, the open system evolution can be computed very efficiently for different time-dependent system Hamiltonians.
To exemplify this issue, we performed an example calculation for the spin boson model with different sub-Ohmic environments. We computed the MPO form for the influence functional obtained with uniTEMPO for a fixed compression tolerance. As can be seen in Fig. 1 (lower panel), the bond dimension of the influence functional MPO increases drastically as the bath temperature increases. This effect is particularly dramatic in the highly sub-Ohmic case (low s) because the tail of the thermal bath correlation function becomes extremely heavy. Note again that the computational demand scales with the third power of the bond dimension.
III. LANDAU–ZENER TRANSITION
The LZ probability P is typically considered as a function of the sweep velocity v. For the particular basis chosen in (2), it describes the probability for reaching the eigenstate at large positive times when starting in the eigenstate at large negative times. In addition to this so-called diabatic basis, the instantaneous eigenstates of the system Hamiltonian H(t), |+(t)⟩, |−(t)⟩ (adiabatic basis), can be considered at any time during the sweep. Note that the lower energy adiabatic state |−(t)⟩ of the system approaches the diabatic state |1⟩ for large negative times, but the state |0⟩ for positive times. Due to the avoided level crossing, the eigenenergies corresponding to the adiabatic basis states are always separated by at least the gap Δ. Having introduced these two bases, the (ground-state) LZ-probability can be expressed either as the probability of the system to stay in the same adiabatic state it was prepared in, , or equivalently as the probability of a transition between the two diabatic basis states, . For a perfectly isolated qubit, this probability is close to one for slow sweeps, as the system is able to evolve adiabatically and always stays in a momentary eigenstate. As the sweep velocity is increased, the LZ probability decreases monotonically, reaching zero for very fast sweeps. Adding the influence of an environment to the system evolution, the probability may deviate strongly from the analytical solution. In the case of dissipative dynamics, the LZ probability is influenced in a highly nontrivial way and can display a nonmonotonic behavior.18,44
As a first example, we revisit the case of an Ohmic bath s = 1 for which previous results using the QUAPI method are available.44 For weak to intermediate coupling, the Ohmic case is computationally rather simple because the bath correlation function decays quickly. Our results are displayed in Fig. 2. One can see a nonmonotonic behavior of the LZ curves, which can be explained in terms of a competition between dissipation and thermal fluctuations. At fast sweep velocities, the effective interaction time with the environment is too short to cause deviations from the unitary behavior. At intermediate sweep velocities, thermal fluctuations are the main influence. These fluctuations allow excitation to the excited state and hence cause the LZ probability to decrease. If the sweep is slowed down even further, the dissipation causes a relaxation to the lower energy state, leading to a subsequent increase in the ground state probability. The dissipation due to the environment decreases with decreasing temperature. In fact, at zero temperature, it was shown analytically that the LZ probability is unaffected by the bath in the case of purely longitudinal coupling .8
In Fig. 2, we compare our simulations with QUAPI results from Ref. 44. We find slight deviations of our converged uniTEMPO calculations with the QUAPI data. We suspect that these deviations arise from the Trotterization used in the QUAPI simulations, where a fairly large time step δ = 0.1 was chosen to cope with the exponential scaling of the numerical effort in QUAPI. We reach close agreement with the QUAPI results using this large time step in uniTEMPO. Clearly, a smaller value for δ is needed for numerical convergence. For the Ohmic bath, it is still feasible to use the bare uniTEMPO method without a stochastic realization of temperature. In this case, the full thermal correlation function is decaying fast and the required uniTEMPO bond dimension remains small even at high temperatures (see Fig. 1).
We now turn to the challenging sub-Ohmic (s < 1) regime where the bath memory becomes large. In this case, we employ the stochastic approach for finite temperature fluctuations in order to still handle the extremely slow decaying thermal part of the bath correlation function. The additional zero-temperature “quantum” part of the bath can be captured efficiently by uniTEMPO. The drawback of this approach is that, at high temperatures, many trajectories are required to reduce the statistical error of the data points. In Figs. 3–5, these statistical errors are smaller than the symbol size and are therefore not displayed. Note that for long sweep times, computing the dynamics of single realizations takes long computation times, even though parallelization is trivial.
As concrete examples, we calculated the dissipative Landau–Zener problem for s = 0.5 in Fig. 3 and s = 0.1 in Fig. 4 at different coupling strengths and temperatures. Qualitatively, we find a similar behavior as in the Ohmic case in Fig. 2. Moreover, we find that the Landau–Zener probability is not very sensitive to the exact value of s, i.e., the low frequency behavior, if the values for the coupling strength α used lead to a comparable influence of the dissipation [as chosen in Figs. 3 and 4, panels (a)–(c)]. For s = 0.1, thermal fluctuations are generally more dominant at low frequencies, leading to a lower probability for reaching the ground state at low sweep velocities compared with s = 0.5.
In the example calculation displayed in Fig. 5, we consider a strongly sub-Ohmic bath at finite temperatures coupled via different mixing angles. Similar to the Ohmic case considered in Ref. 44, we find that transverse coupling leads to a stronger effect on the LZ probability compared with longitudinal coupling. This, in turn, requires a larger bond dimension for uniTEMPO in order to achieve convergence.
Overall, we did not find significant qualitative differences in the behavior of the LZ probability for sub-Ohmic baths compared with Ohmic baths. It should, however, be noted that environments of real quantum devices, such as superconducting flux qubits, are known to have strong low-frequency contributions within a deeply sub-Ohmic regime (the so-called 1/f noise).7,19,45,46 With the combination of stochastic fluctuations and uniTEMPO, we are able to perform exact calculations in this challenging regime over long evolution times with a relatively low computational effort.
In Sec. IV, we provide detailed explanations of the theoretical framework and the numerical methods that were used in the simulations.
IV. METHODS
A. Influence functional theory
B. uniTEMPO representation of the influence functional
For our model (3), the influence functional is well known analytically. Since the operator obeys Gaussian statistics with respect to the Gaussian bath initial state, a cumulant expansion truncates after the second order, which yields a Gaussian form of .32
An important issue that occurs in the case of sub-Ohmic spectral densities is the algebraic decay of the bath correlation function. True algebraic correlations cannot be realized in a MPO with finite bond-dimensions. In order to ensure the stability of the uniTEMPO algorithm, we include a smooth exponential cutoff to the bath correlation function after a large cutoff time. Exponential memory-time cutoffs are standard in common numerical methods for open system evolutions, such as the hierarchical equations of motion (HEOM). QUAPI and the finite PT-TEMPO algorithms use a hard memory cutoff, which, in contrast to the exponential cutoff, can lead to numerical instabilities when used in uniTEMPO. In general, the dynamics of the open system is unaffected by the cutoff if the cutoff time is larger than the considered evolution time. For longer evolution times, convergence of the results with respect to the cutoff should always be checked.
C. Stochastic realization of finite temperature
The derivation above is applicable to any real-valued bath correlation function whenever the coupling to the bath is Hermitian. Moreover, for a non-Hermitian coupling to the bath, an alternative derivation of the stochastic decoupling based on the thermal P-function of the environment can still be applied in the case of both Gaussian34 and spin environments.48,49
V. CONCLUSIONS
In this work, we presented exact numerical simulations of the dissipative Landau–Zener transition probability for strongly non-Markovian sub-Ohmic environments. This setup directly corresponds to solid-state quantum devices, where 1/f noise is the primary source of dissipation.19 Previous simulations of this dissipation over extended periods have been limited by strong memory effects, which preclude the short-time truncation necessary for many numerical open system methods. We overcome these limitations by combining novel MPO-based simulation techniques28 with a long-known but rarely utilized stochastic sampling of thermal fluctuations.32,34 In particular, the new uniTEMPO approach enables the simulation of dynamics with strong dissipation over arbitrarily long times and with very low memory requirements. In addition, it can be easily integrated with the stochastic finite-temperature method, thereby mitigating convergence issues associated with the heavy-tail statistics of thermal noise.
Our framework provides a robust tool for verifying the quantitative agreement between experimental results and theoretical predictions without relying on additional approximations commonly used in deriving master equations. This capability allows a more accurate assessment of the noise characterization in qubit devices over extended periods and may thereby offer deeper insights into the underlying noise processes affecting the qubit dynamics.
The simulation method is, of course, completely independent of the specific driving protocol and can be used to simulate arbitrary dynamics of single or multiple qubits. For instance, one can straightforwardly study transitions in more intricate sweep protocols.52 For many-qubit systems, the numerical effort will scale exponentially with the number of qubits. However, for local dissipation, this exponential scaling could be lifted in favorable cases by using a MPO representation for the qubit state as well, similar to the approaches used in Refs. 25, 53, and 54.
ACKNOWLEDGMENTS
It is a pleasure to thank Xi Dai, Adrian Lupascu, and Johannes Kassel for insightful discussions in connection with this work. All authors are grateful to the Center for Information Services and High Performance Computing [Zentrum für Informationsdienste und Hochleistungsrechnen (ZIH)] at TU Dresden for providing its facilities for high throughput calculations.
AUTHOR DECLARATIONS
Conflict of Interest
The authors have no conflicts to disclose.
Author Contributions
Felix Kahlert: Conceptualization (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal). Valentin Link: Conceptualization (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal). Richard Hartmann: Conceptualization (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal). Walter T. Strunz: Conceptualization (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal).
DATA AVAILABILITY
The data that support the findings of this study are available from the corresponding author upon reasonable request.