Research articlesSimulations of frustrated Ising Hamiltonians using quantum approximate optimization
Abstract
Novel magnetic materials are important for future technological advances. Theoretical and numerical calculations of ground-state properties are essential in understanding these materials, however, computational complexity limits conventional methods for studying these states. Here we investigate an alternative approach to preparing materials ground states using the quantum approximate optimization algorithm (QAOA) on near-term quantum computers. We study classical Ising spin models on unit cells of square, Shastry-Sutherland and triangular lattices, with varying field amplitudes and couplings in the material Hamiltonian. We find relationships between the theoretical QAOA success probability and the structure of the ground state, indicating that only a modest number of measurements () are needed to find the ground state of our nine-spin Hamiltonians, even for parameters leading to frustrated magnetism. We further demonstrate the approach in calculations on a trapped-ion quantum computer and succeed in recovering each ground state of the Shastry-Sutherland unit cell with probabilities close to ideal theoretical values. The results demonstrate the viability of QAOA for materials ground state preparation in the frustrated Ising limit, giving important first steps towards larger sizes and more complex Hamiltonians where quantum computational advantage may prove essential in developing a systematic understanding of novel materials.
This article is part of the theme issue ‘Quantum annealing and computation: challenges and perspectives’.
1. Introduction
Quantum magnetism has been a major focus in condensed matter research, driven by the potential for new disruptive applications ranging from quantum computing to quantum sensing [1]. Quantum material properties are intrinsically related to the structure of the ground states. However, exact ground states are notoriously challenging to calculate classically, requiring the field to resort to using semi-classical limits [2–5] or fully quantum approaches with restricted applicability [6–12]. New, fully quantum computational tools are required to understand current problems including frustrated two-dimensional quantum magnets currently explored by bulk neutron scattering and thin film susceptibility [13,14]. Digital and analogue quantum simulators have emerged as a new tool for the simulation of quantum many-body phenomena towards efficient modelling of exotic quantum phases of matter beyond classical tractability [15,16]. They are naturally suited for magnetic Hamiltonians since spins can be directly mapped to qubits. Non-trivial phases in magnetic systems, such as frustrated phases [17], spin glasses [18], and topologically ordered phases [19,20], have been realized on multiple qubit platforms using a variety of techniques.
In this paper, we investigate an alternative approach to preparing materials ground states using the quantum approximate optimization algorithm (QAOA) [21] on near-term quantum computers. We apply QAOA to lattices of interest in materials science, considering the classical Ising limit (equivalently, ) where standard QAOA is directly applicable. This serves as a stepping stone towards truly quantum problems such as the and Heisenberg models in the fully frustrated limit, which will require further algorithmic research and modifications to the approach presented here. Our results validate that QAOA achieves sufficient accuracy for the simpler classical limit and provides insights into algorithmic behaviour for material lattice problems.
We consider lattice instances with varying degrees of frustration. The smallest building block of a frustrated magnetic Hamiltonian is an anti-ferromagnetic triangular motif of three spins where all the bonds cannot be satisfied simultaneously. In these materials, exchange interactions compete such that it is impossible to satisfy them all simultaneously (figure 1). If all spin configurations are equally favourable, frustration can lead to non-ordered states such as spin liquids [5], spin glasses [22] or plaquette states [8], each with distinct signatures.
Figure 1. Example of frustration on an anti-ferromagnetic triangular motif. Two spins in opposite orientations (black and green) minimize the energy along one bond, however, there is no configuration for the final spin that minimizes energy along both remaining bonds. (Online version in colour.)
We solve three different types of Hamiltonians for unit cells pictured in figure 2. The first is a square unit cell Hamiltonian, which exhibits only simple ferromagnetic and anti-ferromagnetic phases in the infinite size (thermodynamic) limit. The second is the celebrated Shastry-Sutherland lattice. Interestingly, this problem already lends itself to materials applications and experimental data analysis. Among other examples, it is conjectured to describe the class of rare-earth tetraborides (, and ) and allows a direct comparison with several existing results both theoretical [23–27] and experimental [28–33]. The third case is the more complex Ising triangular lattice which represents a maximally frustrated problem with an infinite number of possible ground states in the infinite size limit [34,35]. We compute theoretical probabilities to prepare the ground state for each of these 9-spin Hamiltonians under varying choices of the external field and coupling parameters and compare these theoretical results against computations on a trapped ion quantum computer.
Figure 2. Unit cells of (a) square, (b) Shastry-Sutherland, and (c) triangular geometries. Colours indicate two spin values in examples of ground states with . (d–f) Phase diagrams for each of the unit cells (a–c), respectively, with labels ‘A’, ‘B’, denoting regions with distinct ground states for each lattice. Magnetization at is due to degeneracy in the ground states, where spin-flip-related pairs of states are present in the absence of the field (. (Online version in colour.)
We choose spins as this is the logical minimum number of spins required to construct a unit cell of the Shastry-Sutherland lattice and also it is a feasible size for the Quantinuum quantum computer (with 12 qubits available at the time this work was completed). We focus on layers of the QAOA algorithm; for larger instances more layers of QAOA will be needed to maintain a significant success probability [36–38]. Implementations on quantum computers will also have to overcome predicted limitations due to noise [39–42] including an exponential scaling in the number of measurements with circuit size, depth, and other factors [43]. Noise in modern quantum computers has negative consequences for all quantum algorithms, not only QAOA. Ongoing testing and development of these devices is necessary to assess realistic performance scaling in the presence of noise and to determine whether QAOA, or other quantum algorithms, will ultimately succeed in providing a useful computational advantage over conventional approaches.
2. Ising Hamiltonian and model unit cell lattices
A single unpaired spin on the outermost orbital of a magnetic ion constitutes a state which is implemented straightforwardly on a qubit. In a magnetic material, several such spins in a lattice interact via pairwise superexchange interactions . The nature and strength of these interactions, , depend on several factors. These include the distance between the magnetic ions (typically scales as the inverse cube of the distance between the ions), the shape of the orbitals, the symmetry of the lattice and the local crystal fields. Magnetic frustration can arise in the lattice, for example, if spins arranged on a triangular motif in the lattice experience equal Ising anti-ferromagnetic interactions. Such a magnetic frustration can arise via a combination of straight edge bonds which are either horizontal or vertical, and the diagonal bonds . This is given by the Hamiltonian
Multiple methods have been proposed to solve for ground states of Ising Hamiltonians and related optimization problems, notably Integer Programming method [46], Simulated Annealing [47] and its variants, Large Neighbourhood Search [48] and Quantum or Quantum-inspired physical annealing devices. Among them, the Integer Programming method solves exactly but suffers from exponential scaling of computational time. Simulated Annealing and Large Neighbourhood Search are both herustics methods that promise faster runtime but there is no guarantee of the solution qualities. Quantum annealers, digital annealers and coherent Ising machines are hardwares dedicated to solving Ising models [49,50]. Depending on the connectivity of these various (qu)bits, every backend could be good at a different task—parallel tempering machines could be good at finding the classical phases [26], while others could reveal intricate dynamical behaviour in a transverse field Ising universality [51]. For frustrated lattice problems, QAOA allows us to sample different ground states with certain probabilities due to quantum randomness, whereas classical and deterministic algorithms may generate a ground state efficiently, but fail to explore the states which might arise because of a coherent superposition between all the spins. Additionally, future extensions beyond the Ising limit also become an exciting possibility.
(a) Ground state magnetization phase diagrams
We consider the nine-spin unit cells with geometries in figure 2a–c which represent the number of spins required to simply construct a unit cell of the Shastry-Sutherland lattice. In materials represented by Bravais lattices, these unit cells repeat periodically to realize the very large lattices in a real material. We computed ground states for each unit cell by evaluating (2.1) for each possible spin configuration to identify the lowest energy states, for varying choices of and , with taken as the unit of energy. We plot the magnetization
Starting with the non-frustrated square lattice with —only interactions with simple ferromagnetic or antiferromagnetic ground states, the degree of frustration is tuned progressively by (i) addition of bonds and (ii) bringing . The triangular lattice with uniform coupling parameters represents the maximally frustrated limit with highly degenerate solutions. The Shastry-Sutherland model represents a scenario with the minimum number of bonds required to realize a fully frustrated lattice. The solution of these states represents a problem of polynomial time complexity in two dimensions and without a magnetic field.
The ground state for a given shows a number of magnetization plateaus, where each plateau has a different proportion of spins pointing up. Unsurprisingly, at large , the ground state for each lattice is ferromagnetic, in regions C, G and E in figure 2d–f, respectively. The situation becomes more interesting at small in regions ‘A’, the ground state is anti-ferromagnetic with magnetization , as five spins are aligned with the field while four are anti-aligned. For fields and small , there is a ground state with in which a single spin in the centre of the unit cell is anti-aligned with the field. Besides frustration, these states are also determined by the finite size of the unit cells, where the central spin is distinguished as the only spin with four interactions in the square lattice. As and are varied, frustration leads to a variety of different ground states for the Shastry-Sutherland and triangular lattices, with varying magnetizations in figure 2e,f, with ground states in electronic supplementary material. These are true ground states of the 9-spin Hamiltonians, with boundary spins playing a big role. In the infinite size limit, we expect the ground states to be progressively less dependent on the boundary, and more on the symmetry of the Hamiltonian, which we discuss in the next section.
(b) Finite size effects
The finite sizes of our lattice unit cells, as well as the unusual ground state noted in the previous section, raise questions of how the ground states for our unit cells match with ground states that would be obtained in the large size limit, and the minimum number of spins that are needed to achieve quantitative behaviour consistent with large sizes. To address these questions, we computed the magnetization of triangular and Shastry-Sutherland lattices of spins to analyse the size-dependent behaviour. Owing to the exponential complexity of the problem, we used neal [52], a software implementation of simulated annealing to approximate the ground state configurations with and . Each combination of and was run 50 times and the solution with minimum Ising energy was picked. Examples with are shown in figure 3. We assessed convergence of the global phase diagrams to the large size limit by computing the root mean square error (RMSE) between the target lattice’s and 3030 lattice’s magnetization across points in the phase diagram. We fitted the RMSE to both a power-law as well as an exponential form (see electronic supplementary material, figure S1), and we find that the power-law scaling with exponent and prefactor fits the RMSE better. The equation takes the form:
Figure 3. Ground-state magnetization of triangular spin arrays with a number of spins per dimension (a) ,(b) , (c) , and (d) , computed as described in §2b. (Online version in colour.)
Figure 4. Root mean square errors (2.3) of total magnetization between lattices relative to a lattice, plotted on a log-log scale. Solid lines show fits to the power-law scaling relation (2.3); the slopes indicate the best-fit exponents and , with fit 0.994 and 0.990 for the triangular and Shastry-Sutherland lattices, respectively. Best-fit intercepts are and for the triangular and Shastry-Sutherland lattices, respectively. The fit and the errors in the exponent are based on standard Levenberg–Marquardt routines and assume Poisson statistics at each point. (Online version in colour.)
A more rigorous analysis of finite-size scaling [53,54] around each critical point could yield a careful analysis of the required lattice size for target fidelity for every phase transition. Our result based on an overall RMSE demonstrates that a spin grid is already obtaining results close to the much larger 3030 grid. Based on this trend, we expect that finite size effects in will diminish quickly with the size of the lattice, indicating that lattices of only a few hundred spins may diminish the errors sufficiently to achieve a realistic ‘bulk’, and therefore meaningful results for comparison with experiments which probe bulk properties, such as diffraction and heat-capacity. This suggests that quantum processors with hundreds of qubits, achievable within the noisy intermediate-scale quantum era [55], may be capable of meaningful applications for materials science applications.
3. Quantum approximate optimization algorithm
Quantum computers offer a route to overcoming issues associated with identifying ground states through conventional methods. One approach to address these problems uses the quantum approximate optimization algorithm, which was originally designed to find approximate solutions to difficult combinatorial optimization problems [21] that are often expressed in terms of Ising Hamiltonians [56]. Empirical performance of QAOA has been characterized for a variety of combinatorial problems [36,57–60] and this has also led to generalizations [61–66] that have been applied to preparing chemical ground states [67] as well as ground state preparation for one-dimensional [38,68] and two-dimensional [69] quantum spin models in theory and experiment [70].
To formulate our Ising problems in a structure that is suitable for QAOA, we first express the Ising Hamiltonian (2.1) in terms of a quantum Hamiltonian operator
To find solutions, QAOA uses a quantum state prepared with layers of unitary evolution, where each layer alternates between Hamiltonian evolution under the Ising Hamiltonian and under a ‘mixing’ Hamiltonian
Figure 5. Quantum-classical optimization loop for QAOA. For a given set of parameters , a quantum computer generates and measures states . The measurements are sent to a conventional computer to compute and check its convergence. If is not converged, then an optimization routine selects new and for the quantum computer. If is converged, then the algorithm terminates and the final solution is the measured result that minimizes the energy.
An analytic proof has demonstrated that QAOA can prepare an exact ground state of Ising Hamiltonians as [21,56]. Apart from the formal proof of convergence at large , there has been significant interest in applying QAOA at small , where approximations exceeding conventional lower bounds have been observed in simulations [36,73] and predicted for large problems in certain contexts [74]. Realizing such advantages on devices with hundreds of qubits or more is an important topic of ongoing research as quantum computing technologies continue to develop.
For the materials problems we consider, we are interested in preparing the ground states of the Ising Hamiltonian. We compute exact ground states for our unit cells in figure 2 by evaluating all eigenvalues of the Hamiltonian using equations (2.1) and (3.1) to identify the lowest energy state. For some phases the number of ground states is due to degeneracies, while for other phases there is a single ground state , as pictured in electronic supplementary material. To assess QAOA performance, we compute the average ground-state probability
Ground-state preparation is a goal specific to the materials problem context we are interested in here. This is, importantly, a departure from the standard goal of QAOA in the context of approximate combinatorial optimization, where the goal is to find approximate solutions that are not necessarily the ground states. While QAOA is not expected to efficiently find exact ground states for generic NP-hard optimization problems, it may still prove useful for finding ground states of specific structured problems such as materials problems on a lattice similar to those we explore here [38,68–70].
(a) Numerical simulations of ideal QAOA
We use numerical calculations to assess the theoretical performance of QAOA for ground-state preparation. These demonstrate the ideal performance of QAOA in exact pure state calculations that use matrix multiplication to evaluate (3.3). This gives an ideal baseline for later comparison against results from a noisy quantum computer, where errors lead to mixed states with degraded performance. We use QAOA layers throughout this section and our results.
To identify QAOA states for our Ising problems, we must determine optimized QAOA parameters and . We choose regions to evaluate parameters in determining and as follows. QAOA is periodic when [57], hence we consider , which gives all unique up to symmetries. The periodicity of the parameter is more complicated, as it depends on the Hamiltonian in in (3.3). Here we focus on intervals near the origin and dependent on the magnitude of the Hamiltonian terms, which has been highly successful in previous work [75]. The basic idea is that the QAOA unitary changes at varying speeds with respect to , depending on the Hamiltonian coefficients , and . When the Hamiltonian coefficients increase, then should decrease to obtain a similar unitary. The rate at which the unitary changes with respect to is related to the average magnitude of the Hamiltonian coefficients
We show an example of how the energy expectation value and average ground state probability depend on the choice of parameters in figure 6 for an example with the Shastry-Sutherland unit cell with Hamiltonian coefficients and (similar patterns are observed in sample calculations using other choices of Hamiltonian coefficients and also for the triangular unit cell). There are two regions in figure 6a with optimized in yellow. The ground-state probabilities in figure 6b are also relatively large near the and that optimize .
Figure 6. Numerical simulations of the (a) average energy and (b) average ground-state probability with varying choices of QAOA parameters and , for the Shastry-Sutherland unit cell with and (§3a). Each plot has the same range of and ; the colour scales are reversed in (a) and (b) so that small and large are each represented by bright colours. (Online version in colour.)
We have found in searches over much larger intervals that the local optima for and do not always approximately align as in figure 6, which can lead to poor at optimized in these larger intervals. However, this does not appear to be a prevalent issue for the smaller -dependent intervals in cases we have looked at. The results are somewhat sensitive to the specific choice of interval, however, our choice of gives satisfactory results across the varying lattices.
To identify optimized parameters, we perform a grid search over and for each Hamiltonian considered. We evaluate the QAOA states in (3.3) on 201 evenly spaced intervals with and over 300 evenly spaced intervals with for a total of grid evaluations for each Hamiltonian. This approach gives optimal parameters in our intervals up to coarse graining in the grid search. We select parameters and that optimize .
The optimized parameters and do not necessarily give the optimal ground state probabilities that are possible from QAOA. The reason is that the average energy optimization accounts for energies and probabilities of all states, which together may yield low energies at parameter choices that are not optimal for the ground states alone [36]. To assess performance, we further compare our parameter choices against parameters and that directly optimize . The direct optimization of is used here for benchmarking purposes and is not a realistic approach for large problems where the ground states are unknown. For our small problems, the comparison gives an idea of how the ground state probabilities from a standard optimization of compare against the best ground state probabilities that could be obtained by QAOA in our set-up.
(b) Quantum computations of QAOA
We next investigate the performance of QAOA using the Quantinuum H1-2 quantum computer. H1-2 contains trapped-ion qubits and uses lasers to implement gates on these qubits. Typical error rates are reported as for two-qubit gates and for single-qubit rotation gates [76]. In addition to the device H1-2, we also use the H1-1E device ‘emulator’ to simulate noisy device behaviour. This gives results that approximately correspond to expected device behaviour while avoiding the financial expense and wait times that are associated with running the device. The emulator models a variety of device-specific noise processes for the H1-class computers, including depolarizing noise, leakage errors, crosstalk, dephasing in transport and qubit idling errors [77].
We test QAOA on the H1-2 using the QCOR software stack [78]. The QCOR stack translates the series of unitary operators expressing QAOA into quantum circuits for H1-2; see electronic supplementary material, appendix B, for details. The QCOR program used for submitting jobs to the device as well as our calculations are available online, cf. ‘Data accessibility’.
Furthermore, modern quantum computers are known to be affected by state preparation and measurement (SPAM) errors as well as gate infidelities from a variety of physical sources. We assessed SPAM errors expected in our quantum computations using the device emulator, with details in electronic supplementary material, appendix C. The probability to observe no error in circuits we tested was approximately 96%, with errors distributed approximately uniformly across qubits. We account for these errors using an independent bit-flip model and associated SPAM matrix , which transforms an ideal set of measurement results to the expected noisy set of results. The inverse matrix can then be applied to our noisy measurements from the quantum computer to approximately correct for SPAM errors. A technical issue arises in that the mitigated measurement probabilities can sometimes be negative, due to the approximate nature of the mitigation scheme. This leads to a second mitigation scheme that additionally sets all negative probabilities to zero and renormalizes so the total probability is one. We use each of these approaches to attempt to correct the small SPAM errors we expect from the quantum computer, as described in detail in electronic supplementary material.
4. Results
In this section, we consider the results from QAOA applied to the materials lattices of figure 2. We take as the unit of energy and analyse the success of QAOA in preparing ground states at variable and , first in numerical simulations (§3a) and then in quantum computations on a trapped-ion quantum computer (§3b).
(a) Ground-state measurement probabilities
We first consider theoretical probabilities to measure the ground state with QAOA and how these vary for different parameter choices in the Hamiltonian. We begin with the simple square lattice in figure 2a,d, which does not exhibit frustration as there are no triangles in the interaction graph. The probability to measure the ground state for varying is shown in figure 7. Figure 7a shows the probabilities obtained from optimizing the standard objective , while figure 7b shows the best-case results based on optimizing , as described in §3a. The probabilities in each case are similar, demonstrating that optimizing is nearly as successful in increasing the ground-state probability as a direct optimization.
Figure 7. The square unit cell ground-state probabilities when (a) optimizing and (b) optimizing as described in §3a. The ranges are identical in each figure. Phases A, B, C refer to the anti-ferromagnetic, =7/9, and ferromagnetic phases of figure 2, respectively, with vertical dotted lines showing the phase boundaries.
The average ground-state probability shows distinct behaviours for each of the three ground states at varying , visually separated by dotted lines. In the anti-ferromagnetic ground state at small , the probability approximately oscillates between and , with small probabilities observed near integer values of and larger probabilities near and . The ground state with has a near-constant probability of . At the ground state becomes ferromagnetic and increases significantly, with monotonic increases at larger .
We rationalize the varying success probabilities in the figure as attributable to structures of the ground states at varying and the interplay with the structure of the QAOA state in (3.3). We show in electronic supplementary material, appendix D that QAOA can exactly prepare the ferromagnetic ground state when for arbitrary lattice sizes, based on the fact that in this limit. This is consistent with the behaviour in the figure, where increases monotonically with for the ferromagnetic ground state at . We further show in electronic supplementary material, appendix E how the anti-ferromagnetic ground state probability is maximized at , and we devise large parameters that can further improve these results. (We did not include larger parameters in our numerical searches as this can lead to poor at parameters that optimize for the frustrated lattices, as remarked in §3a.) However, the mechanism for anti-ferromagnetic ground state preparation here depends on the specific choice of the lattice, and it is not clear how QAOA will behave for other lattice sizes. For the M=7/9 phase, the QAOA state is more complicated, as it is a superposition of many basis states that depends on the optimized parameters, and we do not have an analytic account for this behaviour. The optimized parameters that create each QAOA state in the square-lattice phase diagram are shown in electronic supplementary material, appendix F.
Ground-state probabilities for the Shastry-Sutherland and triangular lattices are pictured in figures 8 and 9, respectively. Ground-state probabilities from optimizing are presented in panels (a) while panels (b) show the best case probabilities from a direct optimization of as described in §3a. These probabilities show patterned behaviour, with distinct probabilities observed throughout most of each individual region A,B,, with significant differences in between different regions. At small , there are oscillations in the probability for preparing the anti-ferromagnetic ground states at small , and large success probabilities for the ferromagnetic ground state at large , as foreshadowed by results from the square lattice. On the other hand, as the coupling increases, the triangular and Shastry-Sutherland lattices experience increased frustration, with competing interactions within the triangular motifs in figure 2. The average ground-state probability decreases significantly as increases and frustration becomes dominant. This is especially evident when , for example in the top left of each of figures 8b and 9b. The are mostly uniform across and within each region, qualitatively similar to the nearly uniform probability for the state at varying for the square lattice in figure 7. Ground-state probabilities are typically 0.01, indicating that only measurements are expected to identify a ground state. We now turn to computations on a trapped-ion quantum computer, to benchmark and assess performance of QAOA on a real quantum computing device.
Figure 8. Shastry-Sutherland unit cell ground-state probabilities when (a) optimizing and (b) optimizing as described in §3a. The ranges for and are identical in each figure. (Online version in colour.) Figure 9. Triangular unit cell ground-state probabilities when (a) optimizing and (b) optimizing as described in §3a. The ranges for and are identical in each figure. (Online version in colour.)
(b) QAOA quantum computations
Here we assess QAOA performance in preparing ground states on a trapped-ion quantum computer. Ultimately, our aim is to validate the idea that a current quantum computing technology is capable of preparing each ground state of a frustrated Ising Hamiltonian using QAOA. An important first step is to assess whether optimized parameters from our theoretical calculations are also optimized for the device, or whether further optimization is needed to determine device-specific optimized parameters.
(i) Quantum computational performance with varying parameters
QAOA depends on the choice of parameters, as discussed in connection with figure 6. To test whether our theoretical parameters also yield good performance in the device, we consider QAOA circuits evaluating a point in region E of the Shastry-Sutherland phase diagram figure 2e, with Hamiltonian coefficients and QAOA parameters shown in table 1. The parameters correspond to a local minimum in , similar to the minima observed in figure 6. We use the H1-1E device emulator to evaluate circuits at the optimized and and circuits where either or has been displaced from its optimal value, as shown in figure 10. Black crosses in the figure indicate how increases in pure state simulations as either of these parameters are varied individually. Error bars denote the analytic standard error of the mean (s.e.m.) for measurement shots, with s.e.m. calculated numerically from the pure states. If the quantum computations did not have any noise, then from the central limit theorem we would expect about two-thirds of the from the quantum computer to be within these error bars.
Figure 10. Angle sensitivity analysis for and , with separate variations in (a) and (b) about the ideal values from table 1. Black crosses show results from pure state calculations, with error bars denoting the standard error of the mean at 1000 shots (see text). Data points showing results from the Honeywell emulator are denoted with ‘E’ in H1-1E and results from the trapped-ion quantum computer are labelled H1-2. Data points labelled ‘H1-1E’ and ‘H1-2’ are raw data, labels ‘E.M.’ (error-mitigation) are with basic mitigation, and ‘E.M. ’ are readout error mitigation that forces each probability to be , as described in §3b. (Online version in colour.)
| region | degeneracy | ||||||
|---|---|---|---|---|---|---|---|
| A | 1 | 1/9 | 0.240 | 1.440 | 0.750 | 400 | |
| B | 4 | 1/9 | 3.840 | 0.480 | 0.112 | 1000 | |
| C | 4 | 3/9 | 3.840 | 1.680 | 0.121 | 1000 | |
| D | 2 | 3/9 | 1.680 | 1.920 | 0.131 | 400 | |
| E | 4 | 5/9 | 2.000 | 2.480 | 0.143 | 1000 | |
| F | 1 | 7/9 | 1.680 | 3.600 | 0.182 | 400 | |
| G | 1 | 1.0 | 0.240 | 5.520 | 0.244 | 400 |
The theoretical can be compared against the device emulator, with data point labels in the figure beginning with ‘H1-1E’. There are three sets of data points for the emulator; the first is direct output labelled ‘H1-1E’, the second includes SPAM error mitigation (§3b) in ‘H1-1E, E.M.’, the third includes a variation of the SPAM mitigation that additionally forces the mitigated probabilities to be in ‘H1-1E, E.M., ’. These emulated are larger than the theoretical values and we attribute this to noise in the device emulator, which introduces errors that cause the energy to deviate from its ideal minimum value. Despite these errors, the shape of the landscape is similar to our theoretical calculations, with best performance observed near the theoretical parameters that minimize , and energies that tend to increase away from these parameters.
We further validate that the H1-2 trapped-ion device itself is consistent with the emulator in the data points that begin with ‘H1-2’. These actually yield better energies than the emulator, and are within one standard error of the mean from our theoretical results. The results from the device and emulator indicate that the energy landscape as a function of the QAOA parameters and is consistent between our theoretical calculations, the quantum device, and emulator. We therefore proceed with our theoretically optimized parameters to evaluate success in ground-state preparation using the quantum computer.
(ii) Quantum calculations of ground states
We now perform calculations on the Honeywell H1-2 quantum computer to analyse success probability in ground-state preparation. We consider points in each region of the Shastry-Sutherland lattice, using parameters listed in table 1 that correspond to local minima in , similar to the minima used to evaluate theoretical performance in §4a). We post-process the measurement results using the SPAM mitigation model with probabilities (see §3b), to give a minor correction to the observed results that is designed to counteract this known source of error.
Figure 11 shows the ground-state probabilities from quantum computations in comparison with ideal expectations from pure states. The ground states are separated by regions A,B, with markers a,b, corresponding to the individual ground states pictured in electronic supplementary material, figure S16. The quantum computations succeed in observing each individual ground state in each region of the Shastry-Sutherland lattice, as seen by the positive probabilities in each state ‘a’, ‘b’,
Figure 11. Probabilities to observe each ground state from pure state simulations compared with observed frequencies estimated by quantum computations with the H1-2 device for each different phase (A-G) of the Shastry-Sutherland unit cell. Alphabetical labels ‘a’, ‘b’, etc., identify the different ground states in electronic supplementary material, figure S16. (Online version in colour.)
For a closer comparison of probabilities, we plot error bars denoting the theoretical standard error of the mean s.e.m.=. The s.e.m. defines a range in which we expect about two-thirds of estimated are expected to be found, where is the number of measurement results of a given ground state and is the total number of measurements in table 1. The probabilities from the quantum computation are largely consistent with the pure state results, with the majority of results within one s.e.m. from the ideal , as expected in finite sampling to estimate the ground state probability. There are only two large deviations for states ‘k’ and ‘i’, which may be related to noise in the device. The quantum computations succeed in preparing ground states with probabilities comparable to pure-state expectations.
5. Conclusion
We analysed QAOA as an approach for preparing materials ground states on three types of Ising Hamiltonians with longitudinal magnetic fields, focusing on nine-spin unit cells as a starting size that is amenable to simulations and calculations on a trapped-ion quantum computer.
We applied QAOA to the nine-spin Ising unit cell problems to assess its success in ground state preparation. We found that the theoretical success probability depends significantly on the structure of the ground state, while it is mostly insensitive to the precise Hamiltonian parameters, which can vary within regions consistent with a fixed ground state. Each Hamiltonian yields a ferromagnetic ground state in the presence of large magnetic fields, and QAOA achieved large success probabilities for these relatively simple states. The probabilities for other types of ground states were more variable, and tended to decrease as next-nearest-neighbour couplings became stronger with associated frustration in the lattice. For all of these nine-spin states, we typically find success probabilities indicating that are expected to be necessary from an ideal quantum computer for these problem instances.
To assess QAOA performance under realistic conditions, we implemented the algorithm on a trapped-ion quantum computer. These quantum computations succeeded in observing each of the 17 ground states of the Shastry-Sutherland unit cell. The quantum computations yielded ground state probabilities that were consistent with theoretical expectations based on pure states, indicating that noise was not a significant issue at the sizes and depths tested. This suggests that calculations with current technology can likely be extended to greater QAOA depth parameters , and to larger sizes as greater numbers of qubits become available. At greater depths and sizes we expect higher performance and more realistic results in comparison with the large size limit, respectively.
While the ground states and associated phase diagrams for our nine-spin unit cells were found to have significant finite size effects relative to the large-size limit, the errors from finite size effects on the classically calculated magnetization phase diagrams on lattices up to was found to be suppressed rapidly with , with small errors at indicating that only hundreds of spins may be necessary to reproduce large scale behaviour. This provides a baseline of hundreds of qubits for quantum computational experiments that seek to explain materials science problems, which may be accessible to near-term quantum computers in coming years.
Assessing scaling of the ground-state probability with size will be an essential aspect of extending this approach to larger sizes . This includes numerical simulations to quantify how the ground-state probabilities depend on the number of spins and number of QAOA layers ; previous works have shown maintains a large ground-state probability for simple models in different contexts [37,38], but future work is needed to test scaling in the current model. Benchmarking on quantum computers is also essential to understand how real noise processes effect scalability.
Thus from our results we envision QAOA can be successfully applied to somewhat larger lattice problems as quantum computing technologies develop and larger number of qubits become available. These could be used for optimization of Ising lattice problems as we have here, with increasing sizes that potentially describe real bulk properties of materials in the limit. Additionally, the QAOA algorithm is general, with the application of a magnetic field, and hence could be explored for lattice problems which are NP-Hard [79]. But a more promising future direction leveraging the full benefit of this approach is to extend and modify QAOA to prepare ground states of quantum Hamiltonians such as the and Heisenberg models, which can lead to a variety of quantum phenomena not captured in the Ising model, such as quantum spin glasses [80], spin nematicity [81], Berzinski–Kosterlitz–Thouless states [82,83] and long-range entangled states such as Dirac string excitations [84], the likes of which exist in two-dimensional frustrated quantum spin liquids and spin ice. Many of these topics are fiercely researched and are of considerable interest and importance for future quantum technologies and devices. Conventional numerical methods for understanding these states are hindered by the exponential size of the Hilbert space, making it difficult to generate a theoretical understanding of experimental observations. QAOA or related generalizations [38,61,67–70] offer a potential route to overcome conventional computing bottlenecks. Some successes along these lines have been observed in certain contexts, however, advances in methodology and quantum computing technologies are needed to extend these methods to complicated and larger-scale problems where quantum computational approaches may have a significant impact in understanding and developing materials for technological applications.
Data accessibility
Data and code from this study are available online at https://code.ornl.gov/5ci/dataset-simulations-of-frustrated-ising-hamiltonians-using-qaoa.
The data are provided in electronic supplementary material [85].
Authors' contributions
P.C.L.: conceptualization, formal analysis, investigation, methodology, writing—original draft, writing—review and editing; H.X.: formal analysis, investigation, writing—original draft, writing—review and editing; B.K.: formal analysis, writing—original draft, writing—review and editing; G.B.: formal analysis, supervision, writing—original draft, writing—review and editing; T.S.H.: conceptualization, formal analysis, project administration, supervision, writing—original draft, writing—review and editing; A.B.: conceptualization, formal analysis, project administration, supervision, writing—original draft, writing—review and editing.
All authors gave final approval for publication and agreed to be held accountable for the work performed therein.
Conflict of interest declaration
The authors declare that they have no competing interests.
Funding
This material is based upon work supported by the US Department of Energy, Office of Science, National Quantum Information Science Research Centers, Quantum Science Center.
Acknowledgements
The authors thank Paul Kairys for interesting discussions. This research used resources of the Oak Ridge Leadership Computing Facility at the Oak Ridge National Laboratory, which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC05-00OR22725. P.C.L. and T.S.H. acknowledge support from the Office of Science, Early Career Research Program. A.B., G.B. and B.K. acknowledge funding from the US Department of Energy, Office of Science, National Quantum Information Science Research Centers, Quantum Science Center.
Disclaimer
This manuscript has been authored by UT-Battelle, LLC under Contract No. DE-AC05-00OR22725 with the U.S. Department of Energy. The United States Government retains and the publisher, by accepting the article for publication, acknowledges that the United States Government retains a non-exclusive, paid-up, irrevocable, world-wide license to publish or reproduce the published form of this manuscript, or allow others to do so, for United States Government purposes. The Department of Energy will provide public access to these results of federally sponsored research in accordance with the DOE Public Access Plan. (http://energy.gov/downloads/doe-public-access-plan).
