Theoretical basis for quantum simulation with a planar ionic crystal in a Penning trap using a triangular rotating wall

# Theoretical basis for quantum simulation with a planar ionic crystal in a Penning trap using a triangular rotating wall

A. Khan Department of Physics, Indian Institute of Technology, Roorkee-247667 Uttarakhand, India    B. Yoshimura Department of Physics, Georgetown University, Washington, DC, 20057 USA    J. K. Freericks Department of Physics, Georgetown University, Washington, DC, 20057 USA
July 15, 2019
###### Abstract

One of the challenges with quantum simulation in ion traps is that the effective spin-spin exchange couplings are not uniform across the lattice. This can be particularly important in Penning trap realizations where the presence of an ellipsoidal boundary at the edge of the trap leads to dislocations in the crystal. By adding an additional anharmonic potential to better control interion spacing, and a triangular shaped rotating wall potential to reduce the appearance of dislocations, one can achieve better uniformity of the ionic positions. In this work, we calculate the axial phonon frequencies and the spin-spin interactions driven by a spin-dependent optical dipole force, and discuss what effects the more uniform ion spacing has on the spin simulation properties of Penning trap quantum simulators. Indeed, we find the spin-spin interactions behave more like a power law for a wide range of parameters.

###### pacs:
37.10.Ty, 03.75.âb, 03.67.Ac, 03.67.Lx

## I Introduction

The idea of a quantum simulator, where a complex many-body quantum system is emulated in a controlled analog quantum computer and the results of the simulation are read off of the computer by measuring different properties as a function of time, originates with work from Richard Feynman feynman () in the early 1980s. Cirac and Zoller cirac_zoller () showed how ion traps driven by a spin-dependent optical dipole force could realize quantum computers. Porras and Cirac porras_cirac () further described how one could perform quantum simulations in a Penning trap. One of the issues with these quantum simulations is that the ions are not spaced uniformly. On the one hand, this leads to nonuniform effective spin-spin couplings between the ions, on the other hand, in one-dimensional linear Paul traps, it leads to the linear to zig-zag transition, which limits the number of ions that can be held in the trap in a one-dimensional linear configuration. It was quickly realized that by adding an anharmonic potential, which pushes together the farther out ions preferentially when compared to the central ions, one can achieve a more uniform arrangement, and the precise potential for perfectly uniform trapping is known for the linear Paul trap duan_monroe (). Surprisingly, one can achieve quite uniform crystals by just adding a quartic potential on top of the conventional quadratic trapping potential. This ideology has been extended to the Penning trap by Dubin dubin (), where he also included a triangular-shaped rotating wall potential to reduce dislocation formation, which occurs in two-dimensional Penning traps when the boundary potential does not have the same symmetry of the underlying ionic lattice.

In this paper, we extend the analysis of Dubin to determine the behavior of different numbers of trapped ions, different wall potentials and different rotation rates, to determine the stability of these ionic crystals. We further calculate the axial phonon modes and from them the effective spin-spin interactions induced by a state-dependent optical dipole force. We end by discussing the feasibility of the triangular wall for quantum simulation with ionic crystals in the Penning trap. The organization of this paper is as follows: in Sec. II, we describe the theoretical background for the calculations. In Sec. III, we present the numerical results for the calculations of the spin-spin couplings of the ions. In Sec. IV, we present our conclusions.

## Ii Theoretical Formulation

The Penning trap confines ions by using an electrostatic potential that pushes the ions towards the plane with , and also pushes the ions outwards, radially. A large static magnetic field curves the radial motion into circles, which result in a trapped ion crystal (after taking into account the Coulomb repulsion of each ion). An additional rotating wall potential, with specified angular symmetry, is then applied to control the shape of the crystal and its rotation rate. While many different ions can be employed in a Penning trap, we will focus on the realization with two hyperfine levels of Be ions, and , localized to a single plane. An extensive description of such a set-up can be found elsewhere bollinger_penning_trap (). Cold atoms condensing in such crystals are candidates for building quantum simulators, owing to the ease with which these systems can be prepared for large number of ions, and the precise quantum control of individual ions these systems afford britton ().

In actual experiments, the ionic crystal often acquires an increasing number of impurities as a function of time, such as BeH that form due to collisions of the beryllium ions with hydrogen molecules that are in the dilute background gas. While a theoretical treatment that includes the effect of these impurities is possible mcaneny_yoshimura_freericks (), we consider only the clean limit here, where there are no impurities. This simplifies the analysis below. In addition, experimental protocols to purify the systems may make this effect less important bollinger_photodissociation ().

The theoretical treatment of the equilibrium positions and normal modes of the Penning trap requires a careful analysis employing standard classical mechanics and then an appropriate quantization scheme. Details for how to do this have already appeared russian (); taylor (); wang_keith_freericks (). Here, we provide a quick summary of that formalism to establish our notation and to show how the approach needs to be modified for the more uniform triangular crystals that can be generated in an extra anharmonic potential with a rotating triangular wall. We begin with the ion Lagrangian in the laboratory reference frame which satisfies

 L=N∑j=1[12m˙r2j−eϕj(t)+eA⋅˙rj] (1)

where is the total number of ions, is the (positive) unit charge of an electron and is the mass of a Be ion. The symbol is the position vector for the th ion in Cartesian coordinates, is the total scalar potential acting on the th ion (including the rotating wall potential), and is the vector potential in the symmetric gauge for the uniform axial magnetic field with . The scalar potential includes the potentials that trap the ions and the mutual Coulomb repulsion between the ions. It can be expressed as follows:

 eϕj(t)= V0 [z2j−12ρ2j]+12mω2effC4ρ4j (2) + VWρ3jcos[3(θj+Ωt)]+kee22∑k≠j1rkj

where is the effective trapping frequency in the rotating frame of the crystal with (see below), is the amplitude of the static quadratic potential from the Penning trap electrodes, is the strength of the additional fourth-order anharmonic trapping potential (also coming from the Penning trap electrodes), is the amplitude of the triangular rotating wall potential, is the rotating wall angular frequency (which rotates about the -axis), and is the Coulomb force constant. Here, is the interparticle distance between the th and th ion and is given by , is the polar coordinate radius for the th ion and is given by , and is the polar angular coordinate for the th ion which is given by . The rotating wall potential makes the potential time-dependent in the laboratory frame. The solution for the ion positions is simplified by transforming to the equivalent equilibrium problem in the rotating frame with angular speed (where the effective trapping potential becomes time-independent). Transforming to the rotating frame, we arrive at the following time-independent rotating-frame potential for the th ion (which is confined to the plane with :

 eϕj = 12mω2eff[(ρRj)2+C4(ρRj)4] (3) + VW[(xRj)3−3xRj(yRj)2]+kee22∑k≠j1rRjk

where is the transformed set of coordinates for the rotating frame. The -dependence of the effective trapping frequency is due to the way velocities transform for rotating frames, which now includes the effects of potential terms from the centrifugal and Lorentz forces as well.

Comparing Eq. (3) with a similar expression appearing in Ref. wang_keith_freericks, , we see the following differences: (i) there is an additional anharmonic, fourth-order term that causes the ions in the outer regions of the crystal to be pushed in more strongly, and hence counteracts some of the inhomogeneities that occur due to increasing interion distances as we move outwards; (ii) the angular shape of the rotating wall term has now been adjusted to an angular harmonic as the crystal condenses into a triangular lattice and there is less frustration at the edges if the rotating wall has a symmetry that matches that of the underlying crystal facets dubin (). We note that the additional/modified terms retain their form under a transformation from the laboratory frame to the rotating frame. We also remark that the unconventional choice of (with the coefficient dependent on ) for the anharmonic term was made in anticipation of the simpler and rather standard form we get in Eq. (3).

The ions in a Penning trap crystal do not always crystallize in a two-dimensional plane. The following approximate criterion is usually required to be satisfied for a planar only configuration:

 2eV0(eBzΩ−mΩ2−eV0)≫1 (4)

which basically ensures that the restoring force in the axial direction is several orders of magnitude bigger than in the radial direction, and the crystal lies in the plane only.

The triangular rotating wall also introduces an anisotropy in the radial potential, with deconfinement along certain directions which, if strong enough, can lead to particle loss. For the more familiar case of , this is reflected in the increasing eccentricity of the elliptical equilibrium structures that suggests deconfinement along the ”weak” axis. This is expressed in terms of an approximate criteria in Ref. wang_keith_freericks, for the simplest rotating wall. For , this becomes difficult to express in terms of a single criterion and we resort to a numerical calculation to find the triangular ”separatrix”, the contour lines of the radial part of the gradient of the potential function in Eq. (3), with the Coulomb repulsion terms ignored. These are shown for varying ratios of strengths of the rotating wall term to the effective radial confinement strength, in Fig. 1. We notice that unlike the quadrupole rotating wall potential, deconfinement happens for the triangular wall at large enough distances for all wall amplitude strengths. We see that as the strength of the triangular rotating wall increases, the separatrix moves closer to the center of symmetry of the crystal, with apparent deconfinement centered along the and axes. This has a noticeable effect on the ”shape” of the Penning-trap crystals, which reduces the dislocations in the crystal and helps maintain the uniformity. Of course, a full description of stability of these planar structures requires inclusion of the Coulomb terms. The accurate quantitative description of stability requires solving for the equilibrium positions presupposing a planar arrangement of ions, and then showing that the eigenvalues of the normal modes of oscillation about these equilibrium positions are all positive; that is, the phonon normal mode frequencies are all real.

The full, transformed Lagrangian for the rotating frame is

 LR=N∑j=1[12m|˙rRj|2−eBeff(Ω)2(˙xRjyRj−˙yRjxRj)−eϕRj] (5)

where is the -dependent effective magnetic field in the rotating frame. The modification of the magnetic field is due to velocity dependent terms in the laboratory frame Lagrangian. This affects the oscillating normal modes of the planar motion when the ions are far from their equilibrium positions. However, as we will see soon, this does not have an effect on the axial modes. This observation greatly simplifies the normal-mode analysis for the axial modes.

To find the stable spatial configuration of the ions, we minimize the effective potential energy in the rotating frame of reference. This is a challenging optimization problem to solve for in two (and higher dimensions), especially since different configurations, separated by large potential barriers, can have local minima in the potential energy function with small energy differences to the global minimum. We follow the previous treatments of this problem, where the experimental indication of the fact that the ions condense in a triangular lattice in a single plane is used to construct the optimized solution that lies close to a perfect triangular lattice.

We construct an initial, trial solution based on the “closed-shell” approximation as in Ref. wang_keith_freericks, , but with the important difference that the overall “shape” is triangular and not hexagonal, as it was in the previous solutions. This is to reflect the fact that the overall shape of the crystal is dictated by the equipotential lines of the rotating wall term, which in this case is triangular.

We then proceed to calculate the collective normal mode excitations of the crystal. The ion Lagrangian (in the rotating frame, all superscripts are dropped for clarity) is expanded via a Taylor series about the previously calculated equilibrium positions of the ions up to quadratic order. The ion coordinates are, for the purpose of the expansion, written as , while for the ion velocities we write . Because we are expanding about equilibrium, we can drop the linear terms in the coordinates, and we find

 (6)

where the is due to the equilibrium state and the quadratic terms are due to fluctuations away from equilibrium, which we henceforth call for the phonon Lagrangian. The Lorentz force due to the external magnetic field lies in the plane, and the potential energy is clearly seen to be separable in cylindrical coordinates. This means that the axial phonon Lagrangian can be decoupled from the planar phonon Lagrangian, and there is no harmonic coupling between the planar and axial degrees of freedom. Therefore, we can study the axial and planar modes independently and we can solve just the equations of motion for the axial or the planar modes independent of the other).

We examine only the axial modes in this work. This is due to the fact that the planar modes have a complex structure owing to a coupling of the ion motion in the and directions, the appearance of velocity-dependent forces, as well as complexities introduced by rotation of the ion crystals as observed in the laboratory frame. The fact that there is no harmonic coupling between axial and planar directions of the crystal allows us to exclude the planar modes from our discussion henceforth. Restricting to the axial modes only, is further supported by the fact that the simplest form of quantum simulation works on driving the axial modes with a state-dependent optical dipole force.

The axial Lagrangian is then given by

 Laxialph=12N∑k=1m(δ˙Rzk)2−12N∑j,k=1KzzjkδRzjδRzk (7)

where the spring constants satisfy

 Kzzjk=−∂2L∂Rzj∂Rzk∣∣ ∣∣0 (8)

The absence of cross terms in the velocity part of the Lagrangian can be easily seen from the sum-of-squares structure of the kinetic energy along the -direction and the fact that there are no velocity-dependent forces in the -direction.

An explicit calculation for the matrix elements of the gives the following:

 Kzzjk=⎧⎪ ⎪⎨⎪ ⎪⎩2eV0−kee2∑Nk′,k′≠j1(R0jk′)3j=kkee21(R0jk)3j≠k (9)

where is the distance between ions located at their respective equilibrium positions in the rotating frame. We see that the axial stiffness matrix is Hermitian and symmetric, and is independent of the anharmonic or wall potentials.

To solve for the axial ion normal modes, we apply the Euler-Lagrange equations to the axial phonon Lagrangian in Eq. (7):

 mδ¨Rzj+N∑k=1KzzjkδRzk=0,j=1,2,....N. (10)

which, on substitution of the eigenvector solution ansatz , gives

 N∑k=1[mω2zνδjk−Kzzjk] bzνk=0,j,ν=1,2,…N (11)

where is the normal mode eigenfrequency and is the th axial normal-mode eigenvector. The eigenvalue problem is quadratic, but we can easily map it onto a linear eigenvalue problem by setting the eigenvalue according to . We can then solve for the eigenvalues and eigenvectors numerically in MATLAB. The eigenvectors are real, N-tuples whose norm has been set to unity by convention. The eigenvalues are positive for stable normal modes and negative for unstable normal modes.

The quantization of the normal modes is completely standard: we first identify the positions and momenta associated with each phonon mode as canonically conjugate, and promote the relation given by the Poisson bracket to the commutation relation for the operators and , . To calculate the canonically conjugate variables for the phonon modes, we make the transformation where are the normal coordinates for each phonon mode . We see that the Lagrangian assumes the following diagonal form:

 Laxialph=12N∑ν=1m(˙ξ2ν−ω2zνξ2ν). (12)

Hence, we calculate the conjugate momenta as follows:

 Paxialν=∂Laxialph∂˙ξν=m˙ξν. (13)

The Hamiltonian is then expressed as

 Haxialph=N∑ν=1((Paxialν)22m+12mω2zνξ2ν). (14)

To quantize the normal modes, we identify that the Hamiltonian is a sum of simple harmonic modes with frequencies . We now introduce creation and annihilation operators as follows:

 ^azν=√mωzν2ℏ(ξν+imωzνPaxialν) (15)

and

 ^a†zν=√mωzν2ℏ(ξν−imωzνPaxialν). (16)

Hence, the quantized Hamiltonian operator can be expressed as

 ^Haxialph=N∑ν=1ℏωzν(^nzν+12) (17)

where is the number operator. The operator for displacement along the -direction can be expressed in terms of the creation and annihilation operators as follows:

 δ^Rj=N∑ν=1bzνj√ℏ2mωzν[^a†zν+^azν]. (18)

We also note that the form of the Hamiltonian derived for the axial modes here is invariant when we transform coordinates from the lab to the rotating frame, since the ion oscillations are only along the rotation axis (-direction) and hence are not influenced by rotation of the coordinate axes.

We now need to calculate the effective spin-spin coupling between the ions, generated by the spin-dependent optical dipole force. This analysis has been done in detail elsewhere porras_cirac (); monroe (); wang_freericks (), which we utilize here. The effective spin Hamiltonian is dictated by a time-dependent Ising spin Hamiltonian,

 H(t)=N∑j,j′=1Jjj′(t)σzjσzj′, (19)

where the Ising spin-spin coupling between sites and is given by

 Jjj′(t) = F2O4mN∑ν=1bzνjbzνj′μ2−ω2zν[1+cos(2μt) (20) − 2μωzνsinωzνtsinμt].

Here is the magnitude of the optical dipole force, and is the beat-note frequency corresponding to the frequency difference of the two off-resonant laser beams being applied to the trapped ion crystal. We see that this expression relates the strength of the Ising-like coupling between ions to the phonon mode properties , and , which are calculated from the classical, normal-mode analysis described above. The time-averaged spin-spin couplings are given by the first term in Eq. (20). We can think of the effective spin-spin Hamiltonian as the expression in the parenthesis of Eq. (19) with the time-dependent spin-spin interactions replaced by the time-averaged ones.

## Iii Numerical Results and Analysis

We consider Be ions localized in a plane by the Penning trap potential defined in Eq. (2). We characterize the strength of the end-cap potentials that affect the axial trapping by a characteristic angular frequency , such that . This is fixed for the purpose of all our numerical calculations at the value  kHz, a typical value used in experiments. We normalize subsequent frequencies in terms of . The experiments at NIST typically run at rotational frequencies , and we have concentrated on regions close to this value in our calculations for experimental relevance. The cyclotron frequency associated with the magnetic field is defined as . Fixing  T, we get . The beryllium atom has an atomic mass  a.u. and a positive unit charge  C.

For the strength of the anharmonic term, we use the value in Ref. dubin, , where in the following form of the potential :

 ε = N∑i=1⎡⎣12mω2eff(ρ2i+38~C4ρ4ir2p)+kee22∑j,i≠j1rij (21) + VWallr3pρ3icos(3θi)]

where is the plasma radius parameter whose value is taken to be  cm. We now define a typical length and energy scale,

 l0=(kee2mω2z)13  % and  E=mω2zl20. (22)

Henceforth, we express all lengths in units of and all energies in units of . The potential then becomes

 ~ε=εE = N∑i=1⎡⎣12ω2eff(ρ2i+C4ρ4i)+12∑i≠j1rij (23) + VW(x3i−3xiy2i)]

where the dimensionless parameters and are given by the relations and . All lengths (), energies () and frequencies () appearing in Eq. (23) are expressed in dimensionless units. For the choice of , we have . For the strength of the wall potential, we use two different values, and . These have been chosen to show clearly the effect of variation in the wall strength on key ion-crystal characteristics, while also ensuring that all crystal structures have a stable equilibrium.

We report the rotational frequency in terms of , normalized with respect to . We stick to regions close to the experimental value of , which translates as . We use values for in the range in our analysis.

It is also useful to consider the stability of the crystal under the trap potentials we have used. Because the rotating wall potential varies as the third power of the coordinates, there is no deconfinement frequency as in the case of the quadrupole wall, trap. To probe the stability under deconfining forces, we look at the radial component of the force on the th ion due to the trap potential, excluding the Coulomb potential. which is given by

 Fr=−ω2eff⎡⎣ρi+2C4ω2effρ3i+3VWω2effρi(x3i−3xiy2i)⎤⎦. (24)

If we plot the locus of points where this function becomes zero for various values of , we see that we get three regions, each centered along the axes. For increasing strength of the rotating wall potential, they move closer to the origin and the radius of the separatrix (the smallest distance to these unstable zones along any axis) is seen to decrease. Hence, the deconfinement increases for very high values of the rotating wall potential, which is what we expect. In our analysis, we stay in regions where the extent of crystal is much smaller than the separatrix radius.

We also need to exclude unstable equilibrium configurations (which are indicated by the nonpositivity of the eigenvalues of the stiffness matrices of the axial vibrations). We examine this in detail below, when we discuss our results on the axial phonon modes.

To find the equilibrium configurations, we need to minimize the Hamiltonian of the crystal in the rotating frame, which boils down to finding the best minimum of the transformed potential function of Eq. (23) near a triangular lattice. We only concentrate on the solution we obtain starting from the closed-shell construction. In this case, we start with a seed lattice where we arrange ions in closed, triangular shells while also respecting the triangular lattice symmetry. The shell is closed if we can put all the ions in these complete shells. We relate the number of shells to the number of total ions in the crystal as

 S=⎡⎣√2(N−1)3+14−12⎤⎦. (25)

If we cannot put the ions in an integer number of shells, we put the outermost ions in an incomplete triangular ring according to the minimal potential energy at each of the outer ring sites. For the purpose of our discussion here, we pick and hence . We arrive at the minima guaranteed under such a consideration of the seed lattice using a trust-region algorithm of the MATLAB Optimization Toolbox. The minimization procedure requires us to specify a locally calculated gradient of the potential, which can be input analytically by taking derivatives of the potential. The procedure iterates the minimization steps until the local minimum is found.

The equilibrium configurations we obtain from such a procedure behave as we might expect (see Fig. 2). We obtain structures that form a nearly perfect triangular lattice close to the center, and smoothly transition to the shape of the contour lines of the effective potential as we move radially outward to the edges. The edge effects cause the interionic distances to change as we move outwards. Our strategy to counter these effects are two-fold, as discussed earlier: (i) We introduce a weak anharmonic term whose strength is characterized by and (ii) we match the symmetry of the rotating wall with that of the condensed crystal. The first effect has been incorporated and fixed at a particular value, as discussed earlier. However, we keep the strength of the rotating wall variable as it is only for a certain range of values of the strength of the rotating wall potential that the contour lines of the effective potential are triangular, and hence become least likely to cause edge distortions of the ionic crystal.

This is clearly seen in the Fig. 2, where we show the progression of these structures for different values of the effective radial trapping strength , for and . For the low-strength rotating wall, the structures (a)-(d) are nearly triangular and uniform, whereas the higher value of the rotating wall strength corresponds to the more distorted structures of (e)-(h). For a fixed number of ions and fixed , we see the destabilizing effect of the decrease in separatrix radius with increasing rotating wall strength clearly on these crystal structures. We note here that in the limit of vanishing wall strength or high effective radial trapping frequency, the structures become more isotropic, with uniformly decreasing nearest-neighbor distances as we move towards the edges. The other extreme (very small radial trapping frequency or high strength of rotating wall), the equilibrium configuration we obtain from a closed-shell construction shows that the ions are reduced to (three) pockets of stability and the structure is no longer closed. A detailed normal-mode analysis (see below) shows that these structures are in unstable equilibrium, and hence we can discard them. Another important observation is that a given progression of structures (for differing values of the radial trapping strength and increasing values of ) displays similar structures to those found at smaller trapping strengths, for higher values of the rotating wall strength. This fact will be important to arrive at structures that show the maximum uniformity and also exhibit a stable equilibrium.

Note that in these calculations, we fix the wall potential and then vary . In doing so, we find that we are limited by how many ions we can hold in stable equilibrium. Because the radial deconfinement decreases as the rotational frequency decreases, the effect of the rotating wall will become stronger if it also remains fixed. We have done this here to reduce the parameter space we explore. But, if one wants to examine larger crystals, then one needs to carefully tune the rotating wall strength as the rotational frequency is changed, as well as the strength of the quartic potential, to be able to continue to maintain stable equilibrium. These issues are discussed in Ref. dubin ()

The equilibrium structures we obtain are markedly uniform, and this is born out in Fig. 3 where we plot the distance to the nearest neighbor for each ion in the crystal as a function of the central ions’ distance from the trap symmetry axis. These results are for the most uniform stable structures we could obtain for both the and rotating walls. The nearest-neighbor distances have been calculated based on the Delaunay triangulation algorithm. Each point in the figure represents the distance of the ion in question to an ion in the first nearest-neighbor shell. We focus on the first circle of nearest neighbors only. The values of the relevant parameters, for both the (quadrupole) and the (triangular) wall crystals, are indicated in the caption to Fig. 3. The larger spread of values for the quadrupole wall tells us that, on average, the triangular wall crystal is indeed more uniform spatially than the quadrupole wall crystal. Note that there are the same number of blue squares and red circles in Fig. 3; the uniformity of the triangular lattice has many of the symbols overlap.

We next discuss the features of the normal modes of small oscillations of the ions. We first examine the positivity of the eigenvalue spectrum of the stiffness matrix of axial vibrations as a function of effective trapping strength, for high and low rotating wall potential strengths, in Fig. 4. We see that the eigenvalue spectrum is real and positive for only a “band” of values of (indicated in blue in the figure), and this band shrinks for the higher amplitude rotating wall potential. The positivity of all the eigenvalues indicates stability of the corresponding structures.

This behavior is fundamentally different from that seen in the usual Penning trap crystals where we do not see a lower bound to the strength of the radial trapping, characterized by , where the structures become unstable, although there is an upper limit. Hence, we can only talk of stable structures for certain ranges of for the rotating wall (with an additional quartic potential), and this range gets narrower for increasing wall strength. A similar effect is seen if the number of ions in the trap is increased, in that the band of stability shrinks as we increase the number of ions in the trap. The value of was the maximum number of ions we found that could be trapped with reasonably large bands of stability for the particular ranges of parameters that we chose. Note that more ions can be trapped by carefully choosing the rotation frequency , the rotating wall potential amplitude, and the strength of the quartic potential, but we do not discuss these cases in detail here.

Next, we plot the numerically obtained eigenfrequencies against their mode numbers, for crystal structures corresponding to two values of (stable) , for both high and low rotating wall potentials. Roughly, we can see a trend similar to that of the quadrupole-wall crystal normal modes, where the primary dependence is on and not on . However, there are important distinctions to be made. The “band” structure of the eigenvalue spectrum causes the structure of these curves to change, and quite significantly at that, when we vary the rotating wall from low to high strength. We see this in Fig. 4, where for , the lower edge of the instability of the band shown in Fig.4 shifts to the right when the strength of the rotating wall is increased, and this causes the eigenfrequencies for the higher wall strength to drop abruptly to values very close to zero. In this fashion, we see that the dependence on is now superimposed on a dependence on the strength of the rotating wall due to a -variable bandwidth.

The highest axial mode has a universal eigenfrequency equal to the angular frequency of the trapping strength, for all values of and , and we see that all the branches converge to this point. This behavior is identical to that of the quadrupole-wall rotating crystal. The corresponding eigenmode is the well-known center-of-mass mode, where all ions move equal displacements that are in-phase with each other. This is because the center-of-mass motion does not cost any additional Coulomb energy and all ions have the same axial trapping energy, and hence their motion is independent of the strength of the rotating wall potential applied in the crystal plane. Also, other axial phonon modes will have frequencies lower than the center-of-mass mode, as the average distance between the ions increases when the wave vector is nonzero and there is a reduction in energy due to the Coulomb repulsion. This is also why the axial eigenfrequency branches for higher lie roughly lower in Fig. 5.

Next, we discuss the eigenvectors of the axial mode phonons, which we can calculate immediately from the diagonalization of the stiffness matrix that yields the eigenfrequencies. Modes close to the center-of-mass mode are collective, where ions move with a long-wavelength. It is these long-wavelength modes that are important for purposes of quantum simulation, and we concentrate on these in our discussion here. In Fig. 6, we show maximum axial displacements of the ions, corresponding to the highest three axial modes, with the displacements normalized to unity and color-coded as indicated in the adjacent colorbar. We do this for the two different strengths of the radial trapping strength , for both the high and the low rotating wall potentials . The highest axial mode corresponds to the center-of-mass motion. The other two modes (so-called tilt modes) are nearly degenerate, and are seen to have similar nature even when anisotropy due to the rotating wall is dominant. This is also seen in Fig. 5, where we see that the two modes just below the center-of-mass mode have the same eigenfrequencies. This behavior of the penultimate axial modes is different from the quadrupole-wall crystal, where this degeneracy of the modes is lifted under the corresponding anisotropic rotating wall, and an additional mode is sometimes introduced between them.

Finally, we examine the strength of the effective spin-spin coupling that results from applying a spin-dependent dipole optical force detuned close to the axial phonon modes. We focus only on detuning frequencies to the blue of the center-of-mass mode () where we expect to find a power law dependence of the spin exchange on the distance between the spins in the lattice porras_cirac (). This behavior has been predicted and verified in experiments, and we expect it to hold even for this triangular-wall anharmonic crystal. In fact, the objective for making the lattice distances more uniform, is to also make the spin-spin couplings between adjacent spins more uniform, and this is would follow immediately if there exists a power law relation between the two ().

The parameters of axial phonon modes discussed already are used to calculate the static spin-spin interaction between the spins of ions and , based on Eq. (20). In Fig. 7, we plot this static interaction strength (expressed on a logarithmic scale) as a function of the distance between the ions, for various values of detuning larger than the center-of-mass frequency. We do this for the stable crystal structure exhibiting minimum variance in the nearest-neighbor distances. We see a behavior very similar to the quadrupole-wall potential. A uniform , independent of indicates that the detuning laser excites only the uniform center-of-mass mode. As we move away from the center-of-mass mode, an increasingly large number of eigenmodes participate in the coupling, and we see a clear trend in the values, such that . In the limit of large detuning, we have dipole-dipole interactions where tends to a value of 3. For small detunings, we have the all-to-all case of . There are small departures from the power law behavior for intermediate detunings, while the very small and very large values of detuning show excellent agreement with the power law. Note that these results are indeed much more uniform than what was found for a a quadrupole rotating wall britton ().

In Fig. 8, we plot the fitted power-law exponent versus the strength of the detuning away from the center-of-mass mode. The trend we see here is similar to the one seen in calculations for the quadrupole-wall potential. We see a faster approach to the dipole-dipole limit () for both smaller effective (radial) trapping frequencies, and weaker strengths of the triangular rotating wall potential, just like for the quadrupole rotating wall.

We have already noted that there are deviations from the power law behavior for intermediate values of detuning, and this is apparent in the spread of values away from the linear fit in Fig. 7 (especially for ). To explore these deviations in more detail, we plot the normalized root mean square deviation (RMSD), defined by,

 normalized RMSD=√∑i

as a function of the detuning in Fig. 9. We see an adherence to the power law (characterized by values of the normalized RMSD close to 0) for both small and large detunings . More importantly, we see the largest deviation in the normalized RMSD parameter for strengths of detuning in the intermediate range of to , We can understand this behavior easily if we look at the structure of the static part of Eq. (20), that relates the strength of the spin-spin coupling, , to the normal mode properties of the crystal. Each term in the summation can be understood to correspond to an eigenmode’s contribution to the coupling strength. When the beat-note frequency is very close to (), only the center-of-mass mode (corresponding to uniform motion of the ions) contributes, and the does not depend on distance. This behavior corresponds to a value of . For farther away from the center-of-mass mode, the lower modes begin to contribute increasingly. When only a few modes contribute, we cannot expect the power law behavior to hold mcaneny_freericks (). The structure of the eigenvectors, as we see in Fig. 6, is clearly not compatible with the power law decay of with distance. For many ions sitting at opposite edges, there is a large whereas the coupling is virtually zero for other pairs separated by much smaller distances. This is the origin of the spread of spin-spin couplings in Fig. 7 and the increase in the value of the normalized RMSD in Fig. 9. As we move towards larger values of (and ), we see that all the modes begin to contribute almost equally, and in the limit of dipole-dipole interactions, we get a value of very close to 3, with the normalized RMSD close to 0.

We note that such behavior is independent of the details of the crystal structure itself. The eigenvectors corresponding to the first few modes will have a structure independent of the details of the trap potential, and will cause a similar deviation from the power law behavior as we saw above. This deviation implies that the spin-spin couplings are no longer correlated with the distances between the ions, and hence an increase in spatial uniformity of the crystal is not guaranteed to have a bearing on the uniformity of the spin-spin interactions. This is an important observation, as the increasingly uniform nearest-neighbor distances for the triangular wall crystal would imply a more uniform spin-spin coupling strength only for detuning strength values that are moderately large. For intermediate values of , it is important to consider the nature of modes just below the center-of-mass mode to describe the spin-spin coupling strength between ions corresponding to that strength of detuning.

## Iv Conclusions and Discussion

In this work, we have examined the properties of a Penning trap with an additional anharmonic and triangular rotating wall potential which provide a much more uniform ionic crystal for use in quantum simulation. By performing a detailed analysis of the equilibrium positions, the phonons, and the effective spin-spin interactions, we find that indeed one can generally obtain more uniform spin-spin coupling constants. As one might have predicted, the relationship between ionic spacing in the lattice and the uniformity of the spin-spin interactions is not directly one-to-one. For small values, it is the character of the phonon eigenmodes that lie close to the center-of-mass mode that determine the behavior of the spin-spin couplings more than the interparticle spacing. We hope that the result of this work will be found to be useful in planning future experiments with the Penning trap platform that will employ additional anharmonic trap terms and a triangular rotating wall for a more uniformly spaced triangular lattice.

###### Acknowledgements.
We acknowledge useful conversations with Dan Dubin and John Bollinger. We acknowledge support from National Science Foundation under Grant No. PHY-1314295. A. K. acknowledges support from the Indo-US Science and Technology Forum via the S. N. Bose Scholars Program. J.K.F. further acknowledges support from the McDevitt bequest at Georgetown University. B.Y. acknowledges support from the Achievement Rewards for College Scientists Foundation.

## References

• (1) R. P. Feynman, Int. J. Theor. Phys. 21, 467 (1981).
• (2) J. I. Cirac and P. Zoller, Phys. Rev. Lett. 74, 4091 (1995).
• (3) D. Porras and J. I. Cirac, Phys. Rev. Lett. 92, 207901 (2004).
• (4) G.-D. Lin, S.-L. Zhu, R. Islam, K. Kim, M.-S. Chang, S. Korenblit, C. Monroe and L.-M. Duan, Europhys. Lett. 86, 60004 (2009).
• (5) Daniel H. E. Dubin, Phys. Rev. A 88, 013403 (2013).
• (6) M. J. Biercuk, H. Uys, A. P. Van Devender, N. Shiga, W. M. Itano, and J. J. Bollinger, Quantum Info. Comput. 9, 920 (2009).
• (7) J. W. Britton, B. C. Sawyer, A. C. Keith, C.-C. Joseph Wang, J.K. Freericks, M. J. Biercuk, H.Uys, and J. J. Bollinger, Nature (London) 484, 489 (2012).
• (8) M. McAneny, B. Yoshimura, and J. K. Freericks, Phys. Rev. A 88, 043434 (2013).
• (9) Brian C. Sawyer, Justin G. Bohnet, Joseph W. Britton, and John J. Bollinger, Phys. Rev. A 91, 011401R (2015).
• (10) D. A. Baiko, Phys. Rev. E 80, 046405 (2009).
• (11) J. D. Baltrusch, A. Negretti, J. M. Taylor, and T. Calarco, Phys. Rev. A 83, 042319 (2011).
• (12) C.-C. Joseph Wang, Adam C. Keith, and J. K. Freericks, Phys. Rev. A 87, 013422 (2013).
• (13) Shi-Liang Zhu, C. Monroe, and L.-M. Duan, Phys. Rev. Lett. 97, 050505 (2006).
• (14) C.-C. Joseph Wang and J. K. Freericks, Phys. Rev. A 86, 032329 (2012).
• (15) M. McAneny and J. K. Freericks, Phys. Rev. A 90, 053405 (2014).
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters