Simulating the Shastry-Sutherland Ising Model Using Quantum Annealing

A core concept in condensed matter physics is geometric frustration that leads to emergent spin phases in magnetic materials. These distinct phases, which depart from the conventional ferromagnet or the antiferromagnet, require unique computational techniques to decipher. In this study, we use the canonical Ising Shastry-Sutherland lattice to demonstrate new techniques for solving frustrated Hamiltonians using a quantum annealer of programmable superconducting qubits. This Hamiltonian can be tuned to produce a variety of intriguing ground states ranging from short- and long-range orders and fractional order parameters. We show that a large-scale finite-field quantum annealing experiment is possible on 468 logical spins of this model embedded into the quantum hardware. We determine microscopic spin configurations using an iterative quantum annealing protocol and develop mean-field boundary conditions to attenuate finite-size effects and defects. We not only recover all phases of the Shastry-Sutherland Ising model—including the well-known fractional magnetization plateau in a longitudinal field—but also predict the spin behavior at the critical points with significant ground-state degeneracy and in the presence of defects. The results lead us to establish the connection to the diffuse neutron scattering experiments by calculation of the static structure factors.

Popular Summary

Quantum computing is a rapidly developing discipline, promising to improve our understanding of the microscopic world, as well as impact the wider fields of science and technology. One such promising application is the use of quantum computers to simulate the behavior of materials. However, demonstrating this application is challenging due to the limited size and reliability of current quantum computing hardware. The approach outlined here takes a step toward this goal by demonstrating an example of quantum-enabled calculations that can be validated via experiments on real materials.

In this manuscript, the use of a quantum annealer is outlined, with a scope to model and simulate a class of real-world magnetic materials. In particular, materials known to demonstrate exotic microscopic magnetic structures that emerge when placing these materials in an external magnetic field are considered. Using several novel techniques, quantities that can be probed experimentally in real-world materials using neutron scattering are analyzed. This offers a tantalizing glimpse into the potential for quantum simulation to support the understanding of experimental measurements.

The results from this work open new avenues for exploring the novel physics of magnetism using quantum annealers and support future efforts to use quantum simulation to improve our understanding of materials. These rapidly developing applications of quantum simulations offer new insights into the ways in which quantum computing may impact the broader interdisciplinary fields.

Figure 1

Embedding the Shastry-Sutherland lattice in the D-Wave quantum annealer. (a) The Shastry-Sutherland lattice showing the nearest-neighbor square (solid black lines) and the next-nearest-neighbor dimer (dashed black lines) interactions. Highlighted in gray is the logical unit cell with eight spins labeled. (b) The embedding of the logical unit cell from (a) into a Chimera graph (full embedded graph shown in Fig. 5). This “half-cell” embedding maps each logical spin to a chain of four ferromagnetically connected qubits (solid color lines) and each dimer pair to a single Chimera unit cell. Each logical square bond is mapped to one external coupler (solid black lines) between Chimera unit cells and each dimer bond is mapped to eight internal couplers (dashed black lines) within a Chimera unit cell. This embedding preserves the symmetry of the logical problem and ensures that every physical qubit experiences similar local environments, which prevents chain breaking. (c) The embedding of a 468-spin lattice in the quantum annealer. Missing lattice sites in the underlying Chimera graph reveal defects that originate from defective qubits or bonds in the annealer. Effects from these defects, as well as from finite boundaries, are mitigated by applying self-consistent mean-field boundary conditions (see the text). These boundary conditions are applied on both the external boundaries and the internal boundaries (qubits surrounding defects) simultaneously.

Figure 2

Optimal QEMC annealing parameters. (a) Simulated energy distributions using QEMC chains at $(J_1=1,J_2=1,h_z=2.1)$ with $t_r=1\mu \mathrm{s}$ ramps and $t_p=1998\mu \mathrm{s}$. The energies displayed are rescaled according to the lowest-energy configuration observed in experiment. Two values, $s_p=0.7,s_p=0.8$, are omitted from (a) for clarity because $⟨H(s_p=0.7)⟩\ge 600$ and $⟨H(s_p=0.8)⟩\ge 700$ lie far outside of the domain shown in (a). The inset is an alternative view of (a) that clearly shows the parabolic energy curve as a function of $s_p$ with a minimum at $s_p=0.4$. (b) The convergence of the magnetization from a ferromagnetic state at $(J_1=1,J_2=1,h_z=2.1)$ with $t_r=1\mu \mathrm{s}$ ramps and $t_p=1998\mu \mathrm{s}$ to the $\frac{1}{3}$ magnetization plateau (shown in black). For $s_p<0.5$, convergence to the bulk magnetization occurs in less than $5$ iterations; however, these states are not energetically optimal and convergence to the correct microscopic order takes many tens of iterations. (c) The average returned energy as a function of $s_p$ and $t_p$ with $t_r=1\mu \mathrm{s}$. Unless $s_p$ is sufficiently strong such that the timescale of tunneling dynamics is much shorter than the anneal time, increasing the pause time has a dramatic effect in finding low-energy solutions.

Figure 3

Phase diagram and spin structures. (a) Phase diagram obtained by QEMC chains with $t_r=1\mu \mathrm{s}$, $t_p=1998\mu \mathrm{s}$, and $s=0.4$. The simulated phase boundaries, shown with red crosses, are determined by calculating the point of steepest slope of magnetization. This helps to distinguish the $\frac{1}{3}$ plateau, the dimer, FM, and the Néel states. The white lines outline the exact phase diagram as determined analytically by Dublenych [6]. A slice of the phase diagram at $J_2/J_1=1$ is shown in (b). (b) Slice of phase diagram at ($J_1=1,J_2=1$) obtained by QEMC chains with $t_r=1\mu \mathrm{s}$, $t_p=1998\mu \mathrm{s}$, and $s=0.4$. The observed phase transition is broadened compared to exact results due to a number of compounding effects, the most notable being bond disorder of $J_1,J_2$, residual effects of point defects and finite lattice size, non-negligible persistent transverse field, and temperature effects. Each data point is the average of 100 QEMC chains, each of 100 iterations in which we sample the last 50, providing 5000 spin configurations per data point. We observe that the standard deviation of the magnetization, ${\sigma }_{⟨m⟩}$, lies within the range $0.002\le {\sigma }_{⟨m⟩}<0.03$ and has characteristic peaks within the critical regimes. We find that ${\sigma }_{⟨m⟩}$ is sufficiently small that the round markers obscure the plotted error bars completely and so it is plotted above the magnetization for clarity. Four points are labeled (c)–(f) that correspond to the structures shown in panels (c)–(f). The critical regimes are highlighted in white. (c),(d) The real-space spin motifs determined by QEMC chains in the Néel AFM and critical AFM to plateau phase calculated at the points labeled in (b). In (d) we observe a coexistence between the Néel and plateau orderings separated by domain walls. (e),(f) Two degenerate solutions within the plateau phase. The black and red squares represent spin up and down, respectively.

Figure 4

Static structure factors in equilibrium and at criticality. (a) Structure factor of the Néel AFM phase at $J_1=J_2=1,h_z=0$. The high symmetry points of this model are labeled in white along with a cut path used in (e); these points are $(q_x=0,q_y=0)$, $(q_x=\pi ,q_y=0)$, and $(q_x=\pi ,q_y=\pi )$. (b) Structure factor of the dimer AFM phase at $J_1=1,J_2=3,h_z=0$. (c) Structure factor of the $\frac{1}{3}$ plateau phase at $J_1=J_2=1,h_z=2.1$. (d) Structure factor in the transition between Néel AFM and $\frac{1}{3}$ plateau phases at $J_1=J_2=1,h_z=1.8$. (e) A cut along the symmetry points shown in (a) of the structure factor through the Néel to plateau transition. The inset figure shows the color scale corresponding to the longitudinal field across the phase transition at $J_2/J_1=1$.

Figure 5

Half-cell embedding into the D-Wave quantum annealer with unit cell as depicted in Fig. 1. This is rotated ${45}^{\circ }$ to obtain the lattice shown in Fig. 1.

Figure 6

(a) The annealing parameter $s$ controls amplitudes $A(s)$ and $B(s)$ in the device Hamiltonian, Eq. (2). Here $Q(s)$ quantifies the relative strengths of these terms as a function of $s$. Also shown is the dependence on $s$ of the ratio to amplitude $B(s)$ and the characteristic thermal energy of the system at the operating temperature of 14 mK. (b) A single QEMC anneal schedule that shows the initialization of the annealer with an eigenstate of the classical Hamiltonian [$A(s)\ll B(s)$] at $t=0$, reverse annealing for time $t_r$, pausing for time $t_p$ at $s_p$ [where $Q(s)\ge 1$], and forward annealing for time $t_r$. (c) An iterative protocol where multiple QEMC steps are chained together by reprogramming the next QEMC step with the output from the previous QEMC step. The first half of the samples are discarded as burn-in steps while the chain equilibrates and the last half of the samples are collected and used for analysis.