## Abstract

A digital quantum simulator is an envisioned quantum device that can be programmed to efficiently simulate any other local system. We demonstrate and investigate the digital approach to quantum simulation in a system of trapped ions. With sequences of up to 100 gates and 6 qubits, the full time dynamics of a range of spin systems are digitally simulated. Interactions beyond those naturally present in our simulator are accurately reproduced, and quantitative bounds are provided for the overall simulation quality. Our results demonstrate the key principles of digital quantum simulation and provide evidence that the level of control required for a full-scale device is within reach.

Although many natural phenomena are accurately described by the laws of quantum mechanics, solving the associated equations to calculate properties of physical systems, i.e., simulating quantum physics, is in general thought to be very difficult (*1*). Both the number of parameters and differential equations that describe a quantum state and its dynamics grow exponentially with the number of particles involved. One proposed solution is to build a highly controllable quantum system that can efficiently perform the simulations (*2*). Recently, quantum simulations have been performed in several different systems (*3*–*13*), largely following the analog approach (*2*) whereby an analogous model is built, with a direct mapping between the state and dynamics of the simulated system and those of the simulator. An analog simulator is dedicated to a particular problem, or class of problems.

A digital quantum simulator (*2*, *14*–*16*) is a precisely controllable many-body quantum system on which a universal set of quantum operations (gates) can be performed, i.e., a quantum computer (*17*). The simulated state is encoded in a register of quantum information carriers, and the dynamics are approximated with a stroboscopic sequence of quantum gates. In principle, it can be reprogrammed to efficiently simulate any local quantum system (*14*) and is therefore referred to as a universal quantum simulator. Furthermore, there are known methods to efficiently correct for and quantitatively bound experimental error in large-scale digital simulations (*18*).

We report on digital simulations using a system of trapped ions. We focus on simulating the full time evolution of networks of interacting spin-1/2 particles, which are models of magnetism (*19*) and exhibit rich dynamics. We do not use error correction, which has been demonstrated separately in our system (*20*) and must be included in a full-scale fault-tolerant digital quantum simulator.

The central goal of a quantum simulation is to calculate the time-evolved state of a quantum system ψ(*t*). In the case of a time-independent Hamiltonian *H*, the form of the solution is ψ(*t*) = *e*^{−}^{iHt}^{/}* ^{ħ}*ψ(0) =

*U*ψ(0). A digital quantum simulator can solve this equation efficiently for any local quantum system (

*14*), i.e., where

*H*contains a sum of terms

*H*that operate on a finite number of particles, owing to interaction strengths that fall off with distance, for example. In this case the local evolution operators

_{k}*U*=

_{k}*e*

^{−}

^{iHkt}^{/}

*can be approximated with a fixed number of operations from a universal set. However, these terms do not generally commute*

^{ħ}*U*≠ ∏

_{k}e^{−}

^{iHkt}^{/}

*. This can be overcome with the Trotter approximation (*

^{ħ}*14*,

*21*),

*n*, which is at the heart of the digital quantum simulation algorithm. For finite

*n*, the Trotter error is bounded and can be made arbitrarily small. The global evolution of a quantum system can therefore be approximated by a stroboscopic sequence of many small time-steps of evolution, generated by the local interactions present in the system. The digital algorithm can also be applied to time-dependent Hamiltonians and open quantum systems (

*14*,

*16*,

*17*,

*22*).

Our simulator is based on a string of electrically trapped and laser-cooled calcium ions (*23*). The |*S*_{1/2}〉 = |1〉 and |*D*_{5/2}〉 = |0〉 Zeeman states encode a qubit in each ion. Simulated states are encoded in these qubits and manipulated by laser pulses that implement the operation set: *O*_{1}(θ,*j*) = exp(−*i*θσ* ^{j}_{z}*);

*O*

_{2}(θ) = exp(−

*i*θ∑

*σ*

_{i}*);*

_{z}^{i}*O*

_{3}(θ,ϕ) = exp(−

*i*θ∑

*σ*

_{i}_{ϕ}

*);*

^{i}*O*

_{4}(θ,ϕ) = exp(−

*i*θ∑

_{i}_{<}

*σ*

_{j}

^{i}_{ϕ}σ

^{j}_{ϕ}). Here σ

_{ϕ}= cosϕσ

*+ sinϕσ*

_{x}*and σ*

_{y}*denotes the*

^{j}_{k}*k*-th Pauli matrix acting on the

*j*-th qubit.

*O*

_{4}is an effective qubit-qubit interaction mediated by a common vibrational mode of the ion string (

*24*). Recent advances have seen the quality of these operations increase appreciably (

*25*). For our simulations, we define dimensionless Hamiltonians

*Et*/

*ħ*.

We begin by simulating an Ising system of two interacting spin-1/2 particles, which is an elementary building block of larger and more complex spin models and was one of the first systems to be simulated with trapped ions following an analog approach (*6*, *26*). The Hamiltonian is given by *z* direction and the second an interaction between the spins in an orthogonal direction. The interactions do not commute, giving rise to nontrivial dynamics and entangled eigenstates. Each spin is mapped directly to an ionic qubit (|1〉 = |↑〉, |0〉 = |↓〉). The dynamics are implemented with a stroboscopic sequence of *O*_{2} and *O*_{4} gates, representing the magnetic field and spin-spin evolution operators, respectively. We first simulate a time-independent case *J* = 2*B*, which couples the initial state |↑↑〉 to a maximally entangled superposition of |↑↑〉 and |↓↓〉 (Fig. 1A). The simulated dynamics converge closer to the exact dynamics as the digital resolution is increased. The overall simulation quality is quantified by quantum process tomography (QPT) (*27*), yielding a process fidelity of 91(1)% at the finest digital resolution used. In (*23*), we show how higher-order Trotter decompositions can be used to achieve more accurate digital approximations with fewer operations.

We now consider a time-dependent case where *J* increases linearly from 0 to 4*B* during a total evolution θ* _{t}*. In the limit θ

*→∞, spins initially prepared in the paramagnetic ground state of the magnetic field (|↑↑〉) will evolve adiabatically into the antiferromagnetic ground state of the final Hamiltonian: an entangled superposition of the ∑*

_{t}*σ*

_{j}*eigenstates |←→〉*

_{x}^{j}*and |→←〉*

_{x}*. As a demonstration, we approximate the continuous dynamics, for θ*

_{x}*= π/2, using a stroboscopic sequence of 24*

_{t}*O*

_{2}and

*O*

_{4}gates and measure the populations in the σ

*basis (Fig. 1B). The evolution closely follows the exact case, and an entangled state is created [63(6)% tangle (*

_{x}*28*)]. Full quantum state reconstructions are performed after each digital step, yielding fidelities between the ideal digitized and measured state of at least 91(2)% and overlaps with the instantaneous ground state of no less than 91(2)%. The observed oscillation in expectation values is a diabatic effect, as excited states become populated.

More complex systems with additional spin-spin interactions in the *y* (“XY” model) and *z* (“XYZ” model) directions can be simulated by reprogramming the operation sequence. The dynamics due to an additional spin-spin interaction in the *y* direction is simulated by adding another *O*_{4} operation to each step of the Ising stroboscopic sequence (with ϕ = π/2). A third spin-spin interaction in the *z* direction is realized by adding an *O*_{4} gate sandwiched between a pair of *O*_{3} operations set to rotate the reference frame of the qubits. Dynamics of the initial state |→←〉* _{x}* are simulated for each model, with a fixed digital resolution of θ/

*n*= π/16 and up to 12 Trotter steps (Fig. 2). Up to 24, 48, and 84 gates are used for the Ising, XY, and XYZ simulations, respectively. This particular initial state is chosen because the ideal evolution is different for each model. The results show close agreement with the exact dynamics and results from QPT after four digital steps yield process fidelities, with the exact unitary evolution, of 88(1), 85(1), and 79(1)% for the Ising, XY, and XYZ, respectively. With perfect operations, the Trotter error would be less than 1% in each case. Although analog simulations of Ising models have previously been demonstrated in ion traps (

*6*,

*8*), XY and XYZ models have not.

The digital approach allows arbitrary interaction distributions between spins to be programmed. For three-spin systems, we realize various interactions that give rise to the dynamical evolutions of the initial state |↑↑↑〉 (Fig. 3). Figure 3A shows a system supporting interactions between all spin pairs with equal strength, and between each spin and a transverse field. The initial state couples equally to |↑↓↓〉, |↓↑↓〉, and |↓↓↑〉, while the strength of the field determines the amplitude and frequency of the dynamics. For the case shown (*J* = 2*B*), an equal superposition of the coupled states is periodically created [an entangled W state (*29*)]. Figure 3B shows how nonsymmetric interaction distributions can be programmed, with sequences of *O*_{4} and *O*_{1} to add spin-selective interactions. The interaction between one spin pair is dominant. Owing to this broken symmetry, one coupled state (|↑↓↓〉) is populated faster than the others, yielding more complex dynamics than in the symmetric case. Figure 3C demonstrates the ability to simulate *n*-body interactions; specifically, σ_{z}^{1}σ_{x}^{2}σ_{x}^{3}. A clear signature is observed: direct coupling between ∑* _{j}*σ

*eigenstates |→→→〉*

_{y}^{j}*and |←←←〉*

_{y}*, periodically producing an entangled GHZ state (*

_{y}*29*). Many-body spin interactions of this kind are an important ingredient in the simulation of systems with strict symmetry requirements (

*30*) or spin models exhibiting topological order (

*31*). Measurements in other bases and simulations of nearest-neighbor and many-body interactions with a transverse field using more than 100 gates are presented in (

*23*).

Figure 4A shows the observed dynamics of the four-spin state |↑↑↑↑〉 under a long-range Ising-type interaction. The rich structure of the dynamics reflects the increased complexity of the underlying Hamiltonian: Oscillation frequencies correspond to the energy gaps in the spectrum. This information can be extracted via a Fourier transform of the data (*23*). Specific energy gaps could be targeted by preparing superpositions of eigenstates via an initial quasi-adiabatic digital evolution (*10*). Figure 4B shows the observed dynamics for our largest simulation: a six-spin many-body interaction, which directly couples the states |↑↑↑↑↑↑〉 and |↓↓↓↓↓↓〉, periodically producing a maximally entangled GHZ state.

Direct quantification of simulation quality for more than two qubits is impractical via QPT: For three qubits, expectation values must be measured for 1728 experimental configurations, and this increases exponentially with qubit number (≈3 × 10^{6} for six qubits). However, the average process fidelity (*F*_{p}) can be bounded more efficiently (*32*). We demonstrate this for the three- and six-spin simulations of Figs. 3C and 4B, respectively. Bounds of 85(1)% ≤ *F*_{p} ≤ 91(1)% (three spins) and 56(1)% ≤ *F*_{p} ≤ 77(1)% (six spins) are obtained at ϕ = 0.25 π, with 40 and 512 experimental configurations, respectively (*23*). The largest system for which a process fidelity has previously been measured is three qubits (*33*). A different measure of process quality is given by the worst-case fidelity, over all input states, and may be better used to assess errors in future full-scale fault-tolerant simulations. Regardless of the measure used, the error in large-scale digital simulations built from finite-sized operations can be efficiently estimated. Each operation can be characterized with a finite number of measurements, then the error in any combination can be chained (*34*). To exploit this, the number of ionic qubits on which our many-qubit operators *O*_{2−4} can act must be restricted.

The dominant effect of experimental error can be seen in Figs. 3B and 4B: The dynamics damps as a result of decoherence processes. Laser frequency and ambient magnetic field fluctuations are far from the leading error source: In the absence of coherent operations, qubit lifetimes are more than an order of magnitude longer (coherence times ≈30 ms) than the duration of experiments (≈1 to 2 ms). The current leading sources of error, which limit both the simulation complexity and size, are thought to be laser intensity fluctuations (*23*). This is not currently a fundamental limitation and, once properly addressed, should enable an increase in simulation capabilities.

The digital approach can be combined with existing tools and techniques for analog simulations to expand the range of systems that can be simulated. In light of the present work, and current ion trap development (*35*), digital quantum simulations involving many tens of qubits and hundreds of high-fidelity gates seems feasible in coming years.

## Supporting Online Material

www.sciencemag.org/cgi/content/full/science.1208001/DC1

SOM Text

Figs. S1 to S9

Tables S1 to S6

References

## References and Notes

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
Materials and methods are available as supporting material on
*Science*Online. - ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
**Acknowledgments:**We thank W. Dür, A. Aspuru-Guzik, and M. Brownnutt for discussions. We acknowledge financial support by the Austrian Science Fund (FWF) [SFB F40 FOQUS], the Institut für Quanteninformation GmbH, Intelligence Advanced Research Projects Activity, and the European Commission for support via the integrated project AQUTE, two Marie Curie International Incoming Fellowships, and the ERC advanced grant CRYTERION.