# Tuning heat transport in trapped-ion chains across a structural phase transition

###### Abstract

We explore heat transport across an ion Coulomb crystal beyond the harmonic regime by tuning it across the structural phase transition between the linear and zigzag configurations. This demonstrates that the control of the spatial ion distribution by varying the trapping frequencies renders ion Coulomb crystals an ideal test-bed to study heat transport properties in finite open system of tunable non-linearities.

###### pacs:

05.60.-k,64.60.Ht, 05.70.Fh, 37.10.TyUltracold ion Coulomb crystals represent one of the most promising platforms for the simulation of many-body physics thanks to the high degree of spatial and temporal control of mesoscopic ion crystals they afford us with HRB08 (); Schaetz12 (); BR12 (). Recent years have seen a shift away from the study of ground and thermal state properties, towards the exploration of the potential role of ion traps as a test-bed for models of non-equilibrium statistical mechanics.

In this context it is important to recognize that in addition to the electronic spin degrees of freedom trapped ions also possess motional degrees of freedom that can exhibit highly non-trivial static and dynamical properties including classical and quantum phase transitions. Indeed, ion Coulomb crystals confined in ion traps may support a wide variety of phases including a linear chain and a doubly-degenerate zigzag phase, extending further to increasingly complex configurations in two and three spatial dimensions Walther92a (); Walther92b (); Dubin99 (). The associated structural phase transitions between those configurations are generally of first order, with the exception of the linear-to-zigzag phase transition which is known to be of second order Schiffer93 (); Piacente04 (); Fishman08 (). As a symmetry breaking scenario, it provides a natural testing ground for universal dynamics of phase transitions and topological defect formation Fishman08 (); Retzker08 (); delCampo10 (); DeChiara10 (), recently explored in the laboratory Schaetz13 (); EH13 (); Ulm13 (); Tanja13 ().

Another fundamental setting in non-equilibrium statistical mechanics considers the thermal transport in low-dimensional systems, which exhibits a rich variety of anomalous features, including the breakdown of Fourier’s law of heat conduction, instances in which sub-diffusive and super-diffusive behavior can be observed LWH02 (); Alonso02 (); Alonso04 (), as well as the divergence of the thermal conductivity with the system size LLP03 (); Dhar08 (). Most rigorous theoretical results have been obtained for exactly solvable quasi-free models while systems with non-linearities are typically exceedingly difficult to treat. Equally, the controlled generation of non-linear physics in mesoscopic ion crystals is non-trivial and much recent progress has concerned harmonic models of complex networks and trapped-ion chains MP12 (); Bermudez13 (); Ramm13 (). The richest phenomenology however can be expected in non-integrable models HLZ () which mandates the development of both theoretical and experimental methods for their examination.

In this Letter, we advance further the case for trapped ions as a model system for the examination of challenging problems in mesoscopic physics by considering continuously driven ion chains between two thermal reservoirs as a platform in which to explore heat transport across an ion Coulomb crystal that experience a structural phase transition. As the ion chain crosses the transition from the linear to zigzag phase, the system may exhibit both significant non-linearities and axial-transverse mode coupling which can lead to qualitative changes in both the local temperature profile and the total heat flux through the chain that may be observed experimentally.

The system dynamics.- We consider an effectively -dimensional system composed of ions of mass , charge , positions and momenta , with , which are confined in a trap with axial frequency along the axis and transversal frequency along the axis. The Hamiltonian of the system can be written as

(1) |

The interaction potential accounts for both the harmonic trap and the Coulomb repulsion, and is given by

(2) |

A quasi-linear confinement of the ions along the axis can be achieved by considering a strongly anisotropic trap, with . We note that at variance with lattice systems, none of the ions is pinned, which allows for an intricate interplay between axial and radial modes of motion.

We assume that the dynamics due to the external Doppler cooling lasers acting on the ions can be modeled as Langevin thermostats. This together with the typical separations between the ions (generally of the order of microns) justifies an intrinsically noisy classical description of the dynamics,

(3) | |||

where and are the friction and diffusion coefficients, respectively, denote the Wiener processes resulting from the Gaussian white noise forces associated with the diffusion induced by the interaction with the laser beams, which satisfy and , and .

For small laser intensities the friction and diffusion coefficients can be obtained from the Doppler cooling expressions

(4) |

where is the normalized intensity of the laser beam acting on the -ion along the -direction, is the corresponding laser wavelength, is the detuning of the laser frequency with respect to the frequency of a selected atomic transition in the ions, and is the natural linewidth of the excited state in such transition Phillips92 ().

Heat flux and local kinetic temperature.- A discrete definition of the heat current through the chain can be obtained from the local energy density associated with each ion LLP03 (); Dhar08 (), which can be written as

(5) |

where represents the harmonic trap and the Coulomb term of the interaction potential given in Eq. (2). The time derivative of leads to the discrete continuity equations

(6) |

where

(7) |

can be identified as the energy current from the -ion to the -ion. For the -ion, the first term in Eq. (6) corresponds to the total energy current coming from the ions on the left, whereas the second term is the total energy current going to the ions on the right, see Fig. (1).

The last term

(8) |

is the energy current from the laser reservoirs.

The steady-state average of Eq. (6) implies the balance

(9) |

between the average rate at which each ion receives energy from the ions on the left and the laser beams, and the average rate at which such ion transfers energy to the ions on the right. The average of the energy currents from the reservoirs can be obtained using Novikov’s theorem Novikov (),

(10) |

The total heat current can be derived from a discrete description of the continuity equation LLP03 (); Dhar08 ()

(11) |

by taking the energy and heat flux densities as , respectively, with being the local energy density defined in Eq. (5) and the local flux. A Fourier analysis of Eq. (11) leads to Then the total heat flux, obtained by integration of the flux density over the chain volume, reads

(12) |

In the steady state, the averaged total heat flux is determined by just the local fluxes coming from the reservoirs, and applying Novikov’s theorem Novikov () it follows that

(13) | |||||

A discrete approach can also be considered to define a local temperature through the chain in terms of dynamical variables. According to the virial theorem, the local kinetic temperature of each ion can be defined from its kinetic energy

(14) |

where indicates the average over an ensemble of stochastic trajectories.

Numerical experiments.- We consider a chain composed of ions and analyze the response of the local temperature and the total heat flux to the phase transition from a quasi-linear to the planar zigzag spatial configuration as the transversal frequency of the trap is lowered. We consider a chain of Mg ions, with , and fix an axial frequency of the trap to kHz. We study the dynamics for different transversal frequencies .

An analysis of the static properties of the ion chain in the thermodynamic limit provides an estimate of the local value of the critical ratio of the trap frequencies leading to the phase transition between the linear and the zigzag configurations Fishman08 (), , where is the Riemann-zeta function, is the equilibrium linear density of ions along the trap axis as a function of distance from the chain’s center, and the half-length of the chain Dubin97 (). Due to the axial harmonic confinement the center of the chain experiences a higher axial density and Coulomb repulsion, making the phase transition spatially inhomogeneous.

In our numerical studies we assume that the ions are initially at rest and arranged with random positions in the close vicinity of the linear configuration. Then the reservoir lasers that act on the selected ions are switched on instantaneously. To determine the friction and the diffusion coefficients that characterize the interaction of the Mg ions with the laser beams, we have considered the Doppler cooling expressions (Tuning heat transport in trapped-ion chains across a structural phase transition) applied to the atomic transition with frequency THz NIST () and an excited state natural linewidth MHz Ansbacher89 (). Given these values, the interaction of a laser beam with an ion is a function of the normalized intensity and the detuning of the laser beams. Note that in order to avoid excessively small time steps the numerical model neglects micromotion. While this represents a significant approximation in the zigzag configuration of an ion Coulomb crystal in a rf-Paul trap, it should be noted, however, that micromotion is absent in Penning traps in which analogous structural phase transitions were recently observed and Doppler cooling can be implemented Segal13 ().

To drive a heat current through the chain, the ions at opposite ends of the ion crystal are subjected to different laser beams. In particular, we consider that the three leftmost (rightmost) ions interact with laser beams with normalized intensity () and detuning (), with (). The different detunings and lead to effective reservoirs operating at different temperatures on both ends of the chain, and therefore a stationary non-equilibrium heat current. We assume that no laser beams are acting on the inner ions, corresponding to . If local addressing of individual ions is not available, the use of Mg ions at the opposite ends of the chain that is otherwise composed of Mg can allow for highly localized heat baths by means of frequency selection. Mg ions that are not part of the zigzag structure will not affect the essential features of the structural phase transition. We are principally concerned with the steady state behaviour of the system. As a criterion to determine that the system has reached the steady state we verify expression (9) for each ion, and also apply equality (13) to the whole chain (see Supplementary Material for more details).

As Fig. (2) shows, the temperature profile adopts a gradient along the chain which progressively increases as the transverse frequency of the trap is reduced to drive the transition from the linear to the zigzag spatial distributions. The numerical simulations indicate that the phase transition first emerges for .

For the chain is fully linear and the temperatures of the inner ions tend to settle on a constant value for all , which for this system is close to the mean temperature that is expected in a homogeneous harmonic chain with nearest-neighbor interactions RLL67 (), despite the presence of the axial quadratic potential PC05 (). Indeed, a spectral analysis of the steady evolution of the coordinates of the inner ions in linear chains indicates that their axial dynamics is close to the Brownian motion of a simple harmonic oscillator with characteristic frequency . Figure (3) illustrates that the spectra of these ions are dominated by a well defined single peak corresponding to the axial trap frequency. Also an ordered series of much less intense peaks corresponding to higher order multiples of can be distinguished. The very low intensity of the transversal spectrum evidences the minor role of this mode, due to the strong trap confinement in radial direction.

The simple axial spectrum corroborates the expectation that a harmonic approximation to the system Hamiltonian is valid; the system is effectively integrable, and the heat carriers are freely propagating phonons. The lack of temperature gradient observed in linear chains is characteristic of this ballistic behavior.

This situation is expected to change when the chain approaches the structural phase transition at where the chain buckles with the growth of the zig-zag soft mode. Near this point non-linearities as well as mode coupling between axial and radial modes are expected to play an increasing role and the harmonic chain description is expected to fail. While the non-linearities lead to scattering, the coupling between axial and radial modes lead to an effective dynamics akin to dephasing noise. Both effects, if significant, are known to contribute to the formation of a temperature gradient in the chain. It should be noted that non-linearities tend to be relatively small unless the chain is very close to the phase transition point Marquet02 (). Hence we expect coupling between radial and axial modes to dominate. The deviations from the harmonic picture are also witnessed by the spectra of the axial motion of the ions in the zig-zag configuration, which exhibit more complex features that were absent in the linear chains, see Figure (3), and that can be assigned to the coupling with the transverse motion. The presence of bands of irregular patterns extend across the axial and transverse spectra suggests underlying domains of chaotic dynamics arising from such coupling. As the harmonic chain description ceases to be valid, the heat transport through the chain is modified, as evidenced by the emergence of a temperature gradient, signalling a possible diffusive behavior of the heat carriers.

Further to the emergence of Fourier’s law in the ion chain when approaching the structural phase transition, one may also observe clear signatures of the structural phase transition in the heat flux across the chain as demonstrated by figure (4).

Initially, on approaching the structural phase transition from the linear chain, one observes a small increase of the heat flux. This has two origins. Firstly, in the proximity of the transition the transverse modes will start to contribute to transport. Secondly, the increased thermal motion of the ions as the chain softens upon approach of the phase transition lead to an increased level of fluctuations. This noise may assist transport as it overcome the effects of spatial inhomogeneity in the chain PlenioH08 (). Upon further decrease of the chain buckles, leading to two possible heat conduction paths, while the inter-ion distances increase which in turn leads to a reduction in the interaction between neighboring ions and therefore a reduction of the heat flux. While initially, close to the phase transition, the increase in distance is compensated for by the formation of two independent channels, deeper in the zigzag configuration this is not the case anymore and the heat flux reduces.

Before closing we point out that a thermal conductance could be estimated from the temperature gradient and the axial heat flux measured in a non-equilibrium steady state. The spatial constrains imposed by the finite size and low-dimensionality of trapped-ion chains prompt the analysis of the thermal conductance as a function of the length, instead of a size-independent thermal conductivity of interest in a macroscopic model of thermal conduction. An alternative approach to obtain the heat conductance of a specific ion chain could be based on Green-Kubo type linear response expressions valid for the heat current in finite low-dimensional systems coupled to heat reservoirs Dhar08 (); Kundu09 ().

Conclusions.- Our analysis based on the local temperature profile and the total heat flux indicates that trapped-ion chains exhibit anomalous heat transport. The linearly distributed ions resemble harmonic chains, and therefore an integrable system, in which the free energy transport along the chain by non-interacting axial modes precludes the establishment of a temperature gradient and would lead to a divergent thermal conductivity. The phase transition from the linear to the bi-dimensional zigzag configuration induces a coupling between axial and transverse modes that hinders the transport of energy along the chain, and allows the emergence of a central domain in which a temperature gradient can be set up. Such domain grows as the transversal frequency is lowered and the bidimensional configuration extends towards the ends of the chain, resulting in a significant decrease of the axial heat flux. Heat transport is optimal in the linear configuration in the proximity of the critical point.

Note: After the completion of this work, we learned about reference Freitas13 () devoted to the study of quantum heat conduction in harmonic ion chains.

Acknowledgment.- It is a pleasure to acknowledge discussions with T. E. Mehlstäubler, J. M. Plata and D. Roy. This project was funded by the Spanish MICINN, the European Union (FEDER) (FIS2010-19998), the EU Integrating project SIQS, the EU STREP EQUAM, the Alexander von Humboldt Professorship (MBP) and the U.S. Department of Energy through the LANL/LDRD Program and a LANL J. Robert Oppenheimer Fellowship (AdC).

## References

- (1) H. Häffner, C. F. Roos, and R. Blatt, Phys. Rep. 469, 155 (2008).
- (2) C. Schneider, D. Porras, and T. Schaetz, Rep. Prog. Phys. 75, 024401 (2012).
- (3) R. Blatt and C. F. Roos, Nature Phys. 8, 277 (2012).
- (4) G. Birkl, S. Kassner, and H. Walther, Nature 357, 310 (1992).
- (5) I. Waki, S. Kassner, G. Birkl, and H. Walther, Phys. Rev. Lett. 68, 2007 (1992).
- (6) D. H. E. Dubin and T. M. O’Neil, Rev. Mod. Phys. 71, 87 (1999).
- (7) J. P. Schiffer, Phys. Rev. Lett. 70, 818 (1993).
- (8) G. Piacente, I. V. Schweigert, J. J. Betouras, and F. M. Peeters, Phys. Rev. B 69, 045324 (2004).
- (9) S. Fishman, G. De Chiara, T. Calarco, and G. Morigi, Phys. Rev. B 77, 064111 (2008).
- (10) A. Retzker, R. C. Thompson, D. M. Segal, and M. B. Plenio, Phys. Rev. Lett. 101, 260504 (2008).
- (11) A. del Campo, G. De Chiara, G. Morigi, M. B. Plenio, and A. Retzker, Phys. Rev. Lett. 105, 075701 (2010).
- (12) G. De Chiara, A. del Campo, G. Morigi, M. B. Plenio, and A. Retzker, New J. Phys. 12, 115003 (2010).
- (13) M. Mielenz, J. Brox, S. Kahra, G. Leschhorn, M. Albert, T. Schaetz, H. Landa, and B. Reznik, Phys. Rev. Lett. 110, 133004 (2013).
- (14) S. Ejtemaee and P. C. Haljan, Phys. Rev. A 87, 051401(R) (2013).
- (15) S. Ulm S, J. Roßnagel, G. Jacob, C. Degünther, S. T. Dawkins, U. G. Poschinger, R. Nigmatullin, A. Retzker, M. B. Plenio, F. Schmidt-Kaler, and K. Singer, Nat. Commun. 4, 2290 (2013).
- (16) K. Pyka , J. Keller, H. L. Partner, R. Nigmatullin, T. Burgermeister, D. M. Meier, K. Kuhlmann, A. Retzker, M. B. Plenio, W. H. Zurek, A. del Campo, and T. E. Mehlstäubler, Nat. Commun. 4, 2291 (2013).
- (17) B. Li, L. Wang, and B. Hu, Phys. Rev. Lett. 88 223901 (2002)
- (18) D. Alonso, A. Ruiz, and I. de Vega, Phys. Rev. E 66, 066131 (2002).
- (19) D. Alonso, A. Ruiz, and I. de Vega, Physica D 187, 184 (2004).
- (20) S. Lepri, R. Livi, and A. Politi, Phys. Rep. 337, 1 (2003).
- (21) A. Dhar, Adv. Phys. 57, 457 (2008).
- (22) E. A. Martinez and J. P. Paz, Phys. Rev. Lett. 110, 130406 (2013).
- (23) A. Bermudez, M. Bruderer, and M. B. Plenio, Phys. Rev. Lett. 111, 040601 (2013).
- (24) M. Ramm, T. Pruttivarasin, and H. Häffner, arXiv:1312.5786.
- (25) B. Hu, B. Li, and H. Zhao, Phys. Rev. E 57, 2992 (1998); 61, 3828 (2000).
- (26) See for example W. D. Phillips, Laser cooling and trapping of neutral atoms, in Laser Manipulation of Atoms and Ions, edited by E. Arimondo, W. D. Phillips, and F. Strumia, Proceedings of the International School of Physics “Enrico Fermi”, Course CXVIII (North-Holland, Amsterdam, 1992), p. 289.
- (27) E. A. Novikov, Soviet Phys. JETP 20, 1290 (1965).
- (28) D. H. E. Dubin, Phys. Rev. E, 55, 4017 (1997).
- (29) The frequency was obtained from the energy transition eV taken from the NIST Atomic Spectra Database 1999. http://physics.nist.gov.
- (30) W. Ansbacher, Y. Li, and E. H. Pinnington, Phys. Lett. A 139, 165 (1989).
- (31) S. Mavadia, J. F. Goddwin, G. Stutter, S. Bharadia, D. R. Crick, D. M. Segal, and R. C. Thompson, Nat. Comm. 4, 2571 (2013).
- (32) Z. Rieder, J. L. Lebowitz, and E. Lieb, J. Math. Phys. 8, 1073 (1967).
- (33) T. Prosen and D. K. Campbell, Chaos 15, 015117 (2005). See Subsection II A. 2.
- (34) C. Marquet, F. Schmidt-Kaler, and D. F. V. James, Appl. Phys. B 76, 199 (2003).
- (35) A. Kundu, A. Dhar, and O. Narayan, J. Stat. Mech. L03001 (2009)
- (36) M. B. Plenio and S. F. Huelga, New J. Phys. 10, 113019 (2008).
- (37) N. Freitas, E. Martinez, and J. P. Paz, arXiv:1312.6644.

## Appendix A Dimensionless variables

The analysis of the system dynamics can be simplified considering the equations of motion (3) in terms of dimensionless variables. To define such variables we start by introducing a characteristic system length , given by the relation

(15) |

Defining the dimensionless ion coordinate and momentum vectors as

(16) |

and

(17) |

and the one-parameter dimensionless interaction potential

(18) |

where is the aspect ratio of the trap frequencies, the equations of motions take the form

(19) | |||

with the dimensionless time

(20) |

the Wiener processes

(21) |

the friction coefficients

(22) |

and the diffusion coefficients

(23) |

In terms of the dimensionless variables, the total heat flux (12) is expressed as

(24) |

with the dimensionless local energy densities and currents,

(25) |

and

(26) |

respectively. The dimensionless local temperatures are given by

(27) |

In the next section we will continue the analysis in terms of the dimensionless variables. We will remove the tilde symbol from all the variables and parameters to simplify the notation.

## Appendix B The numerical integration

The dimensional stochastic differential equations (19) can be expressed in a compact form as

(28) |

where the components of the variable vector have been ordered in the form

(29) |

The components of the vector containing the deterministic terms in the equations of motion are

The matrix contains the diffusion coefficients . In our model it can be expressed by a diagonal matrix with the elements

(30) |

The vector denotes the dimensional Wiener process, whose elements are

(31) |

To integrate the stochastic differential equations (28) we consider the multi-dimensional explicit order 2.0 weak scheme proposed by Platen Kloeden (). Since in our model the matrix does not depend explicitly on the variable , such scheme is particularly simple. Given the variable at a time step , its value at the following time step is given by

(32) |

where

(33) |

is the constant time interval between two consecutive time steps, and is the vector with elements

(34) |

with (standard Gaussian) a normally distributed random variable selected at time step for the -ion, along the -direction.

In order to get stable local temperatures and a total heat flux from the solutions of the equations (28) numerical simulations were performed considering time intervals , and integrated up to a final time in which the steady conditions (9) and (13) were satisfied. The process became computationally expensive as the multiple interaction potential (18) had to be evaluated more than times for each stochastic, and averages that included over of these trajectories were considered. It took over two months to carry out the simulations and get stable results using a CPU machine with AMD Opteron (tm) Processors 6134.

## Appendix C Characteristic timescales in the non-equilibrium dynamics towards the steady-state

So far we have presented results obtained in terms of dynamical variables collected from long enough simulations for the system to reach a non-equilibrium steady state. In this section we indicate the characteristic timescales that are needed to achieve such state, with associated time-independent local temperatures and heat fluxes.

Figures (5) and (6) show some representative time evolutions of the local temperatures and the total heat fluxes towards the steady state configuration, obtained from an average over more than stochastic trajectories. The ions at both ends of the chain that are directly connected to the laser beams reach steady values faster than the inner ions. The slowest convergence occurs for the central ions of the linear chain. In general, it can be assumed that the system stabilizes after approximately ms. To ensure accurate temperature profiles and heat fluxes, the simulations have been performed up to ms, and the final results have been obtained from a time average within the interval ms.

Hence we are considering short-time-scale experiments in which effects such the motional heating of the trapped ions confined in rf-traps due to fluctuating electric fields from the trap electrodes should not be relevant.

## References

- (38) P.E. Kloeden and E. Platen, Numerical Solution of Stochastic Differential Equations (Springer, 1999).