Real-time operator evolution in two and three dimensions via sparse Pauli dynamics
Tomislav Begušić Garnet Kin-Lic Chan




gkc1000@gmail.com
Division of Chemistry and Chemical Engineering, California Institute of Technology, Pasadena, California 91125, USA


§ ABSTRACT 
We study real-time operator evolution using sparse Pauli dynamics, a recently developed method for simulating expectation values of quantum circuits. On the examples of energy and charge diffusion in 1D spin chains and sudden quench dynamics in the 2D transverse-field Ising model, it is shown that this approach can compete with state-of-the-art tensor network methods. We further demonstrate the flexibility of the approach by studying quench dynamics in the 3D transverse-field Ising model which is highly challenging for tensor network methods.
For the simulation of expectation value dynamics starting in a computational basis state, we introduce an extension of sparse Pauli dynamics that truncates the growing sum of Pauli operators by discarding terms with a large number of X and Y matrices. This is validated by our 2D and 3D simulations. Finally, we argue that sparse Pauli dynamics is not only capable of converging challenging observables to high accuracy, but can also serve as a reliable approximate approach even when given only limited computational resources.


Real-time operator evolution in two and three dimensions via sparse Pauli dynamics
    Garnet Kin-Lic Chan
    September 9, 2024
==================================================================================



./figures/



§ INTRODUCTION


Numerical simulations of quantum dynamics are essential for our understanding of strongly correlated many-body physics. For large systems and long-time dynamics, exact state-vector or tensor network contraction <cit.> methods become computationally intractable and approximate numerical methods are required. These can be formulated in different pictures, for example, the Schrödinger picture, where the state is evolved, or the Heisenberg picture, where the state is unchanged and time evolution is applied to the observable.

Operator time evolution appears naturally in the dynamics of high-temperature systems <cit.> and the theory of operator spreading <cit.>, and has proven useful in the computation of time-correlation functions <cit.> and out-of-time-order correlations <cit.>. Compared to the Schrödinger picture, working in the Heisenberg picture has the benefit that one can take advantage of the dynamical light-cone structure, i.e., the fact that the support of some relevant observables is initially local and grows in time. Operator time evolution has been the subject of tensor network studies, e.g., those based on matrix-product operators (MPO) <cit.> or projected entangled-pair operators (PEPO) <cit.>, whose performance depends on the degree of operator entanglement in the evolved observable.

Alternative Heisenberg-picture methods have been formulated to take advantage of the sparsity of the observable in the Pauli operator representation. 

Under time evolution, local observables spread to a superposition of an increasing number of Pauli operators, and 
various strategies to curb the (at worst) exponential growth of the number of Pauli operators have been employed, including stochastic sampling of Pauli paths <cit.> or discarding Pauli paths based on Pauli weights <cit.>, perturbation order with respect to the nearest Clifford transformation <cit.>, Fourier order <cit.>, or combinations thereof <cit.>. These methods have mainly been developed in the context of simulating quantum circuit expectation values, where the number of Pauli operators can be low due to the use of Clifford gates or noise <cit.>.

Here, we consider one of these methods, namely the sparse Pauli dynamics (SPD) <cit.>, which was applied successfully recently to simulate the kicked transverse-field Ising model quantum simulation of Ref. <cit.>. There, we showed that it can simulate the expectation values faster than the quantum device to reach comparable accuracy and that, given more time, can produce more accurate results than the quantum processor.  These simulations also demonstrated that SPD performs well compared to other classical simulation strategies based on state <cit.> and operator evolution <cit.>. In this work, we analyze its performance in simulating real-time dynamics on 1D spin chains, where we compare it to MPO dynamics including its recent dissipation assisted variant <cit.>, as well as in the transverse-field Ising model on square (2D) and cubic (3D) lattices, where we use available 2D tensor network simulation benchmarks. For the latter, we introduce a modification of the original method that further reduces the computational cost of simulating time-dependent expectation values when the initial state is a computational basis state (e.g., | 0^⊗ n⟩).



§ REAL-TIME SPARSE PAULI DYNAMICS


In SPD we write the observable operator

    O = ∑_P ∈𝒫 a_P P

as a sum of Pauli operators P with complex coefficients a_P, where 𝒫 is a set of Pauli operators that contribute to O. The central part of our algorithm is the action of a Pauli rotation operator U_σ(θ) = exp(-i θσ / 2) on the observable (<ref>):

    Õ = U_σ(θ)^† O U_σ(θ) = ∑_P ∈𝒫_C a_P P + ∑_P ∈𝒫_A (a_P cos(θ) P + i a_P sin(θ) σ· P),

which follows from

    U_σ(θ)^† P U_σ(θ) = 
            P,                        [P, σ] = 0, 
    cos(θ) P + i sin(θ) σ· P,   {P, σ} = 0,

where [·, ·] denotes a commutator and {·, ·} an anticommutator. In Eq. (<ref>), 𝒫_C (𝒫_A) is a set of Pauli operators in 𝒫 that commute (anticommute) with σ. The right-hand side of Eq. (<ref>) can be brought into the form of Eq. (<ref>) by identifying which σ· P already exist in 𝒫 and which have to be added to represent the rotated observable Õ. In general, the number of Pauli operators that we have to store, N = |𝒫|, grows exponentially with the number of unitary Pauli rotation operators applied. To limit the growth of N, we replace the exactly rotated observable Õ by

    Õ_δ = Π_δ (Õ),

where Π_δ acts by discarding all Pauli operators P with |a_P| < δ. Here, the threshold δ defines the approximation error, i.e., δ = 0 corresponds to exact dynamics. In practice, the goal is to converge the simulation with respect to this tunable parameter. Additional implementation details can be found in Appendix <ref>.

For real-time dynamics under a Hamiltonian

    H = ∑_j c_j H_j,

where c_j are real coefficients and H_j are Pauli operators, we replace the exact time-evolution operator U_Δ t = exp(- i H Δ t), corresponding to a time step Δ t, by the first-order Trotter splitting formula

    U_Δ t≈∏_j exp(- i c_j Δ t H_j).

Now the real-time evolution takes the form of applying multiple Pauli rotation gates, which allows us to use the SPD method as defined above. The size of the time step determines not only the Trotter error in Eq. (<ref>) but also the truncation error in SPD. Specifically, for small values of threshold δ and time step Δ t, the error is a function of their ratio δ / Δ t. To show this, let us consider one substep of the dynamics in which we apply one Pauli rotation operator exp(- i c_j Δ t H_j) to the observable (<ref>), the result of which is shown in Eq. (<ref>) with θ = 2 c_j Δ t. For a sufficiently small time step, we can expand Eq. (<ref>) up to first order in Δ t,

    Õ≈∑_P ∈𝒫_C a_P P + ∑_P ∈𝒫_A (a_P  P + 2 i a_P c_j Δ t σ· P) = O + 2 i c_j Δ t ∑_P ∈𝒫_A a_P H_j · P.

For small threshold δ, we can assume that the threshold-based truncation will mainly discard Pauli operators H_j · P that do not already exist in O. The condition under which they are discarded reads |2 c_j Δ t a_P| < δ, i.e., |a_P| < (δ / Δ t) / 2|c_j|. Consequently, the truncation error and the number of Pauli operators depend only on the ratio δ / Δ t. Since the number of Pauli operators determines the computational cost per time step, it is generally preferred to use a larger time step—at fixed δ / Δ t, the cost per time step is constant, but a larger time step requires less steps to reach a given total time. In turn, the size of the time step is limited by the Trotter error, which must be validated before converging the calculation with respect to truncation threshold δ.



§ RESULTS




 §.§ Spin and energy diffusion constants in 1D chains


We begin with the computation of diffusion constants of conserved densities in spin chains of length L at infinite temperature. The diffusion constant <cit.>

    D = 1/2∂/∂ t d^2(t)

is defined through the time-derivative of the mean-square displacement

    d^2(t) = ∑_j C_j(t) j^2 - ( ∑_j C_j(t) j )^2,

where

    C_j(t) = Tr[q_j q_L+1/2(t)]/{∑_j Tr[q_j q_L+1/2(0)] }

are the dynamical correlations between the operators q localized at the central site (L+1)/2 and sites j. These operators represent conserved densities in the sense that ∑_j q_j(t) = ∑_j q_j(0), which leads to ∑_j C_j(t) = ∑_j C_j(0) = 1.

This problem is naturally formulated in the Heisenberg picture and we employ SPD to numerically evolve q_L+1/2(t) in time. Within the sparse Pauli representation of the operator, it is also easy to evaluate the overlap with another Pauli operator (or a sum of Pauli operators) q_j (see Appendix <ref>). Since the conservation laws are not strictly obeyed by the non-unitary truncation scheme employed in SPD, we replace the denominator in Eq. (<ref>) by ∑_j Tr[q_j q_L+1/2(t)], i.e., we renormalize the correlations at post-processing.



  §.§.§ Models


Below, we introduce the examples of Ref. <cit.>, which were studied using dissipation assisted operator evolution (DAOE) combined with a MPO representation of the operator. There, the authors introduced an artificial dissipator that reduces the entanglement of the said MPO and demonstrated that this computational strategy is well founded for the simulation of diffusive transport at high temperature.

The first example is the one-dimensional tilted-field Ising model <cit.>

    H = ∑_j=1^L-1 Z_j Z_j+1 + ∑_j=1^L (1.4 X_j + 0.9045 Z_j)

 with open boundary conditions, while the conserved densities are the local energies
 
    q_j = 1/2 Z_j-1 Z_j + 1/2 Z_j Z_j+1 + 1.4 X_j + 0.9054 Z_j.

L=51 is the number of sites in the chain.

The second model is the XX-ladder Hamiltonian <cit.>

    H = 1/4∑_j=1^L-1∑_a=1,2 (X_j,a X_j+1, a + Y_j,a Y_j+1,a) + 1/4∑_j=1^L (X_j,1 X_j, 2 + Y_j,1 Y_j,2),

where the total spin is a conserved property and we consider the diffusion of q_j = (Z_j,1 + Z_j,2)/2 along the chain of length L=41 (number of sites is n=2L=82).

For these two models, we simulated the mean-square displacement (<ref>) for a total time of t = 20, with a time step of Δ t = 0.02, unless stated otherwise. The diffusion constant was computed by linear regression of d^2(t) between t=10 and t=20.



  §.§.§ Numerical results


Figure <ref> shows the mean-square displacements for the tilted-field Ising (panel A) and XX-ladder models (panel C), simulated using SPD with different values of the threshold δ ranging from 2^-18 to 2^-13. We can observe how d^2(t) converges onto a straight line at large t as we reduce the threshold and make the simulations more accurate. Figures <ref>B, D plot the corresponding diffusion constant as a function δ/Δ t for two sets of data, one using Δ t = 0.01 and the other using Δ t = 0.02. For this property, the results are amenable to an extrapolation in δ→ 0, and we find D ≈ 1.4 for the tilted-field Ising chain (same value was reported in Ref. <cit.>), while D ≈ 0.94 for the XX-ladder (compared to D ≈ 0.95 of Ref. <cit.> and D ≈ 0.96-0.98 reported in Ref. <cit.>). The plots of the diffusion constant also reveal the scaling relationship discussed in Sec. <ref>. Namely, as the threshold is reduced, the simulations using different δ and Δ t but the same δ / Δ t become closer to each other. In addition, we can numerically verify that using a larger time step reduces the computational cost. For example, the two points in Fig. <ref>D at fixed δ / Δ t = 2^-18/0.02 ≈ 0.00019 take around 84h (Δ t = 0.01, blue circle) and 43h (Δ t = 0.02, orange square) to simulate on 6 cores.



To illustrate that the above problems are non-trivial to simulate, we present a comparison between SPD and matrix-product operator (MPO) dynamics without dissipation assistance (Fig. <ref>, see details in Appendix <ref>), a standard method for 1D dynamics.
Using a MPO bond dimension up to χ=2^9 we find that the MPO simulations are still far from the exact results even for shorter chains (L=9-21) of the tilted-field Ising model, whereas SPD is more accurate already at rather large values of δ. Similarly, Ref. <cit.> reported that the MPO dynamics of the L=51 chain with bond dimension χ = 2^9 diverges from the exact result already at t=8. For the same system, SPD is visually well converged up to t=15 already at δ = 2^-15.




In the remainder of this Section, we analyze the operator evolution in terms of contributions from Pauli operators of different weights. Here, the Pauli weight is defined as the number of non-identity Pauli matrices in a Pauli operator. To this end, for any operator O=∑_P ∈𝒫 a_P P, we introduce F_m = ∑_P ∈𝒫_m |a_P|^2, where 𝒫_m ⊆𝒫 is a subset of Pauli operators with Pauli weight m. Then the sum F = ∑_m=1^n F_m is constant for unitary dynamics and equal to the square of the Frobenius norm, i.e., F = Tr(O^†O)/2^n. Figure <ref>A shows the breakdown of F into individual components F_m (up to m=12) for the example of tilted-field Ising model [Eqs. (<ref>) and (<ref>)], demonstrating how the dynamics evolves initially low-weight Pauli operators into a sum of operators with higher weights. After some time, F of the operator evolved with SPD deviates from its initial value because of threshold-based truncation. This truncation appears to affect high-weight Pauli operators more than the low-weight Paulis. This is confirmed in Fig. <ref>B, where we show F_m contributions for m up to 5. As shown in previous works <cit.>, Pauli operators corresponding to the local conserved quantity (here m=1 and m=2 for the local energy (<ref>)) obey the long-time scaling F_m ∼ t^-1/2, whereas other Pauli operators (here m>2) obey F_m ∼ t^-3/2 in the long-time limit.



This observation motivates methods that truncate the time-evolved operator based on Pauli weights, including the dissipation assisted operator evolution (DAOE) approach <cit.>, low-weight simulation algorithm (LOWESA) <cit.> and observable's back-propagation on Pauli paths (OBPPP) <cit.>, two approaches that were originally developed for the simulation of noisy quantum circuits, and the restricted state space approximation <cit.>, introduced in the context of simulating nuclear magnetic resonance experiments. Our results indicate that SPD can take advantage of such behavior without explicitly truncating the sum based on Pauli weights because high-weight Pauli operators appear with small coefficients that do not meet the threshold criterion.



 §.§ Time-dependent expectation values




  §.§.§ X-truncated sparse Pauli dynamics


In the examples above, we focused on infinite-temperature time-correlation functions of a time-evolved operator with a local, low-weight Pauli operator. Here, we consider expectation values of the form ⟨ O ⟩_t = ⟨ 0^⊗ n | O(t) | 0^⊗ n⟩, where high-weight Pauli operators can contribute as long as they are composed only of identity and Z matrices (i.e., if they are diagonal in the computational basis). Default Heisenberg evolution does not account for the information about the state over which we take the expectation value but rather treats all Pauli operators equally. Within the framework of SPD, this means that we keep a large number of Pauli operators, of which only a fraction contribute to the observable.

To further truncate the number of Pauli operators without introducing large errors in the expectation value ⟨ O ⟩_t, we propose to discard Pauli operators composed of more than M X or Y Pauli matrices. We refer to the number of X/Y matrices as the X-weight of a Pauli operator and we call this additional truncation scheme X-truncated SPD (xSPD). The truncation introduces the additional assumption that for certain (short) times, there is limited operator backflow from high X-weight Paulis to the manifold of Z-type Pauli operators. For each calculation, we test the value of M to ensure that the error introduced by the X-truncation scheme is sufficiently small (for example, smaller than the target convergence criterion). The X-truncation is applied only every T steps of the dynamics to limit the impact on the accuracy. In our calculations, we fixed T to 5 steps. Finally, we note here that alternatives to our hard M cutoff could also be considered, such as introducing an artificial dissipator based on the Pauli's X-weight that would be similar in spirit to DAOE <cit.>.

In the following, we apply this modification of the original SPD method to dynamics in the 2D (square lattice) and 3D (cubic lattice) transverse-field Ising models described by the Hamiltonian

    H = -∑_⟨ j k ⟩ X_j X_k - h ∑_j=1^L  Z_j,

where the first sum runs over nearest neighbors on the lattice with open boundary conditions and h controls the magnitude of the field. We consider the time-dependent magnetization ⟨ Z⟩_t = ⟨ 0^⊗ n| Z_0(t) |0^⊗ n⟩, where Z_0 denotes the Z Pauli operator on the central site. Physically this corresponds to the magnetization induced after a sudden quench from infinite h to a finite value of h.

In this setting, the first-order Trotter splitting

    U^(1)(t) = [e^ i Δ t ∑_⟨ j k ⟩ X_j X_k e^ i h Δ t ∑_j=1^L  Z_j]^K,

where K = t/Δ t is the number of time steps, is equivalent to the second-order splitting

    U^(2)(t)    =   [e^ i 1/2 h Δ t ∑_j=1^L  Z_j e^ i Δ t ∑_⟨ j k ⟩ X_j X_k e^ i 1/2 h Δ t ∑_j=1^L  Z_j]^K 
       =    e^ i 1/2 h Δ t ∑_j=1^L  Z_j U^(1)(t) e^ - i 1/2 h Δ t ∑_j=1^L  Z_j

because Z-Pauli rotations commute with the observable and apply only a phase to the initial state.



  §.§.§ 2D transverse-field Ising model


The quantum quench dynamics of magnetization in the 2D transverse-field Ising model has been studied by means of infinite projected entangled pair state (iPEPS) <cit.> and neural network quantum state <cit.> calculations. While the iPEPS simulations correspond to dynamics in the thermodynamic limit, neural network simulations were performed on a 10 × 10 lattice, which was shown to be sufficiently large to exhibit negligible finite-size effects <cit.>. In our simulations, we used an 11 × 11 square lattice. We set h=h_c, where h_c = 3.04438(2) corresponds to the quantum critical point <cit.>, and simulate the magnetization up to t = 0.92, where we can compare our results to different update schemes used in iPEPS simulations, namely the full update (FU) <cit.>, neighborhood tensor update (NTU) <cit.>, and gradient tensor update (GTU) <cit.>.

Figure <ref> shows the convergence of xSPD with respect to the threshold δ. As expected, the method converges quickly at short times but requires small values of δ to converge the values at longer times. Our most accurate simulation agrees well with FU and GTU iPEPS results, and shows some deviation from the NTU scheme at the end of the simulation time. 
We note that although the NTU iPEPS calculation corresponds to the largest bond dimension amongst the iPEPS data, the accuracy of truncation is also believed to be less than for the FU iPEPS simulation, thus the relative accuracy of the different reference iPEPS data is unclear.
The disagreement between the two smallest δ xSPD simulations at the end of the simulation time is only 0.002 (0.007 over the 3 smallest δs) which provides an estimate of the threshold error. This threshold error is comparable to the estimated Trotter and X-truncation errors (discussed below).



In this example, we employed a time step of Δ t = 0.04 and set the X-truncation parameter to M=5. To validate this choice of parameters, we analyze the associated errors in Fig. <ref>. Specifically, for the time step (Fig. <ref>A), we set δ / Δ t = 2^-19/0.01  = 1.90734 × 10^-4 and compute the observable using three different time steps. We estimate that the Trotter error is below 0.003 within the simulated time of t=0.92. Similarly, we perform SPD and xSPD simulations with varying values of M, using δ=2^-20. We estimate that the error due to M=5 X-truncation is about 0.003. In contrast, employing M=7 would lead to almost no error but also limited computational savings, while M=3 produces an error that is greater than our convergence target of less than 0.01. Note that due to symmetry, all Pauli operators appear with an even X-weight, which is why we only consider odd values of M.



Regarding the computational cost, the most accurate simulation at δ = 2^-23 generated up to 8.5 billion Pauli operators, used over 1TB of memory, and took around 36 hours on 16 CPU cores. For comparison, the least accurate simulation shown in Fig. <ref>, corresponding to δ = 2^-18, generated at most 84 million Pauli operators in 32 minutes on 16 CPU cores.



  §.§.§ 3D transverse-field Ising model


As our final example, we present the quench dynamics of magnetization for spins on a simple cubic lattice with L=11 and n=L^3=1331 sites. Here we consider two values of h, h=1 (weak field) and h=5.15813(6) (critical point). 

In the case of h=1, using Δ t = 0.04 and M=7, we could run the dynamics up to t=1 with thresholds as low as 2^-19 (see Fig. <ref>). The most accurate simulations (δ≤ 2^-17) agree to within ≈ 0.02, which is comparable to the estimated time step (Trotter) and X-truncation errors (see Fig. <ref>A, C in Appendix <ref>). Interestingly, even the fastest, least-accurate simulation (δ=2^-14) exhibits “qualitative accuracy”, i.e., recovers the general trend of the most accurate available result. We ascribe this to the fact that sparse Pauli representation can easily reproduce the dynamics dominated by few Pauli operators. However, the method struggles to include small contributions from many Pauli operators generated by the time evolution. For example, while only about 10^6 Pauli operators are generated at threshold δ=2^-14 after 1 time unit, around 7 × 10^8 Pauli operators are generated during the same dynamics with a threshold of δ=2^-19. Yet, the difference in the observable appears limited to around 0.04.

The system with h=5.158136 (Fig. <ref>B) was propagated up to t=0.6 using Δ t= 0.02 and M=5. For this choice of parameters, the Trotter and X-truncation errors are estimated to be below 0.005 (see Fig. <ref>B, D in Appendix <ref>). Because of the shorter dynamics we could use a smaller value of the X-truncation parameter M compared to the weak-field case. The results, using thresholds as low as δ=2^-20, are converged to below 0.01 for times t<0.5, after which our most accurate simulations begin to deviate from each other.

In general, these 3D calculations are expected to pose challenges for tensor network techniques, which for a fixed bond dimension show an exponential scaling with the connectivity of the lattice (assuming site tensors with a number of bonds equal to the number of neighbours). Although SPD does not show this exponential scaling, the 3D simulations for a given side-length L are still more costly than the 2D ones, primarily because the number of sites n=L^3 is a factor of L larger than in the 2D case.
For these reasons, we cannot simulate as many Pauli operators as in the 2D case. Specifically, our memory budget of about 1.5TB is reached already with less than 10^9 Pauli operators, an order of magnitude less than in our 2D square lattice simulations. With this number of Pauli operators, the computation with h=1 takes around 73 h on 16 CPU cores, about two times longer than our most accurate 2D calculation. Nonetheless, the relative feasibility of these simulations illustrate the kinds of systems that can be (approximately) studied by SPD dynamics that would otherwise be challenging to consider.





  §.§.§ Computational savings from X-truncation


We now analyze the savings due to the X-truncation scheme. Figures <ref>A, B show the number of Pauli operators as a function of time for 2D and 3D simulations. While the number of Pauli operators in xSPD is suppressed by the X-truncation, the number of Z-type operators (Pauli operators composed only of Z and identity matrices) is almost the same as in SPD. Since the total cost of the computation is proportional to the number of all Pauli operators N, we can use the ratio N_ SPD/N_ xSPD to quantify the computational savings (Fig. <ref>C). We observe a factor of 6 decrease in the number of Pauli operators in the 2D case, and up to a factor of 5 for the 3D simulation. Such savings allow a factor of 2–4 lower thresholds compared to the original SPD. For example, within a budget of a few days and 1.5 TB of memory, the most accurate SPD calculation we could run for the 2D system would be limited to δ=2^-21, which is not converged with respect to threshold at the longest simulated times (see Fig. <ref>). In contrast, with xSPD we could afford to run the same simulation with δ=2^-23.





§ CONCLUSION


To conclude, we have presented numerical simulations of real-time operator evolution with SPD and its modified version, xSPD. We have shown that, in systems for which reference data is available, the performance of these methods is at the level of state-of-the-art tensor network approaches, while the flexibility and simplicity of SPD allows for applications to dynamics problems where tensor network approaches have yet to make an impact, such as in 3D lattices. For the studied systems, we found that we can obtain very accurate results either by converging the observables in the case of short-time, 2D and 3D transverse field Ising model dynamics, or by extrapolating to zero threshold for the long-time diffusion coefficients in 1D chains. Apart from reaching quantitatively converged results for these challenging examples of time-dependent observables, we also found that SPD could serve as a practical method for computing qualitatively accurate results. In this respect, our work motivates further research into using sparse Pauli representations for real-time quantum dynamics.


We thank Jacek Dziarmaga for sharing the iPEPS simulation data presented in Fig. <ref> and Huanchen Zhai for help with setting up MPO simulations presented in Fig. <ref>.
The authors were supported by the US Department of Energy, Office of Science, Office of Advanced Scientific Computing Research and Office of Basic Energy Sciences, Scientific Discovery through Advanced Computing (SciDAC) program under Award Number DE-SC0022088. TB acknowledges financial support from the Swiss National Science Foundation through the Postdoc Mobility Fellowship (grant number P500PN-214214). GKC is a Simons Investigator in Physics. Computations presented here were conducted in the Resnick High Performance Computing Center, a facility supported by Resnick Sustainability Institute at the California Institute of Technology.






§ ADDITIONAL DETAILS OF SPD IMPLEMENTATION AND WORKING EQUATIONS


SPD is implemented as described in Ref. <cit.>. Briefly, a sum of N Pauli operators is stored in the form of three arrays: an array of N complex coefficients a, an array of N integer phases φ, and a N × 2n_64 array of 64-bit unsigned integers that stores two bitstrings x and z for each Pauli operator:

    O = ∑_j=0^N-1 a_j (-i)^φ∏_k=0^n-1 Z_k^z_jk X_k^x_jk.

The number of unsigned integers needed to store n bits is n_64 = ⌈ n ⌉. Pauli operators are sorted using lexicographic ordering on the bitstrings. In this way, we can find the position j of a given Pauli operator in the sum—or the position at which to insert a new Pauli operator so that the ordering is preserved—in 𝒪(log N) time. Similarly, deleting Pauli operators preserves the order trivially.

Apart from searching, inserting, and deleting Pauli operators, other key operations on this sparse representation of a sum of Pauli operators include identifying which Pauli operators in the sum anticommute with a given Pauli operator and multiplying the sum of Pauli operators by a Pauli. For the anticommutation of Pauli operators A and B, we have

    A  anticommutes with  B = z_A · x_B - x_A · z_B.
 
Here the multiplications on the right-hand side correspond to the AND logical operator between bits and additions to the XOR logical operator (addition in ℤ_2). The product of two Pauli operators C = A B is given by

    (z_C, x_C)    =    (z_A + z_B, x_A + x_B), 
    φ_C    =   φ_A + φ_B + 2 z_A · x_B, 
    
        a_C    =    a_A a_B.


In Sec. <ref>, we also introduced an overlap (inner product) between sums of Pauli operators Tr[O_1^† O_2]/2^n represented in the sparse Pauli format described above. Let us assume that N_1 < N_2, i.e., that O_1 has fewer Pauli operators than O_2. Then, we can search for all Pauli operators of O_1 in O_2 in 𝒪(N_1 log N_2) time and the overlap is 

    1/2^nTr[O_1^† O_2] = ∑_j (-i)^φ_2, k[j] - φ_1, j a_1, j^∗ a_2, k[j],

where the sum runs over Pauli operators in O_1 that were found in O_2 and k[j] is the index of the j-th found Pauli in O_2. The expectation over the all-zero state, as needed in Sec. <ref>, can be computed as

    ⟨ 0^⊗ n | O | 0^⊗ n⟩ = ∑_j a_j (-i)^φ_j,

where the sum runs over Z-type Pauli operators (for which all bits in x_j are 0). Finally, the X-weight of an operator is computed as the number of set bits in the corresponding x array.

For convenience, our implementation interfaces to Qiskit <cit.> for setting up the calculations. Specifically, it converts Qiskit's SparsePauliOp into our representation described above that is then used in the simulations.



§ DETAILS OF MPO SIMULATIONS


MPO simulations were performed using the time-dependent density matrix renormalization group (TD-DMRG) method, as implemented in the Block2 code <cit.>. After constructing the MPO of the observable at time zero, we convert it to an MPS | O with a doubled number of sites, i.e., ⟨ i_1 i_2 … | O | j_1 j_2 …⟩ / 2^n/2 =  i_1 j_1 i_2 j_2 … | O. The sites of the MPS are ordered so that the two physical legs of a single site in the MPO appear on neighboring sites in the MPS. The Liouvillian superoperator L |O≡ [H, O] that governs the dynamics of the observable is then represented as an MPO in the extended Hilbert space with twice as many sites. Specifically, each Pauli operator in the Hamiltonian corresponds to a sum of two Pauli operators in the Liouvillian. For single-site Pauli operators σ_i ∈{I, X, Y, Z} at site i, we have σ_2i - σ_2i+1, while the nearest-neighbor interaction terms σ_i σ_i+1 correspond to σ_2iσ_2i+2 - σ_2i+1σ_2i+3 in the superoperator picture.

The initial observable MPS is propagated with the Liouvillian MPO using the time-step-targeting method <cit.> and a time step of Δ t = 0.04. The correlation functions of the form Tr[O_1^† O_2]/2^n, used in Eq. (<ref>), were evaluated as the inner product O_1 | O_2 of the two matrix product states.



§ TIME STEP AND X-TRUNCATION PARAMETERS IN 3D SIMULATIONS














46
fxundefined [1]
 ifx#1

fnum [1]
 #1firstoftwo
 secondoftwo
 
fx [1]
 #1firstoftwo
 secondoftwo
 





noop [0]secondoftwo

ref[1]@startlink#1@href
href[1]#1@endlink
anitize@url [0]`
12`$12`&12`#12`1̂2`_12`%12
startlink[1]
endlink[0]

rl [1]href #1












@bib@innerbibempty

[Gray and Kourtis(2021)]gray2021hyper
  author author J. Gray and author S. Kourtis,  title title Hyper-optimized tensor network
  contraction, @noop  journal journal
  Quantum volume 5, pages 410
  (year 2021)NoStop
[Huang et al.(2021)Huang,
  Zhang, Newman, Ni,
  Ding, Cai, Gao, Wang, Wu, Zhang et al.]huang2021efficient
  author author C. Huang, author F. Zhang,
  author M. Newman, author X. Ni, author
  D. Ding, author J. Cai, author X. Gao, author T. Wang, author F. Wu, author
  G. Zhang, et al., title title Efficient parallelization of tensor network
  contraction for simulating quantum computation, @noop  journal journal Nature Computational Science volume 1, pages 578 (year
  2021)NoStop
[Gray and Chan(2024)]gray2022hyper
  author author J. Gray and author G. K.-L. Chan, title title Hyperoptimized Approximate
  Contraction of Tensor Networks with Arbitrary Geometry, https://doi.org/10.1103/PhysRevX.14.011009 journal journal Phys. Rev. X volume 14, pages 011009 (year 2024)NoStop
[Rakovszky et al.(2022)Rakovszky, von Keyserlingk, and Pollmann]rakovszky2022dissipation
  author author T. Rakovszky, author C. W. von
  Keyserlingk, and author F. Pollmann, title title Dissipation-assisted
  operator evolution method for capturing hydrodynamic transport, https://doi.org/10.1103/PhysRevB.105.075131 journal journal Phys. Rev. B volume 105, pages 075131 (year 2022)NoStop
[Nahum et al.(2018)Nahum,
  Vijay, and Haah]Nahum2018
  author author A. Nahum, author S. Vijay, and  author J. Haah, title title Operator Spreading in Random Unitary Circuits,  https://doi.org/10.1103/PhysRevX.8.021014 journal
  journal Phys. Rev. X volume 8,  pages 021014 (year 2018)NoStop
[von Keyserlingk et al.(2018)von Keyserlingk, Rakovszky, Pollmann, and  Sondhi]VonKeyserlingk2018
  author author C. W. von Keyserlingk, author T. Rakovszky, author F. Pollmann, and author S. L. Sondhi, title title Operator Hydrodynamics,
  OTOCs, and Entanglement Growth in Systems without Conservation Laws,  https://doi.org/10.1103/PhysRevX.8.021013 journal
  journal Phys. Rev. X volume 8,  pages 021013 (year 2018)NoStop
[Khemani et al.(2018)Khemani, Vishwanath, and Huse]Khemani2018
  author author V. Khemani, author A. Vishwanath, and author D. A. Huse, title title Operator Spreading and the
  Emergence of Dissipative Hydrodynamics under Unitary Evolution with
  Conservation Laws, https://doi.org/10.1103/PhysRevX.8.031057
  journal journal Phys. Rev. X volume 8, pages 031057 (year
  2018)NoStop
[von Keyserlingk et al.(2022)von Keyserlingk, Pollmann, and Rakovszky]VonKeyserlingk2022
  author author C. von
  Keyserlingk, author F. Pollmann, and author T. Rakovszky, title title Operator backflow and
  the classical simulation of quantum transport, https://doi.org/10.1103/PhysRevB.105.245101 journal journal Phys. Rev. B volume 105, pages 245101 (year 2022)NoStop
[Mi et al.(2021)Mi,
  Roushan, Quintana, Mandrà, Marshall, Neill,
  Arute, Arya, Atalaya,
  Babbush, Bardin, Barends,
  Basso, Bengtsson, Boixo,
  Bourassa, Broughton, Buckley,
  Buell, Burkett, Bushnell,
  Chen, Chiaro, Collins,
  Courtney, Demura, Derk,
  Dunsworth, Eppens, Erickson,
  Farhi, Fowler, Foxen,
  Gidney, Giustina, Gross,
  Harrigan, Harrington, Hilton,
  Ho, Hong, Huang,
  Huggins, Ioffe, Isakov,
  Jeffrey, Jiang, Jones,
  Kafri, Kelly, Kim,
  Kitaev, Klimov, Korotkov,
  Kostritsa, Landhuis, Laptev,
  Lucero, Martin, McClean,
  McCourt, McEwen, Megrant,
  Miao, Mohseni, Montazeri,
  Mruczkiewicz, Mutus, Naaman,
  Neeley, Newman, Niu,
  O'Brien, Opremcak, Ostby,
  Pato, Petukhov, Redd,
  Rubin, Sank, Satzinger,
  Shvarts, Strain, Szalay,
  Trevithick, Villalonga, White, Yao, Yeh, Zalcman,
  Neven, Aleiner, Kechedzhi,
  Smelyanskiy, and Chen]Mi2021
  author author X. Mi, author P. Roushan,
  author C. Quintana, author S. Mandrà, author
  J. Marshall, author
  C. Neill, author F. Arute, author K. Arya, author J. Atalaya, author R. Babbush,
  author J. C. Bardin, author R. Barends, author
  J. Basso, author A. Bengtsson, author S. Boixo, author A. Bourassa, author M. Broughton, author B. B.  Buckley, author D. A.  Buell, author B. Burkett, author N. Bushnell,
  author Z. Chen, author
  B. Chiaro, author R. Collins, author W. Courtney, author S. Demura, author A. R.  Derk, author A. Dunsworth, author D. Eppens,
  author C. Erickson, author E. Farhi, author
  A. G. Fowler, author
  B. Foxen, author C. Gidney, author M. Giustina, author J. A.  Gross, author M. P. Harrigan, author S. D. Harrington, author J. Hilton,
  author A. Ho, author
  S. Hong, author T. Huang, author W. J.  Huggins, author L. B.  Ioffe, author S. V. Isakov, author E. Jeffrey,
  author Z. Jiang, author C. Jones, author
  D. Kafri, author J. Kelly, author S. Kim, author A. Kitaev, author P. V. Klimov,
  author A. N. Korotkov, author F. Kostritsa, author
  D. Landhuis, author
  P. Laptev, author E. Lucero, author O. Martin, author J. R.  McClean, author T. McCourt, author M. McEwen,
  author A. Megrant, author K. C. Miao, author
  M. Mohseni, author S. Montazeri, author W. Mruczkiewicz, author J. Mutus, author O. Naaman, author M. Neeley, author M. Newman, author M. Y.  Niu, author T. E. O'Brien, author A. Opremcak,
  author E. Ostby, author B. Pato, author
  A. Petukhov, author
  N. Redd, author N. C. Rubin, author D. Sank, author K. J.  Satzinger, author V. Shvarts, author D. Strain,
  author M. Szalay, author M. D. Trevithick, author B. Villalonga, author
  T. White, author Z. J. Yao, author P. Yeh, author A. Zalcman, author H. Neven,
  author I. Aleiner, author K. Kechedzhi, author
  V. Smelyanskiy, and author
  Y. Chen, title title Information scrambling in quantum circuits, https://doi.org/10.1126/science.abg5029 journal journal Science volume 374, pages
  1479 (year 2021)NoStop
[Verstraete et al.(2004)Verstraete, García-Ripoll, and Cirac]Verstraete2004
  author author F. Verstraete, author J. J. García-Ripoll, and author J. I.  Cirac, title title Matrix
  Product Density Operators: Simulation of Finite-Temperature and Dissipative
  Systems, https://doi.org/10.1103/PhysRevLett.93.207204 journal journal Phys. Rev. Lett. volume 93, pages 207204 (year
  2004)NoStop
[Zwolak and Vidal(2004)]Zwolak2004
  author author M. Zwolak and author G. Vidal,  title title Mixed-State Dynamics in
  One-Dimensional Quantum Lattice Systems: A Time-Dependent Superoperator
  Renormalization Algorithm, https://doi.org/10.1103/PhysRevLett.93.207205 journal
  journal Phys. Rev. Lett. volume 93,  pages 207205 (year 2004)NoStop
[Anand et al.(2023)Anand,
  Temme, Kandala, and Zaletel]Anand2023
  author author S. Anand, author K. Temme,
  author A. Kandala, and author M. Zaletel, http://arxiv.org/abs/2306.17839 title Classical benchmarking
  of zero noise extrapolation beyond the exactly-verifiable regime (year 2023), https://arxiv.org/abs/2306.17839
  arXiv:2306.17839 NoStop
[Alhambra and Cirac(2021)]Alhambra2021
  author author Á. M.  Alhambra and author J. I.  Cirac, title title Locally
  Accurate Tensor Networks for Thermal States and Time Evolution, https://doi.org/10.1103/PRXQuantum.2.040331 journal journal PRX Quantum volume 2, pages
  040331 (year 2021)NoStop
[Liao et al.(2023)Liao,
  Wang, Zhou, Zhang, and  Xiang]liao2023simulation
  author author H.-J. Liao, author K. Wang, author Z.-S. Zhou, author
  P. Zhang, and author
  T. Xiang, @noop  title
  Simulation of IBM's kicked Ising experiment with Projected Entangled Pair
  Operator (year 2023), https://arxiv.org/abs/2308.03082 arXiv:2308.03082 [quant-ph] NoStop
[Rall et al.(2019)Rall,
  Liang, Cook, and Kretschmer]rall2019simulation
  author author P. Rall, author D. Liang,
  author J. Cook, and author W. Kretschmer, title
  title Simulation of qubit quantum circuits via pauli
  propagation, https://doi.org/10.1103/PhysRevA.99.062337 journal journal Phys. Rev. A volume
  99, pages 062337 (year 2019)NoStop
[Shao et al.(2023)Shao,
  Wei, Cheng, and Liu]Shao2023
  author author Y. Shao, author F. Wei, author S. Cheng, and author
  Z. Liu, http://arxiv.org/abs/2306.05804 title Simulating Quantum Mean
  Values in Noisy Variational Quantum Algorithms: A Polynomial-Scale
  Approach (year 2023), https://arxiv.org/abs/2306.05804 arXiv:2306.05804 NoStop
[Begušić et al.(2023)Begušić, Hejazi, and Chan]begusic2023simulating
  author author T. Begušić, author K. Hejazi, and author G. K.-L. Chan, @noop  title Simulating quantum circuit
  expectation values by Clifford perturbation theory (year
  2023), https://arxiv.org/abs/2306.04797 arXiv:2306.04797
  [quant-ph] NoStop
[Nemkov et al.(2023)Nemkov,
  Kiktenko, and Fedorov]Nemkov_Fedorov:2023
  author author N. A. Nemkov, author E. O. Kiktenko, and author A. K. Fedorov, title title Fourier expansion in
  variational quantum algorithms, https://doi.org/10.1103/PhysRevA.108.032406 journal journal Phys. Rev. A volume 108, pages 032406 (year 2023)NoStop
[Fontana et al.(2023)Fontana, Rudolph, Duncan, Rungger, and Cîrstoiu]fontana2023classical
  author author E. Fontana, author M. S. Rudolph, author R. Duncan,
  author I. Rungger, and author C. Cîrstoiu, @noop 
  title Classical simulations of noisy variational quantum
  circuits (year 2023), https://arxiv.org/abs/2306.05400 arXiv:2306.05400 [quant-ph] NoStop
[Rudolph et al.(2023)Rudolph, Fontana, Holmes, and  Cincio]rudolph2023classical
  author author M. S. Rudolph, author E. Fontana,
  author Z. Holmes, and author L. Cincio, @noop 
  title Classical surrogate simulation of quantum systems with
  LOWESA (year 2023), https://arxiv.org/abs/2308.09109 arXiv:2308.09109 [quant-ph] NoStop
[Aharonov et al.(2023)Aharonov, Gao, Landau, Liu, and Vazirani]Aharonov2023
  author author D. Aharonov, author X. Gao,
  author Z. Landau, author Y. Liu, and author
  U. Vazirani, title title A Polynomial-Time Classical Algorithm for Noisy Random Circuit
  Sampling, in https://doi.org/10.1145/3564246.3585234 booktitle Proc. 55th Annu. ACM Symp. Theory Comput.  (publisher ACM, address New York, NY, USA,  year 2023) pp. pages 945–957, https://arxiv.org/abs/2211.03999 arXiv:2211.03999 NoStop
[Schuster et al.(2024)Schuster, Yin, Gao, and Yao]Schuster2024
  author author T. Schuster, author C. Yin,
  author X. Gao, and author N. Y. Yao, https://doi.org/10.48550/arXiv.2407.12768 title A
  polynomial-time classical algorithm for noisy quantum circuits (year 2024), https://arxiv.org/abs/2407.12768
  arXiv:2407.12768 NoStop
[Begušić et al.(2024)Begušić, Gray,  and Chan]Begusic2024
  author author T. Begušić, author J. Gray, and author G. K.-L. Chan, title title Fast and converged
  classical simulations of evidence for the utility of quantum computing before
  fault tolerance, https://doi.org/10.1126/sciadv.adk4321
  journal journal Sci. Adv. volume 10, pages eadk4321 (year
  2024)NoStop
[Kim et al.(2023)Kim,
  Eddins, Anand, Wei,
  Van Den Berg, Rosenblatt, Nayfeh, Wu, Zaletel, Temme, and Kandala]kim2023evidence
  author author Y. Kim, author A. Eddins,
  author S. Anand, author K. X. Wei, author
  E. Van Den Berg, author
  S. Rosenblatt, author
  H. Nayfeh, author Y. Wu, author M. Zaletel, author K. Temme,  and author A. Kandala,  title title Evidence for the utility of quantum
  computing before fault tolerance, https://doi.org/10.1038/s41586-023-06096-3 journal journal Nature volume 618, pages
  500 (year 2023)NoStop
[Kechedzhi et al.(2024)Kechedzhi, Isakov, Mandrà, Villalonga, Mi, Boixo, and Smelyanskiy]Kechedzhi2024
  author author K. Kechedzhi, author S. Isakov,
  author S. Mandrà, author B. Villalonga, author
  X. Mi, author S. Boixo, and author V. Smelyanskiy, title title
  Effective quantum volume, fidelity and computational cost of noisy quantum
  processing experiments, https://doi.org/https://doi.org/10.1016/j.future.2023.12.002 journal journal Future Gener. Comput. Syst. volume 153, pages 431 (year
  2024)NoStop
[Tindall et al.(2024)Tindall, Fishman, Stoudenmire, and  Sels]Tindall2024
  author author J. Tindall, author M. Fishman,
  author E. M. Stoudenmire, and  author D. Sels, title title Efficient Tensor Network Simulation of IBM's
  Eagle Kicked Ising Experiment, https://doi.org/10.1103/PRXQuantum.5.010308 journal journal PRX Quantum volume 5, pages
  010308 (year 2024)NoStop
[Patra et al.(2024)Patra,
  Jahromi, Singh, and Orús]Patra2024
  author author S. Patra, author S. S. Jahromi,
  author S. Singh, and author R. Orús, title
  title Efficient tensor network simulation of IBM's largest
  quantum processors, https://doi.org/10.1103/PhysRevResearch.6.013326 journal
  journal Phys. Rev. Res. volume 6,  pages 013326 (year 2024)NoStop
[Bertini et al.(2021)Bertini, Heidrich-Meisner, Karrasch,
  Prosen, Steinigeweg, and Žnidarič]Bertini2021
  author author B. Bertini, author F. Heidrich-Meisner, author C. Karrasch, author T. Prosen,
  author R. Steinigeweg, and  author M. Žnidarič,  title title Finite-temperature transport in
  one-dimensional quantum lattice models, https://doi.org/10.1103/RevModPhys.93.025003 journal
  journal Rev. Mod. Phys. volume 93,  pages 025003 (year 2021)NoStop
[Kim and Huse(2013)]Kim2013
  author author H. Kim and author D. A. Huse,  title title Ballistic Spreading of Entanglement
  in a Diffusive Nonintegrable System, https://doi.org/10.1103/PhysRevLett.111.127205 journal
  journal Phys. Rev. Lett. volume 111,  pages 127205 (year 2013)NoStop
[Steinigeweg et al.(2014)Steinigeweg, Heidrich-Meisner, Gemmer,
  Michielsen, and De
  Raedt]Steinigeweg2014
  author author R. Steinigeweg, author F. Heidrich-Meisner, author J. Gemmer, author K. Michielsen,  and author H. De Raedt,  title title Scaling of diffusion constants in the
  spin-1/2 XX ladder, https://doi.org/10.1103/PhysRevB.90.094417 journal journal Phys. Rev. B volume 90, pages 094417 (year 2014)NoStop
[Karrasch et al.(2015)Karrasch, Kennes, and Heidrich-Meisner]Karrasch2015
  author author C. Karrasch, author D. M. Kennes, and author F. Heidrich-Meisner, title title Spin and
  thermal conductivity of quantum spin chains and ladders, https://doi.org/10.1103/PhysRevB.91.115130 journal journal Phys. Rev. B volume 91, pages 115130 (year 2015)NoStop
[Kloss et al.(2018)Kloss,
  Lev, and Reichman]Kloss2018
  author author B. Kloss, author Y. B. Lev,  and author D. Reichman,  title title Time-dependent variational principle
  in matrix-product state manifolds: Pitfalls and potential, https://doi.org/10.1103/PhysRevB.97.024307 journal journal Phys. Rev. B volume 97, pages 024307 (year 2018)NoStop
[Rakovszky et al.(2018)Rakovszky, Pollmann, and von
  Keyserlingk]Rakovszky2018
  author author T. Rakovszky, author F. Pollmann, and author C. W. von Keyserlingk, title title Diffusive
  Hydrodynamics of Out-of-Time-Ordered Correlators with Charge Conservation,  https://doi.org/10.1103/PhysRevX.8.031058 journal
  journal Phys. Rev. X volume 8,  pages 031058 (year 2018)NoStop
[Kuprov et al.(2007)Kuprov,
  Wagner-Rundell, and Hore]kuprov2007polynomially
  author author I. Kuprov, author N. Wagner-Rundell, and author P. Hore, title title Polynomially scaling spin
  dynamics simulation algorithm based on adaptive state-space restriction,  https://doi.org/10.1016/j.jmr.2007.09.014 journal
  journal J. Magn. Reson. volume 189,  pages 241 (year 2007)NoStop
[Kuprov(2008)]Kuprov2008
  author author I. Kuprov, title title Polynomially scaling spin
  dynamics II: Further state-space compression using Krylov subspace techniques
  and zero track elimination, https://doi.org/10.1016/j.jmr.2008.08.008 journal journal J. Magn. Reson. volume 195, pages 45 (year 2008)NoStop
[Karabanov et al.(2011)Karabanov, Kuprov, Charnock, van der Drift, Edwards, and Köckenberger]Karabanov2011
  author author A. Karabanov, author I. Kuprov,
  author G. T. P. Charnock,
  author A. van der Drift, author L. J. Edwards, and author W. Köckenberger, title title On the accuracy of the state space restriction
  approximation for spin dynamics simulations, https://doi.org/10.1063/1.3624564 journal journal J. Chem. Phys. volume 135, pages 084106 (year 2011)NoStop
[Hogben et al.(2011)Hogben,
  Krzystyniak, Charnock, Hore,  and Kuprov]Hogben2011
  author author H. Hogben, author M. Krzystyniak,
  author G. Charnock, author P. Hore, and author
  I. Kuprov, title title Spinach – A software library for simulation of spin dynamics in
  large spin systems, https://doi.org/10.1016/j.jmr.2010.11.008
  journal journal J. Magn. Reson. volume 208, pages 179 (year
  2011)NoStop
[Czarnik et al.(2019)Czarnik, Dziarmaga, and Corboz]Czarnik2019
  author author P. Czarnik, author J. Dziarmaga,  and author P. Corboz,  title title Time evolution of an infinite
  projected entangled pair state: An efficient algorithm, https://doi.org/10.1103/PhysRevB.99.035115 journal journal Phys. Rev. B volume 99, pages 035115 (year 2019)NoStop
[Dziarmaga(2021)]Dziarmaga2021
  author author J. Dziarmaga, title title Time evolution of an
  infinite projected entangled pair state: Neighborhood tensor update,  https://doi.org/10.1103/PhysRevB.104.094411 journal
  journal Phys. Rev. B volume 104,  pages 094411 (year 2021)NoStop
[Dziarmaga(2022)]Dziarmaga2022
  author author J. Dziarmaga, title title Time evolution of an
  infinite projected entangled pair state: A gradient tensor update in the
  tangent space, https://doi.org/10.1103/PhysRevB.106.014304
  journal journal Phys. Rev. B volume 106, pages 014304 (year
  2022)NoStop
[Schmitt and Heyl(2020)]Schmitt2020
  author author M. Schmitt and author M. Heyl,  title title Quantum Many-Body Dynamics in Two
  Dimensions with Artificial Neural Networks, https://doi.org/10.1103/PhysRevLett.125.100503 journal
  journal Phys. Rev. Lett. volume 125,  pages 100503 (year 2020)NoStop
[Blöte and Deng(2002)]Blote2002
  author author H. W. J.  Blöte and author Y. Deng, title title Cluster
  Monte Carlo simulation of the transverse Ising model, https://doi.org/10.1103/PhysRevE.66.066110 journal journal Phys. Rev. E volume 66, pages 066110 (year 2002)NoStop
[Javadi-Abhari et al.(2024)Javadi-Abhari, Treinish, Krsulich,
  Wood, Lishman, Gacon,
  Martiel, Nation, Bishop,
  Cross, Johnson, and Gambetta]Qiskit
  author author A. Javadi-Abhari, author M. Treinish, author K. Krsulich,
  author C. J. Wood, author J. Lishman, author
  J. Gacon, author S. Martiel, author P. D.  Nation, author L. S.  Bishop, author A. W.  Cross, author B. R. Johnson, and author J. M. Gambetta, http://arxiv.org/abs/2405.08810 title
  Quantum computing with Qiskit (year 2024), https://arxiv.org/abs/2405.08810 arXiv:2405.08810 NoStop
[Zhai et al.(2023)Zhai,
  Larsson, Lee, Cui,
  Zhu, Sun, Peng, Peng, Liao, Tölle, Yang, Li, and Chan]Zhai2023
  author author H. Zhai, author H. R. Larsson,
  author S. Lee, author
  Z.-H. Cui, author T. Zhu, author C. Sun, author L. Peng, author R. Peng, author
  K. Liao, author J. Tölle, author J. Yang, author S. Li, and author G. K.-L. Chan, title title Block2 : A
  comprehensive open source framework to develop and apply state-of-the-art
  DMRG algorithms in electronic structure and beyond, https://doi.org/10.1063/5.0180424 journal journal J. Chem. Phys. volume 159, pages 234801 (year 2023)NoStop
[Feiguin and White(2005)]Feiguin2005
  author author A. E. Feiguin and author S. R. White, title title Time-step targeting
  methods for real-time dynamics using the density matrix renormalization
  group, https://doi.org/10.1103/PhysRevB.72.020404 journal journal Phys. Rev. B volume
  72, pages 020404 (year 2005)NoStop
[Ronca et al.(2017)Ronca,
  Li, Jimenez-Hoyos, and Chan]Ronca2017
  author author E. Ronca, author Z. Li, author C. A. Jimenez-Hoyos, and  author G. K. L. Chan,  title title Time-Step Targeting Time-Dependent
  and Dynamical Density Matrix Renormalization Group Algorithms with ab Initio
  Hamiltonians, https://doi.org/10.1021/acs.jctc.7b00682 journal journal J. Chem. Theory Comput. volume 13, pages 5560 (year
  2017)NoStop


