Quantum Phases of Trapped Ions in an Optical Lattice

We propose loading trapped ions into microtraps formed by an optical lattice. For harmonic microtraps, the Coulomb coupling of the spatial motions of neighboring ions can be used to construct a broad class of effective short-range Hamiltonians acting on an internal degree of freedom of the ions. For large anharmonicities, on the other hand, the spatial motion of the ions itself represents a spin-1/2 model with frustrated dipolar XY interactions. We illustrate the latter setup with three systems: the linear chain, the zig-zag ladder, and the triangular lattice. In the frustrated zig-zag ladder with dipolar interactions we find chiral ordering beyond what was predicted previously for a next-nearest-neighbor model. In the frustrated anisotropic triangular lattice with nearest-neighbor interactions we find that the transition from the one-dimensional gapless spin-liquid phase to the two-dimensional spiraling ordered phase passes through a gapped spin-liquid phase, similar to what has been predicted for the same model with Heisenberg interactions. Further, a second gapped spin-liquid phase marks the transition to the two-dimensional Neel-ordered phase.


Introduction
The merging of atomic physics and condensed matter physics has opened exciting new perspectives for the creation and manipulation of quantum states of matter. For example, cold atoms in optical lattices are a setup which allows experimentalists to simulate quantum many-body lattice systems [1,2]. In recent years, trapped ions have also been proposed as an experimental system where a rich variety of quantum many-body models can be implemented [3,4,5,6,7]. The advantage of trapped ions is that we profit from the technology that has been developed for quantum information processing [8,9,10,11,12,13]. In particular, it is possible to measure and manipulate the system at the single-particle level, due to the large inter-particle distances and the optical accessibility of many electromagnetic traps.
So far, previous proposals fall into two categories: (i) effective-spin models, in which the motion of the ions induces an effective interaction between internal levels [4,5,14], and (ii) interacting-Boson models, in which the vibrations (phonons) play the role of interacting particles [6,7]. In this work we show that these proposals may be scaled up to a wide range of spatial dimensions, geometries, and particle interactions by placing the ions in a regular structure. This structure may be induced by optical lattices or by an array of ion microtraps.
The presence of the trapping structure offers us a way to increase the dimensionality of the system beyond the traditional linear Paul trap scheme. Further, it also opens up new possibilities for controlling particle interactions. On one hand, since ions are tightly confined by individual trapping potentials, their motion in any spatial direction can be used to induce spin-spin interactions which decay with the third power of the distance. This situation is to be compared with the case of Coulomb crystals confined in an overall trapping potential, where vibrational phonons with small wavevector, corresponding to highly collective vibrations, have a low-energy spectrum and induce infinite-range spin-spin interactions only [5]. On the other hand, in the optical lattice trapping scheme, we can control the anharmonicities of the trapping potentials and engineer effective phonon-phonon interactions. Thus, vibrational degrees of freedom may follow a Bose-Hubbard model along the lines of [6].
The realization of a periodic trapping potential with an optical lattice is a standard experimental method. In this work we focus on the study of the quantum many-body phases which arise in such experimental setups. We discuss two limiting cases. First, the case of small anharmonicities, and the effective quantum spin models which can be implemented by using the ion motion as a medium for inducing effective spin interactions among internal degrees of freedom. Second, we discuss in more detail the ground vibrational state in trapping schemes in which large anharmonicities induce phonon-phonon interactions and allow mapping hard-core phonons onto S = 1 2 spins. We focus on the case of periodic trapping potentials with interesting geometries, e.g., the triangular lattice. This case provides us with a playground for studying frustrated quantum magnets in a clean experimental setup. In this way, trapped ions allow us to implement models which have attracted much attention in condensed matter physics.

Microtraps
The starting point for all of the models in this work is a set of ions trapped in individual microtraps. Such traps have been proposed in the literature to be microfabricated, e.g., on surfaces [15], which allows for arbitrary geometries of the resulting ion "crystal" and, in principle, gives the experimentalist great control over the effective anharmonicity of the spatial ion motion near the microtrap minima.
We propose an alternative scheme for creating microtrap potentials, postponing detailed calculations to a further publication [16]. Optical lattices have been used extensively in experiments with cold atoms, which interact mostly on short length scales and can populate d-dimensional lattices densely. If a d-dimensional Wigner crystal, assembled in an electromagnetic trap, is subjected to a sufficiently strong optical lattice potential, the spatial motion of each ion is restricted to a region near the minimum of the lattice well closest to the crystal equilibrium position. However, the effective amplitude of a realistic optical lattice is much smaller than the repulsion of two ions in the same or even neighboring lattice wells. Therefore, as opposed to experiments with cold atoms, the optical lattice will be very sparsely populated with ions, but their strong interactions over large distances are nevertheless sufficient for experimentally exploring the models presented in section 4. Depending on the stiffness of the Wigner crystal used as a starting point, the optical lattice does not need to be so strong as to trap the ions by itself; it is sufficient that the zero-point motion of the Wigner crystal is reduced by the optical lattice (since the optical lattice increases the effective mass of the ions) to the point where each ion is constrained to moving within a single lattice minimum. In a one-dimensional optical lattice potential V (z) = αE R sin 2 (qz), with the recoil energy E R = 2 q 2 /(2m), the harmonic frequency of the ions moving in their local microtraps is ω mt = 2 √ αE R . Recoil energies are typically tens to hundreds of kHz; with an optical lattice depth of hundreds of recoil energies (α 100) it is thus feasible to reach microtrap frequencies ω mt far in excess of the ion trap frequency Ω, the latter being typically on the same order as E R but possibly as low as necessary. Lowering the ion trap frequency significantly below the recoil energy of the optical lattice is an alternative to large values of α for achieving well-defined microtraps. The energy scale of the couplings between the vibrational motion of different ions, separated by a typical distance d 0 , is E t = E C (qd 0 ) −2 α −1/2 , with the Coulomb energy scale E C = e 2 /(4πε 0 d 0 ) [see also (2)]. The low filling fraction of the optical lattice (qd 0 ≫ 1) reduces the very large Coulomb energy scale, typically tens or hundreds of GHz, to below the microtrap energy scale ω mt : in a harmonic trap, the typical distance between ions is of order Keeping this coupling energy scale E t below the microtrap energy scale is necessary in order to prevent processes that change the number of vibrational excitations, which are assumed to be absent throughout the rest of this work.
As we have already mentioned, many interesting models require significant phonon-phonon interactions (see section 3.3). The anharmonicity of the Coulomb interaction in a Wigner crystal is insufficient for observing these phenomena, and a harmonic normal-mode analysis captures the full physics of low-energy ion motion. The reason for this is that the length scale ζ of ion vibrations is only a fraction of the optical-lattice wavelength and thus much smaller than the typical distance d 0 between neighboring ions; quartic terms of the Coulomb interaction, which give rise to phonon-phonon interactions, are smaller than the harmonic terms (energy scale E t ) by a factor (ζ/d 0 ) 2 ≪ 1 [17]. However, we can make use of the natural anharmonicity of the sinusoidal optical lattice potential to introduce significant phonon-phonon interactions. This anharmonicity is close to the recoil energy E R for α 10. The ratio of the anharmonicity to the coupling energy scale is E R /E t = ω mt /E t /(2 √ α); with the above assumption of E t ≪ ω mt , this ratio is typically larger than one, giving the effective hard-core repulsion between phonons which is the basis of the models presented in section 4. But it must be noted that the sign of the anharmonicity is negative for a strong monochromatic optical lattice. If the lattice is sufficiently strong, the ions are forced to the lattice minima, where the quartic anharmonicity is negative and leads to an effective phonon-phonon attraction; effective phonon models with attractive interactions have very simple ground states and do not merit experimental study on a quantum simulator. Much more interesting phases result from repulsive phonon-phonon interactions, which can be achieved by using an optical superlattice [18]. Briefly, if a second optical lattice with half the wavelength of the first one is superposed with the same phase and a reduced amplitude α/64 < α ′ < α/16, then the anharmonicity at the superlattice minima can be made positive while still retaining the same spatial arrangement of microtrap minima; further, the amplitude ratio α ′ /α is a sensitive parameter for adjusting the ratio between the energy scales of phonon tunneling (E t ) and repulsion (E R ). In what follows we assume repulsive phonon-phonon interactions, irrespective of how the microtraps are put in place.

Varying the Anharmonicity
This section describes three regimes of anharmonicities of the microtraps. If the anharmonicities are negligible, the internal degrees of freedom of trapped ions can be made to interact via dipolar phonon-mediated couplings. At the opposite end, for strong anharmonicities the quantized vibrational modes themselves describe a S = 1 2 spin model with XY interactions. In the intermediate regime, the vibrational motion of the ions can be effectively described by an interacting Bose-Einstein condensate.

No anharmonicity: short-range lattice Hamiltonians
In [5] a model is elaborated where internal states of trapped ions represent effective S = 1 2 spin degrees of freedom. In this model, the vibrational modes of a self-assembled ion crystal in a harmonic trap mediate variable interactions between the effective spins, allowing for controllable interactions in all spin directions. A drawback of the proposed architecture is that the longitudinal vibrational modes, where the ion crystal moves as a whole in the ion trap, are naturally slow ("floppy" modes) and lead to long-range ferromagnetic interactions. This limits the generality of the Hamiltonians that can be constructed with such crystals of trapped ions.
If, however, local microtraps are added to the confining potential, then this problem can be overcome by increasing the potential energy curvature in all directions of motion of the individual trapped ions. In the simplest case, adding the same isotropic harmonic microtrap potential to each ion, the dispersion and shape of normal modes of the ion crystal are unchanged except for an overall shift in frequency corresponding to the curvature of the microtraps. Such a frequency shift renders the vibrational modes "stiff" for all three spatial directions of each ion. Consequently the effective spin-1 2 degrees of freedom of the trapped ions can be made to interact through Ising, XY, or XYZ couplings, all of which have adjustable amplitudes and decay as r −3 since they stem from the dipole-dipole interactions of the ionic vibrations. [5] further shows that if the phonon modes have only finite stiffness, then next-nearest-neighbor interactions can be suppressed, approaching effective nearest-neighbor Hamiltonians on arbitrary lattice geometries.

Small anharmonicity: interacting Bose-Einstein condensate
As discussed in section 2, the vibrational motion of a self-assembled ion crystal in a global trap is very well approximated by a harmonic model, as was used in Sec. 3.1. At low energies, cubic terms of the Coulomb interaction are unimportant because the three-phonon processes they describe are highly nonresonant; quartic terms, on the other hand, can be neglected because of their small amplitudes.
This picture changes when microtraps are added to the local potentials of the individual ions. In addition to the harmonic confinement, these microtraps can add significant anharmonicities U α to the ionic vibrations. As discussed in [6,7], in the absence of fast phonon decay channels the spatial motion of trapped ions thus behaves like a Bose-Hubbard model of phonons, with Hamiltonian where the indices α, β run over the 3N eigen-directions of motion of the N ions.
Here we assume that the quartic phonon-phonon interaction terms (anharmonicities) do not couple different modes, even if the modes are localized on the same ion. This decoupling occurs naturally in a cubic optical lattice; but, irrespective of the microtrap geometry, we will assume throughout this work that of the three directions of motion of each ion, only one is "active" and contributing to the models we wish to study, while the other two modes are off-resonant and do not influence the lowenergy dynamics of the active mode. This is achieved with strongly anisotropic prolate microtraps, where one slow direction of motion is decoupled from the two fast directions of motion. (Extensions of the present models to two or three ionic degrees of freedom are straightforward and do not introduce new physics if the degrees of freedom are separable.) For example, in a three-dimensional cubic optical lattice where one standing wave is much weaker than the other two, every lattice minimum constitutes such a prolate microtrap.
Calling m α the direction of motion of the α th ion, the effective tunneling matrix elements of the phonons are of the dipolar form where e 2 /(4πε 0 ) is the electrostatic coupling strength, r α,β = r are the equilibrium inter-particle spacings, and n α,β = r α,β /r α,β . The length scales ζ α are those of the harmonic-oscillator motion in the quadratic component of the local microtrap potentials, related to the local potentials as ζ α = / √ mV α , where m is the mass of the ions. In a homogeneous lattice of microtraps, where all m α are equal, the angular dependence of (2) is proportional to the Legendre polynomial t α,β ∝ −P 2 ( n α,β · m). In the following sections, we will focus our attention on the central region of the trap, where the tunneling coefficients t α,β and local potentials V α are sufficiently isotropic to be approximated by a translationally invariant model.
In the absence of anharmonicities (U α = 0 ∀α), the Hamiltonian (1) can be easily diagonalized. Irrespective of the signs of the tunneling matrix elements t α,β , there will be a lowest-energy normal mode z (possibly degenerate by symmetry). The ground state in the canonical ensemble of fixed phonon number n is simply (n!) −1/2 (b † z ) n |0 , which can be interpreted as a Bose-Einstein condensate (BEC) of phonons in mode z.
The BEC analogy can be extended to small but nonzero phonon-phonon interactions U α . With the canonical transformation a α = k Γ α,k b k that diagonalizes the harmonic part of the Hamiltonian (i.e., describes the normal modes), the Hamiltonian is re-expressed as with the coefficients U k3,k4 describing interactions between normal modes which deplete the BEC. When approximating the central region of the ion crystal by a quasi-infinite periodic d-dimensional lattice of N ions, the normal modes become Fourier modes b k = N −1/2 j a j exp(−i k · r j ). The interaction of these modes conserves linear momentum: wheret( k) = δ t( δ)e i k· δ are the Fourier transformed tunneling couplings, and In this approximation, a small interaction U partly depletes the BEC at k z and broadens the structure factor around this momentum component. The condensate depletion fraction is then [19] n ex n 0 where n 0 = n/N is the undepleted BEC density, and the dispersion relation is shifted such thatt( k z ) = 0.

3.3.
Large anharmonicity: spin-1 2 XY models In the limit of strong repulsion |U α | ≫ | β t α,β | ∀α and for low phonon filling n < 1, double phononic occupancy is strongly suppressed, and the system of phonons attains the hardcore limit, where phonons can be conveniently mapped onto S = 1 2 spins via the Holstein-Primakoff transformation For each ion, the vibrational ground state |0 is thus mapped onto the spin state |↓ , and the lowest excited state |1 onto the state |↑ . The resulting spin Hamiltonian has the form of an XY model with long-range dipolar couplings and a site-dependent field: Any spin-spin couplings for the z components must derive from off-site densitydensity interactions between the phonons, namely from terms of the form n i n j . The only source of such interactions is the fourth-order terms of the Coulomb interaction between the ions; as discussed in section 2, they are smaller than the XY interaction by a factor of (ζ/d 0 ) 2 and therefore negligible. Nonetheless the implementation of the Hamiltonian of (6) is of significant interest, because it offers the possibility of exploring the rich physics of frustrated XY spin models in various dimensions, with a fully tunable frustration. Detailed discussions of several such spin models and associated quantum phases and phase transitions are provided in the next section.

Frustrated XY Spin-1 2 Models
In this section we discuss three frustrated XY lattice spin models belonging to the family described by the Hamiltonian of (6) and realizable in linearly or planarly trapped ions. We henceforth assume for simplicity uniform local potentials V α = V in the trapping region; under the condition of preparing the system with a well-defined number of phonons, we can then discard the field term in (6). Moreover we focus on a half-filled system of hardcore bosons, corresponding to the zero-magnetization sector of the Hilbert space for the spin system. In the case of bipartite lattices, the ground state of an XY antiferromagnet in this zero-magnetization sector can be rigorously shown to coincide with that of a (more common) XY antiferromagnet without magnetization constraints and in zero field [20]. As for frustrated lattices, this remains certainly true for finite-size systems, where no spontaneous breaking of the Z 2 symmetry of the Hamiltonian can occur. We have focused on spin models with XY interactions, featuring full frustration for two equivalent (and non-commuting) spin components, a paradigm of quantum frustration. When using internal states of the ions to encode the spin variable (section 3.1), the effective spin-spin interaction which arises is typically of Ising-like 00 00 00 11 11   Linear chain. The Hamiltonian of (6) on a linear chain, obtained by strong transverse confinement of the ions, realizes an S = 1 2 one-dimensional XY antiferromagnet with dipolar interactions. If we neglect all interactions except nearest-neighbor ones, we recover an exactly solvable model [21] with power-law decaying correlations in the ground state, S z i S z i+r ∼ r −2 and S x i S x i+r ∼ r −1/2 . Long-range dipolar tunneling simply modifies the decay exponents [14], leading to a faster decay of S x i S x i+r and to a slower decay of S z i S z i+r .
Zig-zag ladder. A richer physics emerges in the case of a zig-zag ladder of ions, which develops naturally for weaker transverse confinement of the ions [22,23]. Assuming longitudinal inter-ion spacing d 0 and a zig-zag amplitude ξ×d 0 , inter-chain neighboring ions are distant by r 1 = d 0 1 + ξ 2 and intra-chain neighbors by r 2 = 2d 0 (see figure 1). The relevant phonon direction of motion m in (2) is taken perpendicular to the main trap axis and at an angle φ from the ion plane normal. Hence the ratio between the dominant intra-chain (t 2 ) and inter-chain (t 1 ) couplings is In trapped-ion experiments with cylindrically symmetric traps, at ξ ≈ 0.965 spontaneously generated zig-zag ladders deform into helices [23]. For φ = 0 (outof-plane vibrational motion) this limits the range of easily accessible coupling ratios to R 0.335. In practice, larger amplitudes ξ (and hence a larger R) can be engineered by breaking the cylindrical symmetry of the trap and constraining the ions to two dimension, or by microtrap stabilization [24]. However, choosing φ = π instead (in-plane vibrational motion) gives access to all ratios R > 1/8 for ξ < 1/ √ 2. In what follows we therefore assume in-plane vibrational motion. We use ξ and R interchangeably through (7). Figure 2 shows the normal-mode dispersion relation in the zig-zag ladder. For small zig-zag amplitudes ξ <ξ classical dipole = 0.461 (R < 0.353) the normal-mode dispersion t(k) has a unique minimum at k z d 0 = π, corresponding to Néel order; but for ξ >ξ classical dipole this minimum bifurcates into two symmetric minima at k ± z d 0 = π ± q(ξ), corresponding to spiral order. This is to be compared with a similar result for the system with intra-chain and inter-chain nearest-neighbor couplings only ({t 1 , t 2 } system), where the transition from Néel to spiral order occurs at a smaller value of R = 1/4 in the limit of classical XY spins [25]. A spiraling ground state is endowed both with spontaneous magnetic order, namely the breaking of the rotational symmetry in spin space, and with chiral order, corresponding to the choice of helicity of the spiraling state. In a quasi-1D system we define the chirality as [26] where the sites i and i + 1 belong to different chains. A striking result in quantum systems is that chiral order can survive even when magnetic order disappears. In fact, the {t 1 , t 2 } XY zig-zag ladder with S = 1 2 is found in [26] to exhibit power-law decaying spin-spin correlations for all values of R, but it is also found to develop longrange chiral order for R > 1.26. We extend the calculation of [26] to an XY ladder with dipolar interactions, making use of exact diagonalization up to sizes L = 20. We define an order parameter for the chiral phase as the averaged chiral correlation over all pairs of the system: (9b) Figure 3 shows this chirality order parameter as a function of the zig-zag amplitude. For the {t 1 , t 2 } model, we confirm a relatively sharp transition to a chirally ordered state for ξ ≈ 0.59 (corresponding to R ≈ 0.9), consistent with the transition of [26]. The full long-range dipolar model essentially retains this transition, shifted to slightly higher zig-zag amplitudes (ξ ≈ 0.63), along with what seems like a reorientation transition around ξ ≈ 0.8 with significantly higher chiral order, which is not observed in the classical limit. At present we do not provide a description of the phase 0.8 < ξ < 1, but we point out that despite the perturbative appearance of the dipolar interaction tails, they change the phase diagram of this zig-zag ladder drastically. This effect is currently the subject of further study.
Triangular lattice. When put in a planar trap, ions self-assemble into a triangularlattice Wigner crystal [27]. Adding an optical lattice to this setup, the degrees of freedom associated with the strongly non-linear vibrational motion in the lattice minima can realize a large family of spatially anisotropic XY triangular antiferromagnets. Unlike the case of the zig-zag ladder, we consider fixed average ion positions here, and the model parameters are varied by tilting the direction of motion of the ions m with respect to the lattice normal. We further assume that m remains perpendicular to one of the lattice directions, indicated by τ 1 , and makes an angle ϑ with the lattice normal (see figure 4). For ϑ = 0 the model of (6) realizes an isotropic triangular XY antiferromagnet with dipolar couplings, and it is maximally frustrated. For this system the minima of the normal-mode dispersion are on the corners of the hexagonal first Brillouin zone, corresponding in the classical limit to the well-known three-sublattice ordered state with 120 • angles between the spin directions on the sublattices [see figure 5(a)].
A nonzero angle ϑ produces a spatial anisotropy in the couplings, leaving the intrachain couplings along the τ 1 -axis unchanged while modifying the transverse (interchain) couplings. The value ϑ ≈ 31 • produces the minimum effective interaction between the τ 1 -chains, in the sense that the normal-mode dispersiont( k) is least a) Classical spin ordering: decoupled Néel-ordered horizontal (τ 1 ) chains for t 2 /t 1 = 0; three-sublattice ordering for t 2 /t 1 = 1; Néel order on the t 2 lattice for t 2 /t 1 ≥ 2. b) Proposed quantum phase diagram for S = 1 2 spins, including three spin-liquid (SL) phases.
corrugated in the direction perpendicular to τ 1 . We note that ϑ = 31 • is also the point where the BEC depletion of (5) grows fastest with U , approximating the 1D limit where the BEC is fully depleted for any nonzero interaction U . At ϑ ≈ 42 • the nearestneighbor inter-chain coupling vanishes, while residual weaker couplings to further neighbors survive due to the long-range nature of the dipolar interaction. On the opposite end, a value of ϑ = π 2 maximizes the inter-chain nearest-neighbor interaction. Thus by rotating the direction of vibration m with respect to the lattice plane we produce a similar effect to that produced in the zig-zag ladder by deformation of the Coulomb crystal. This different method is necessary because tuning the interaction coefficients in the effective spin Hamiltonian by mechanical deformation of the ion crystal is problematic, due to the significant stiffness of a two-dimensional Wigner crystal. The two approaches could possibly be combined in order to cover a larger range of system parameters.
At the classical level, varying ϑ from zero continuously deforms the threesublattice structure of the ground state of the isotropic triangular antiferromagnet, shifting the peaks in the structure factor away from the corners of the first Brillouin zone. For both zero and nonzero ϑ, the classical ground state has a finite vector chirality, defined on a plaquette with counter-clockwise labeled corners (i, j, k) as [28] and long-range chirality correlations, defined on plaquette pairs as [29]   where the plaquette pairs (△, ▽) and (△ ′ , ▽ ′ ) each share two spins. We find that even in the quantum model it is numerically irrelevant which edge is being shared (along the τ 1 or a τ 2 direction), as shown in figure 6.
In the quantum limit S = 1 2 we neglect the long-range dipolar tail of interactions, restricting our attention to the case of nearest-neighbor intra-chain (t 1 ) and interchain (t 2 ) couplings. According to the previous discussion on the dependence of inter-chain couplings upon the angle ϑ, one can experimentally span the range of parameter ratios 0 < |t 2 /t 1 | < 5/4; the sign of t 2 /t 1 is of no relevance since it can be reversed by a π rotation of every other τ 1 -row of spins [with a corresponding modification to (10)]. We have done calculations up to t 2 /t 1 = 3/2, which may be experimentally accessible by minor deformations of the Wigner crystal. Numerical treatment of the two-dimensional model with exact diagonalization is restricted to very small lattices, while other approaches such as quantum Monte Carlo are hindered by the sign problem due to frustration. We hence resort to a variational calculation based on projected-entangled pair states (PEPS) [30,31], which allows us to look at much larger lattices than with exact diagonalization. The calculations are performed on rhombic 10x10 and 20x20 lattices with open boundary conditions; such boundary conditions, while being essentially inherent to the PEPS method, are particularly welcome in a system developing incommensurate spiraling order due to frustration: in figure 6, exact results on a 4x4 lattice with periodic boundary conditions are shown to exhibit much reduced chiral order. Moreover, they mimic the natural boundary conditions realized experimentally in a finite Coulomb crystal. Due to the significant scaling of the computational effort with the bond dimension D of the variational PEPS basis, we limit ourselves to the smallest bond dimensions (D = 2, 3) going beyond the mean-field limit D = 1. Despite this aspect we still capture dramatic quantum features due to the interplay between frustration and quantum fluctuations.
In figure 6 we report the evolution of chiral order upon increasing the spatial anisotropy in the couplings. We observe strong indications for the breakdown of chiral order around t 2 /t 1 ≈ 0.6 and t 2 /t 1 ≈ 1.3 via quantum phase transitions. Indeed, at the classical level long-range chiral order is expected for all t 2 /t 1 < 2, as it is associated with the a spiral state with a well defined helicity of the spirals; chiral order is expected to vanish only for t 2 /t 1 ≥ 2 (square lattice XY antiferromagnet). Quantum mechanically, the limits t 2 = 0 and t 2 /t 1 → ∞ are well known to give a spin liquid with algebraically decaying correlations, and long-range Néel order [32], respectively. The way these two limits are attained is nonetheless highly non-trivial, as clearly shown by the evolution of spin-spin correlations upon changing t 2 /t 1 . In the following we will focus on the spin-spin correlations along the τ 1 direction (of the t 1 couplings) and along the τ 2 direction (of the t 2 couplings). In particular, we indicate with (i, j) the position of a point in the (τ 1 , τ 2 ) reference system, and for any inter-site separation δ we average over all pairs of sites at a distance of δ, along τ 1 and τ 2 . For instance, for correlations along τ 1 , we define in analogy to (9a) where N δ is the number of sites (i, j) δ in the sum for which (i + δ, j) is also within the simulation cell. Figure 7 shows the spin-spin correlations for different values of t 2 /t 1 . The correlations along the τ 2 axis build up very slowly upon increasing t 2 /t 1 , and they become significant only for t 2 /t 1 ≈ 0.7, consistently with the appearance of chiral order as observed in figure 6. Hence in the interval 0 < t 2 /t 1 0.6 the interplay between quantum fluctuations and frustration leads to an effective decoupling between the τ 1 chains. What is then the state along each τ 1 -chain? Its evolution with increasing t 2 /t 1 is absolutely non-trivial and completely dominated by quantum effects. In fact, classically a staggered ordering appears along each chain for t 2 = 0 and it turns continuously into (generally incommensurate) spiraling order for any finite t 2 . The spiral order reproduces the three-sublattice structure of an isotropic triangular lattice for t 2 = t 1 , and it realizes Néel order on the t 2 -lattice when t 2 /t 1 ≥ 2. This evolution is sketched in figure 5(a). Quantum mechanically, on the other hand, for t 2 /t 1 ≪ 1 the dominant correlations along each chain remain antiferromagnetic, and the decay of correlations becomes stronger for larger t 2 : as shown in figure 8, the decay is initially of the type 1/r K/2 with K = 1 at t 2 = 0. K is seen to grow upon increasing t 2 : the two-dimensional system hosts a one-dimensional Luttinger-liquid state on each of the τ 1 -chains, with a non-universal K exponent. Even more strikingly, the algebraically decaying staggered correlations turn to exponential decay for 0.4 t 2 /t 1 0.6, where the state on each chain evolves from a gapless to a gapped spin liquid. Only upon  increasing t 2 /t 1 further we observe that long-range correlations build up again, this time with a spiraling structure.
When increasing t 2 /t 1 above 1, a similarly non-trivial evolution of correlations takes place. We observe that the loss of chiral order, occurring at t 2 /t 1 ≈ 1.3, is accompanied by a loss in correlations along the τ 2 direction, rendering the τ 1chains again effectively decoupled. The state along each τ 1 -chain appears to be again a gapped spin liquid, with exponentially decaying correlations. It is only for t 2 /t 1 ≈ 1.5 that the system exhibits long-range correlations in both the τ 1 and τ 2 directions, the latter of which already exhibits a full Néel structure. Figure 5(b) summarizes our quantum phase diagram based on PEPS calculations, which exhibit three spin-liquid phases: a gapless one and two gapped ones. The inclusion of longrange dipolar interactions, which physically arise in the ion system, is at present technically challenging; however, we expect the long tail of interactions to give only small corrections to phases which develop a finite gap in the spectrum or finite longrange order in the model limited to short-range couplings, as is the case for most phases in figure 5(b).
The study of quantum antiferromagnetism on the anisotropic triangular lattice has recently attracted a deep interest in the case of Heisenberg interactions, in connection with the physics of quasi-2D antiferromagnets such as Cs 2 CuCl 4 [33]. In particular, a similar phase diagram to the one presented in figure 5(b) has been obtained in [34] by variational calculations. It is extremely intriguing to observe that one of the most striking features of the phase diagram predicted for the Heisenberg case [34], namely the emergence of intermediate gapped spin-liquid phases in the system, is also present in the case of XY interactions. Hence the evolution from a 1D gapless spin liquid to a spiraling ordered state, and from spiral to Néel order, apparently acquires a "universal" discontinuous structure: instead of a continuous deformation of correlations in the ground state, exhibited by the classical system, the quantum system first shows a complete loss of (quasi-)long-range correlations in favor of a short-range spin-liquid state, and then a revival of correlations at a different wavevector.

Preparation of the half-filled system
The results for the spin models in this section assume a ground state with zero total magnetization, which corresponds to a phonon model with population N/2 assuming N active degrees of freedom. Such a state may be more easily prepared in the absence of anharmonicities, as a BEC of phonons described in section 3.2. The optical lattice constituting the microtraps is subsequently switched on adiabatically, which transforms the BEC into the half-filled ground state of the Bose-Hubbard Hamiltionian (1) and, for very large anharmonicities, into that of the S = 1 2 spin model.
The phonon BEC at the starting point of this adiabatic passage cannot be prepared by cooling a thermal population of phonons, since their number is not conserved during the cooling process. Instead, we propose that the harmonic trap holding the ion crystal is opened to the point where just the lowest-frequency vibrational mode becomes soft. In a linear chain, for example, we open the radial confinement to the point where the chain acquires a small zig-zag amplitude [22,23]. The ion crystal is then cooled to its vibrational ground state in the new symmetrybroken configuration, after which the trap is rapidly closed to its original state. Depending on the amplitude of this symmetry breaking, the final state of the system is a coherent state of phonons in the lowest-frequency mode. In a sufficiently large system, adjusting the mean number of phonons to the half-filling point results in a coherent state which is not very far from a Fock state of exactly half filling. A systematic discussion of this procedure is forthcoming [16].

Discussion and Conclusions
In this paper we have presented a variety of quantum many-body Hamiltonians realized by the anharmonic vibrations of ions trapped in the wells of a deep optical lattice. While the dominant Coulomb interaction and the overall trapping potential determine the geometry of the ionic Wigner crystal, at the same time the optical lattice allows to engineer the phononic states for the ions, and it endows the system with a broad range of tunability. On the one hand one can control the non-linearity of phonons (their effective on-site repulsion/attraction) by modifying the optical lattice or superlattice potential. In this way one can explore a regime of essentially harmonic phonons, which can mediate anisotropic spin-spin interactions between internal degrees of freedom of the ions [5]; a regime of weakly interacting phonons, which, in the case of a two-dimensional triangular-lattice Wigner crystal, condense in a non-trivial finitemomentum state containing vortex-antivortex pairs [figure 5(a)]; and a regime of essentially hard-core phonons, which allows a mapping of the phononic Hamiltonian onto a frustrated S = 1 2 XY model with spatially anisotropic couplings, tuned by changing the orientation of the optical lattice with respect to the ionic crystal. In this latter case we show that an extremely rich physical picture can be obtained in the case of a half-filled system of phonons: in particular, one can realize a gapless spin-liquid phase with long-range chiral order in zig-zag ion ladders, and a variety of two-dimensional spin-liquid phases, either gapped or gapless, in the case of a triangular Wigner crystal. Hence this setup offers a platform for the exploration of highly non-trivial quantum phases and quantum critical phenomena in tunable frustrated quantum magnets.
As an outlook, it is important to realize that an optical lattice potential is just a simple instance of a bigger family of microtraps, in which each ion is trapped individually in a tunable local potential. In the case of an optical lattice all individual trapping potentials are typically identical, so that only global tunability is possible. On the other hand it is possible to fabricate independent microtraps for the ions [24,8,12,13,35,15,36], which guarantees at the same time individual tunability of the phononic states of each ion independently, and the possibility of realizing ion-lattice geometries which go beyond the natural Wigner crystals, and even beyond the paradigm of translationally invariant arrays. Work is in progress [16] to identify the quantum many-body Hamiltonians that can be realized in the context of microfabricated traps.