State Preparation in the Heisenberg Model through Adiabatic Spiraling

An adiabatic state preparation technique, called the adiabatic spiral, is proposed for the Heisenberg model. This technique is suitable for implementation on a number of quantum simulation platforms such as Rydberg atoms, trapped ions, or superconducting qubits. Classical simulations of small systems suggest that it can be successfully implemented in the near future. A comparison to Trotterized time evolution is performed and it is shown that the adiabatic spiral is able to outperform Trotterized adiabatics.


Introduction
Quantum simulations of lattice gauge theories, such as quantum chromodynamics (QCD) and quantum electrodynamics (QED), are anticipated, in the future, to enable reliable predictions for non-equilibrium dynamical processes, ranging from fragmentation in high-energy hadronic collisions, through to transport in extreme astrophysical environments. While a quantum advantage has yet to be established for a scientific application, including quantum field theories, there are substantial efforts underway to perform quantum simulations that can be compared with experiment, or impact future experiments, and impressive progress has been made toward these objectives in the last decade. This includes the development of techniques to simulate abelian gauge theories , non-abelian gauge theories , fermionic field theories [98][99][100][101], and scalar field theories [102][103][104][105][106][107][108]. There has also been the development of techniques to extract observables of interest to nuclear physics [109][110][111][112][113][114][115], scattering processes in high energy physics [116][117][118][119][120][121][122][123] and methods to mitigate errors on noisy quantum hardware [124][125][126][127][128][129]. Currently available hardware and our present understanding of quantum algorithms has so far limited quantum simulations of lattice gauge theories to one-and two-dimensions with only a small number of lattice sites [66,71,73,83,85,[91][92][93]. Qualitative insights can be gained from simulations of spin models that share one or more features of QCD, QED or low-energy effective field theories (EFTs) relevant to nuclear and particle physics. These include models that are in the same universality class as these theories, that can be fruitfully digitized onto qubit registers or mapped to analog quantum simulators, such as arrays of Rydberg atoms.
The Heisenberg model with arbitrary couplings is computationally universal in the sense that all other lattice models can be simulated in arbitrary dimensions, particle content and interactions by simulations of Heisenberg models [130]. Therefore, detailed understandings of quantum simulations of the Heisenberg model inform the simulations of quantum field theories describing the forces of nature. Translating results obtained from lattice field theories to predictions that can be compared with experiment requires that all relevant physical length scales are much larger than the scale of discretization of spacetime, and universality guarantees that low-energy continuum physics can be reproduced from simulations that are tuned near a 2nd order critical point [131][132][133]. As an example, it has been proposed that universality assures that the continuum physics of the 1 + 1d O(3) NLσM, which has been studied as a toy model of quantum chromodynamics (QCD) due to sharing a number of qualitative features such as asymptotic freedom, dynamical transmutation, the generation of a non-perturbative mass gap and non-trivial θ vacua, can be recovered from simulations of an anti-ferromagnetic Heisenberg model [134][135][136][137][138][139][140][141]. Thus, quantum simulations of the low-energy dynamics of the anti-ferromagnetic Heisenberg model, which requires preparing a low-energy state and evolving it forward in time, are expected to provide key insights into strategies for simulating QCD, including state preparation.
To enable practical quantum simulations of physical systems, preparation of states that have energies much less than the inverse lattice spacing is required. One proposal for preparing low energy states in both digital and analog quantum simulation is adiabatic switching. This works by beginning in the ground state of a known Hamiltonian and slowly varying the Hamiltonian through a path in parameter space where the energy gap does not close. Implementation of adiabatic switching in a quantum simulation requires the ability to prepare the eigenstate of the initial Hamiltonian and simulate time evolution. Schemes for simulating the Heisenberg model's time evolution have been proposed using digital quantum simulation [142,143], hybrid digital-analog simulation [144], periodically driven trapped ions [145], global microwave pulses on Rydberg atoms [146], nuclear spins [147,148], and adding strong single qubit terms to systems described by Ising interactions [149].
In this work, we present an argument that the eigenstates of the Heisenberg model are approximate eigenstates of the Ising model with a strong external field pointed in an appropriate direction. This is used to develop an analogue quantum simulation technique, called the adiabatic spiral, to adiabatically prepare the ground state of the Heisenberg model by adiabatically varying the direction of the external field in the Ising model. The feasibility of implementing the adiabatic spiral on Rydberg atoms and D-Wave's quantum annealer is investigated. It is found that current Rydberg atom experiments have sufficient coherence time and drive field strengths to implement the adiabatic spiral, while the D-Wave quantum annealer suffers from some limitations.  t T and t varies from t = 0 to t = T . During the evolution, e(t) is precessing around h(t) while the opening angle θ changes adiabatically, resulting in a spiral path on the unit sphere in the interaction picture. Vectors from the origin indicate the direction at the end of the spiral evolution.

Adiabatic Spirals: Spiraling Toward Ground States
We recently showed that the Ising model with large external fields in the transverse and longitudinal directions can approximate the time evolution of the Heisenberg model at discrete time intervals [149]. This generalizes previous work showing the Ising model with a large transverse field approximates the dynamics of the XY model [150][151][152][153][154], and is related to techniques used to study pre-thermalization [155][156][157][158]. Explicitly, if a quantum simulator evolves under the Hamiltonian then the time evolution ofĤ will be approximated at times that are integer multiples of t = 2π Ω , up to a change of basis and corrections that are O 1 Ω . In this work,X,Ŷ andẐ refer to the respective Pauli operators. We now observe that if the time evolution was reproduced exactly at these time intervals, it would guarantee thatĤ Heis. andĤ Ising share the same eigenstates and that their energy levels agree up to integer multiples of Ω. This would suggest that by beginning with θ = 0 and adiabatically increasing θ, it should be possible to prepare an eigenstate of the Heisenberg model from an eigenstate of the PauliẐ operators. When viewed in the interaction picture where the free part of the Hamiltonian is taken to be the single spin driving terms the localẐ operator becomesẐ where σ j is a vector of Pauli matrices and e(t) is a unit vector. For the Hamiltonian in Eq. (1), e(t) rotates in a circle about a vector pointing in the cos θẑ + sin θx direction. As θ is adiabatically varied to prepare the ground state of a Heisenberg model, e(t) will move in a spiral motion along the surface of the sphere as shown in Fig. 1. For this reason we call this method of state preparation adiabatic spiraling. It is important to note that the implementation of the adiabatic spiral does not require the switching time to be an integer multiple of 2π Ω . Studying the evolution of the Ising model at periodic times was only necessary to argue that the eigenstates of the Heisenberg model are also eigenstates of the Ising model up to O 1 Ω corrections. Typically in quantum simulation, adiabatic switching is done between ground states of gapped Hamiltonians. In contrast, the adiabatic spiral can be understood as performing adiabatic switching in the middle of the spectrum of the Ising model to prepare eigenstates of the Heisenberg model.
In practice, the initial eigenstate of the Ising model will often be degenerate, and this degeneracy will need to be split for the adiabatic approximation to be valid. This can be done by modifying the single-spin terms in the Hamiltonian. Explicitly, at times that are integer multiples of t = 2π Ω , the time evolution generated bŷ can be approximated by the Floquet operator, whereÛ B is a local change of basis given byÛ B = j e i θ 2Ŷ j . The additional single-qubit term can be tuned to create an energy penalty that breaks the degeneracy of the initial Hamiltonian which enables the application of the adiabatic approximation. With these additional single-qubit terms present, the same arguments at large Ω can be used to justify the applicability of the adiabatic spiral.
As an explicit demonstration, we consider the preparation of the ground state of the anti-ferromagnetic Heisenberg model on a 1D chain, with a Hamiltonian of the form Preparing the ground state of this system with the adiabatic spiral will require beginning in an eigenstate of the Ising model that is adiabatically connected to the ground state of the Heisenberg model. Equation 2 reproduces Eq. 7 when θ = arccos 1 √ 3 and would suggest that the ground state of Eq. 7 is connected to the ground state of Eq. 2 at other values of θ. The ground state of the Heisenberg Hamiltonian with θ = 0 is a state with spins alternating up and down in theẑ direction (a Néel state), e.g., |↑↓↑↓↑↓ ... . This ground state is degenerate and the degeneracy can be split by adding a single qubit term to the Hamiltonian with alternating signs. Therefore, the ground state of the full Heisenberg model can be prepared by beginning in a Néel state and applying a time-dependent Hamiltonian of the form Explicitly, if |V ac is the ground state of the anti-ferromagnetic Heisenberg model given in Eq. (7) and |Néel is a Néel state, then up to O 1 Ω and finite time corrections whereĤ(t) is the time-dependent Hamiltonian defined in Eq. (8). Typically, analog quantum simulators are initialized with all qubits in their ground state, e.g., | ↓ ⊗n . However, to apply the adiabatic spiral to the anti-ferromagnetic Heisenberg model, relevant for simulations of the O(3) NLσM, the initial state should be a Néel state, which can be accomplished by applying a rotation on every other qubit. Alternately, computations could be performed in a different basis. If X is defined to be a product of Pauli X's on every other site such that |Néel = X |↓ ⊗n where |↓ ⊗n is the state with all spins down, then Eq. (9) can be written as Multiplying both sides of this equation by X yields If the quantum simulator can directly implementH(t), then an adiabatic spiral can be used to adiabatically prepare the Heisenberg ground state from the |↓ ⊗n state up to a basis transformation. However, on some analog quantum simulators such as Rydberg atom systems, the sign of the two-spin interaction is fixed to be positive. This does not present an issue in using an adiabatic spiral, as it can be implemented with −H(t). This is because a change of overall sign does not change the eigenstates or presence of an energy gap between eigenstates, indicating that the adiabatic approximation remains valid.

A Numerical Example
To implement the adiabatic spiral in Eqs. (8) and (12), specific choices have to be made for h P (t) and f (t). As an example, we consider the Heisenberg comb that has recently been shown to reproduce the O(3) NLσM in the continuum and infinite-volume limits. It has modest qubit requirements, making it a good candidate for near term quantum simulations [137,138]. The Hamiltonian is given bŷ where S x,y = 1 2 σ x,y is the vector of Pauli matrices divided by 2 at position (x, y) on the lattice. An adiabatic spiral can be used to prepare the ground state of this model from a Néel state with the Hamiltonian The simplest choices for f (t) and h P (t) are linear functions of t. Once the functional forms of f (t) and h P (t) have been chosen, implementing the adiabatic spiral further requires choices of Ω, T , and h P (0). To ensure the eigenstates of the Ising model are as close to eigenstates of the Heisenberg model as possible, and the conditions of the adiabatic theorem are satisfied, Ω and T should be taken to be as large as possible. As an example, Fig. 2 shows the energy of a state obtained using the adiabatic spiral as a function of Ω for a Heisenberg comb of length 4 with J = J P = 1 and a switching time of T = 25. In this calculation, h P (t) = 0 for all times and f (t) was taken to be a linear function. Evolution underĤ(t) was evaluated numerically by computing and increasing N until convergence. At small values of Ω, the eigenstates of the Ising model and Heisenberg model are not close and the adiabatic spiral fails. At large Ω, the adiabatic spiral is able to prepare a state with an energy below the energy of the first excited state. The reason the energy of the state in Fig. 2 saturates above the ground-state energy is due to a combination of the amount of time used in the switching and the initial degeneracy in the ground state. The energy of the state prepared by the adiabatic spiral can be lowered by taking a non-zero value of h P (0). Unlike with Ω and T , the energy of the state prepared by an adiabatic spiral is not monotonic in h P (0). If h P (0) is taken to be too large, the process of switching off the initial energy penalty in a finite amount of time can break adiabaticity, leading to a state with larger energy being prepared. An optimal value of h P (0) can be found variationally. Fig. 3 shows the energy obtained from the adiabatic spiral with Ω = 8 and h P (0) = 0.18. This value for h P (0) = 0.18 was selected by minimizing the energy obtained by performing the adiabatic spiral with Ω = 8 and T = 25.
The adiabatic spiral's performance can be improved further by optimizing the path taken through parameter space. If h p (t) is taken to be a linear function of time, and f (t) taken to be an optimal path through parameter space can be found by minimizing the energy obtained as a function of the β n s [159]. Note that if the maximum driving field strength is limited on the analog simulator, this will constrain the values that the β n s can take. For long switching times, assuming that the maximum value of the driving field is obtained at the end of the spiral and truncating at β 1 , it was found that the energy obtained by the adiabatic spiral was monotonic as a function of β 1 . This means that in practice, the largest value of β 1 allowed by the constraint should be used so that the maximum value of the driving field is obtained at the end of the spiral, which is β 1 = 1 π . The energy of the state obtained by performing the adiabatic spiral with this path is also shown in Fig. 3. At long switching times, this path can offer an improved performance over linear switching. It was also found that including more terms in the expansion in Eq. (15) offered minimal improvements in the performance of the adiabatic spiral.

Rydberg Atoms
Arrays of Rydberg atoms are an experimental platform that could be used to implement the adiabatic spiral. If the atoms in the Rydberg state interact through the Van der Waals interaction, the Hamiltonian describing their dynamics isĤ wheren i is the Rydberg state occupation at site i, x i is the position of site i andX i couples the ground state of the atom at site i to its Rydberg state. This native implementation of Ising interactions has enabled the use of Rydberg atoms to perform analog simulations of the Ising model in 1D and 2D [160][161][162][163][164][165][166][167][168][169]. Their native Ising interactions allows for the implementation of adiabatic spirals. This is not the only technique that can be used to simulate the Heisenberg interaction on Rydberg atoms. In recent work, a method of simulating Heisenberg interactions on Rydberg atoms making use of dipole-dipole interactions and an external microwave field was proposed [146]. This approach uses global microwave pulses to rotate the dipole-dipole interaction to generate time evolution approximating an XXZ Heisenberg Hamiltonian, analogous to how one could Trotterize the Heisenberg interaction on a digital quantum computer. While this proposal makes use of dipole-dipole interactions, the same approach could be used for Rydberg atom simulations that use the Van der Waals interaction to couple atoms. In that experimental setup, an effectiveẐ iẐj is available. In principle, the ground state of a Heisenberg model could be prepared on this hardware using the adiabatic spiral or by using global pulses to perform a Trotterized adiabatic switching. The choice of which method to use will depend on the specific analog simulator being used. As a case study, we can consider hardware parameters used in a recent experiment to study phases of the 2D Ising model [169]. This experimental setup used a driving field with a max amplitude of Ω max = 2π × 4.3MHz and had a max coherence time of 3µs. With these hardware parameters, the preparation of the ground state of a Heisenberg comb with the adiabatic spiral parameters from the previous section can be performed by placing atoms in an array as shown in Fig. 4 with a lattice spacing of a ≈ 10.5µm. Note that the Hamiltonian simulated will not quite be the same as the Hamiltonian from the previous section due to the 1 r 6 coupling between the atoms, but the alternating vertical positions of atoms will ensure the next-to-nearest neighbor interactions will be suppressed. This same arrangement of atoms could be used to implemented Trotterized adiabatics. In this approach, π 2 pulses are alternated with periods of the external field turned off to construct a first order Trotterized approximation of time evolution according toĤ where (−1) i represents an alternating pattern of signs that breaks the initial degeneracy of the Hamiltonian. The errors in this state preparation method come from a finite switching time, standard Trotterization errors and the amount of time required to implement the π 2 pulses. This is similar in spirit to previous proposals to simulate time evolution using global controls [170]. The explicit pulse sequences for first order Trotterization can be found in Appendix A. Figure 5 shows the energy that can be obtained using the path optimized adiabatic spiral and Trotterized adiabatics with different numbers of Trotter steps as a function of total runtime on the analog simulator. At short times, a small number of Trotter steps are able to outperform the spiral. However, if the number of Trotter steps is fixed and the switching time is increased, the Trotterization error will grow and the energy of the state obtained eventually will increase. This can be mitigated by performing more Trotter steps, but as the number of Trotter steps increases the contribution to the error from the finite pulse time will accumulate. As shown in Fig. 5, this effectively puts a limit on the total number of Trotter steps that can be performed before adiabatics is no longer being effectively simulated. As shown in Appendix B, some of these errors can be cancelled at leading order with a modified pulse sequence. The dashed curves in Fig. 5 show the energy of the state that can be obtained by mitigating these errors at leading order in Trotterized adiabatics. The adiabatic spiral does not suffer from either of these issues. Increasing the amount of time used in the adiabatic spiral only increases the overlap with the ground state. It should also be noted that in the Trotterization simulations it was assumed that the drive field can be instantly ramped to its maximum value. In practice, this is not possible and the shape of the pulse while ramping on the drive field will contribute to errors in the time evolution. Given these limitations, it appears that the adiabatic spiral is a better choice for preparing the ground state of a Heisenberg model than Trotterized adiabatics on a Rydberg system. It should also be noted that the pulse sequences required to implement the adiabatic spiral are quite similar to other pulse sequences that have been implemented on Rydberg atoms [161,166,168,169]. As such, it is expected that the adiabatic spiral will be robust against experimental errors when implemented on Rydberg atoms.

D-Wave
Another experimental platform that could potentially implement the adiabatic spiral is D-Wave's quantum annealer. The D-Wave quantum annealer implements time evolution according to the Hamiltonian, where A(s) and B(s) are fixed functions, but s(t), h i and J ij can be programmed by the user, up to some constraints. Although it is designed to solve optimization problems, D-Wave's hardware has been used to perform analog simulations of the Ising model [171][172][173][174][175][176][177][178][179][180][181]. In previous work, D-Wave's quantum annealer has also been used to simulate time evolution and perform state preparation for lattice gauge theories, however these mappings to quantum simulation have more in common with digital quantum simulation than the analog protocol we are proposing [72,85,92]. To implement the adiabatic spiral on the D-Wave annealer, a set of qubits with a comb pattern coupling must be selected. Given the connectivity of the hardware, this is straightforward to do. J ij and h i should be chosen so that the implemented Hamiltonian iŝ H DW (t) = − A(s(t)) 2 x,yX x,y + B(s(t)) 2 x,y The spiral can be implemented by using D-Wave's reverse annealing function to begin in a Néel state with A(s(0)) = 0, and adiabatically evolve to s * such that A(s * ) = √ 2B(s * )h. As long as J, J p h, the adiabatic spiral will be able to prepare an equal superposition of the ground state and the first excited state of the Heisenberg comb. The pure ground state will not be prepared because an energy penalty cannot be implemented to break the degeneracy within the initial state. Note that to perform a measurement on the D-Wave QPU, the pulse schedule must end with A(s = 1) = 0. The process of turning off A(s) will create errors in the state preparation whose size are controlled by the slew rate of A(s) and the choice of s * . One choice of time dependence of A(t) and B(t) to implement the adiabatic spiral is shown in Fig. 6. As an explicit example, we once again consider the Heisenberg comb of length 4 with J = J p . The energy of states prepared with classical simulations (Mathematica and python) of the D-Wave Advantage 6.1 hardware for different choices of J and h are shown in Fig. 7. For small values of J and h, simulations of the D-Wave hardware are able to produce states with a lower energy than the initial state. However, the quantum device has calibration errors on the order of δJ ≈ 4 * 10 −3 , making it unlikely that these parameters can be used in practice. For larger values of J, the simulations cannot reach states with energy lower than the initial state. As an example, the energy of states as a function of time for J = 0.1 and h = 2 found in simulations of Advantage are shown in Fig. 8. As seen, Advantage is expected to be able to prepare a low energy state with respect to the Heisenberg Hamiltonian part-way through a simulation. However, in the process of turning off A(s) as required to perform a measurement, the time evolution is still approximately adiabatic, and the system ends up in a final state with large overlap with the initial state. If the slew rate of Advantage could be increased, or if measurements could be performed without fully turning off A(t), it would be possible to perform measurements of the lowenergy state prepared by the adiabatic spiral and potentially use Advantage in a scalable manner for quantum simulation of the Heisenberg model. Implementing the adiabatic spiral on Advantage also suffers from the complication that the single qubit terms must be tuned to a specific value to successfully implement the spiral. While a pulse schedule is provided by D-wave for this hardware, the pulses actually implemented can vary from this schedule by up to 30%. This is not an issue for the classical optimization problems the Advantage quantum annealer is designed to solve, but represents a significant issue for reliable implementation of the adiabatic spiral on this hardware.

Discussion
In this work, a method of simulating the Heisenberg model on quantum hardware with native Ising interactions was extended to perform adiabatic state preparation of the Heisenberg ground state. Due to the number of quantum simulation platforms able to natively implement the Ising model, this technique is likely to have experimental applicability. It was shown that it should be feasible to implement this technique on Rydberg atoms, and that it can out-perform a Trotterized approach to adiabatic state preparation. The feasibility of implementing the adiabatic spiral on D-Wave quantum annealers was also investigated, and it was found that the requirement of turning off theX driving field at the end of the calculation limits its efficacy. The ability of the adiabatic spiral to prepare low-energy states of the Heisenberg model on these hardware platforms suggests that near-term quantum simulators can be used to perform analog quantum simulations of the O(3) NLσM, and other such theories of importance to nuclear and particle physics. While the case of the fully symmetric XXX Heisenberg model was focused on in this work, the adiabatic spiral can also be applied to the XXZ Heisenberg model as well. As the approximation of the dynamics of the Heisenberg model by an Ising model with strong drive fields is a special case of the more generic phenomena of prethermalization [155][156][157][158], similar ideas may find applicability for other Hamiltonians on analog quantum simulators with different native couplings.

A First Order Trotterization with Rydberg Atoms
In this section, it will be shown how to engineer a Trotterized approximation to an adiabatic turning on of the Heisenberg interaction on an array of Rydberg atoms. Assuming access to a global and staggered driving field, the Hamiltonian describing the time evolution of an array of Rydberg atoms is given bŷ It will also be assumed that the strength of the global drive field is limited, | ω(t)| ≤ Ω, as is the case in actual Rydberg atom experiments. The approximation to Heisenberg time evolution will be assembled from the following set of global analog gates, R ± X ( ) and R ± Y ( ) correspond to global π 2 -pulses about the x and y axes respectively and can be generated using the global drive field. R Z (t, κ) can be generated by only turning on the staggered driving field. Due to the maximum driving field strength, Ω, the π 2 -pulses, R ± X,Y , need a minimum device time = π 2Ω during which theẐ iẐj interaction cannot be "turned off." This leads to O( ) cross-talk errors that will be studied in further detail in Appendix B. With this analog gate set, a first order Trotter approximation to a generic XXZ Heisenberg evolution with a staggered single qubit term is given by Note that since Ẑ j commutes with the XXZ Heisenberg Hamiltonian, the product overẐ j 's can be neglected when performing adiabatic state preparation as its only contribution is to change the overall phase of the prepared eigenstate. With this pulse sequence, the time evolution generated by the Hamiltonian in Eq. 17 can be approximated by

B Cross-talk Mitigated Trotter Sequence
The accuracy of the Trotterized approximation to the time-evolution operator can be improved by using shorter Trotter steps or higher order formulas. However, for Trotter step sizes comparable to the π 2 pulse length, it becomes important to compensate for the O( ) contributions to the error. In this section, it will be shown how a modification of the pulse sequence in Appendix A can be performed to cancel this error at leading order. To leading order in , the π 2 rotations are given by Using these expressions, it can be shown that at leading order in the pulse sequence from Appendix A is By shifting the length of pulses in the sequence from Appendix A, the O( ) terms can be absorbed into the time evolution. This motivates the improved adiabatic Trotter sequence given by where θ(x) is the Heaviside theta function. When T N 2 is larger than , this pulse sequence cancels the crosstalk errors at O( ). When the Trotter step size is taken to be smaller, this pulse sequence will only cancel the cross-talk errors for the later steps in the sequence.