Spontaneous nucleation of structural defects in inhomogeneous ion chains

Structural defects in ion crystals can be formed during a linear quench of the transverse trapping frequency across the mechanical instability from a linear chain to the zigzag structure. The density of defects after the sweep can be conveniently described by the Kibble-Zurek mechanism. In particular, the number of kinks in the zigzag ordering can be derived from a time-dependent Ginzburg-Landau equation for the order parameter, here the zigzag transverse size, under the assumption that the ions are continuously laser cooled. In a linear Paul trap the transition becomes inhomogeneous, being the charge density larger in the center and more rarefied at the edges. During the linear quench the mechanical instability is first crossed in the center of the chain, and a front, at which the mechanical instability is crossed during the quench, is identified which propagates along the chain from the center to the edges. If the velocity of this front is smaller than the sound velocity, the dynamics becomes adiabatic even in the thermodynamic limit and no defect is produced. Otherwise, the nucleation of kinks is reduced with respect to the case in which the charges are homogeneously distributed, leading to a new scaling of the density of kinks with the quenching rate. The analytical predictions are verified numerically by integrating the Langevin equations of motion of the ions, in presence of a time-dependent transverse confinement. We argue that the non-equilibrium dynamics of an ion chain in a Paul trap constitutes an ideal scenario to test the inhomogeneous extension of the Kibble-Zurek mechanism, which lacks experimental evidence to date.

Abstract. Structural defects in ion crystals can be formed during a linear quench of the transverse trapping frequency across the mechanical instability from a linear chain to the zigzag structure. The density of defects after the sweep can be conveniently described by the Kibble-Zurek mechanism. In particular, the number of kinks in the zigzag ordering can be derived from a time-dependent Ginzburg-Landau equation for the order parameter, here the zigzag transverse size, under the assumption that the ions are continuously laser cooled. In a linear Paul trap the transition becomes inhomogeneous, being the charge density larger in the center and more rarefied at the edges. During the linear quench the mechanical instability is first crossed in the center of the chain, and a front, at which the mechanical instability is crossed during the quench, is identified which propagates along the chain from the center to the edges. If the velocity of this front is smaller than the sound velocity, the dynamics becomes adiabatic even in the thermodynamic limit and no defect is produced. Otherwise, the nucleation of kinks is reduced with respect to the case in which the charges are homogeneously distributed, leading to a new scaling of the density of kinks with the quenching rate. The analytical predictions are verified numerically by integrating the Langevin equations of motion of the ions, in presence of a time-dependent transverse confinement. We argue that the non-equilibrium dynamics of an ion chain in a Paul trap constitutes an ideal scenario to test the inhomogeneous extension of the Kibble-Zurek mechanism, which lacks experimental evidence to date. The high-symmetry phase corresponds to a linear chain, while the broken-symmetry phase is characterised by a doubly degenerate zig-zag configuration.

Introduction
Ion crystals in Paul or Penning traps represent a prominent example of a self organized system that is amenable to accurate experimental characterization and manipulation [1]. Ions confined by means of static and radio-frequency electromagnetic potentials reach crystallization when laser cooled. Crystals made up from tens to millions of ions have been observed both in Paul traps [2,3,4,5] and in Penning traps [6]. Different structures can be realized by varying the particle density or the trap anisotropy. In Ref. [2] quasi one-dimensional structures have been experimentally characterized for first time. Here, the first three structures encountered when decreasing the transverse confinement are a linear chain, a planar zigzag structure, and a helicoidal arrangement of the ions about the trap axis. Phase transitions separating two structures are usually discontinuous. In this respect, the transition from the linear chain to the zigzag structure constitutes an exception. In fact, it has been demonstrated numerically [7,8] and later analytically using Landau theory [9] that the transition from a linear chain to the zigzag is of second order. The transition can be induced either by increasing the ion density or by decreasing the transverse confinement ν t so that their value exceeds or is below, respectively, a certain critical value determining the mechanical instability. Recently, the transition linear-zigzag chain has been suggested as a test-bed for manybody quantum effects and the creation of double-well potentials [11], and is well suited for the distinction between the nucleation of defects and quasiparticles in a symmetry breaking scenario [12].
The aim of this paper is to study the out-of-equilibrium dynamics of a linear chain of trapped ions when ν t is lowered in time from above to below the mechanical instability separating the linear from the two-dimensional zigzag configuration, illustrated in Fig. 1. Due to the finiteness of the sound velocity, space-like separated regions may develop one of the two possible zigzag orderings: odd ions up and even Figure 2. Ion chain in a zigzag configuration exhibiting structural defects. The ions are confined in a ring trap, and we report here the distribution of charges in two dimensions along the ring mapped to a line for clarity. The solid line joining the ions serves as a guide to the eye. In the example reported there are 4 defects. According to the Kibble-Zurek mechanism, the average size of the domains is given by the correlation length at the freeze-out time,ξ (see text for details).
ions down or odd ions down and even ions up as depicted in Fig. 2. These regions are the analogs of magnetic domains in a ferromagnetic material and the interface between these domains is a structural defect. The classical and quantum properties of these defects and their possible use for quantum information processing were studied in [13]. Defects formation can be understood from simple statistical mechanics considerations. If the rate of change of ν t is larger than the relaxation rate of the crystal, the latter does not have enough time to relax to the minimum energy configuration, which is characterized by a perfectly ordered zigzag structure, thus exhibiting no defects and corresponding to the state that would have been obtained after a perfect adiabatic transition. In principle, for an infinite system whose phase transition is described by Landau theory [14,15], the dynamical relaxation time scale diverges at the critical point and therefore no matter how slowly ν t changes, there will always be proliferation of defects. This is the celebrated Kibble-Zurek mechanism (KZM) which can generally account for the nucleation of topological defects in scenarios of spontaneous symmetry breaking [16,17]. Moreover, Zurek predicted the scaling of the number of defects as a function of the rate of passing through the phase transition [17]. The KZM scenario and Zurek's prediction for the scaling of the number of defects have been verified in a variety of systems numerically [18] and experimentally [19]. Recently, the KZ scaling prediction has also been extended to quantum phase transitions [20].
In this work we extend the results in our proposal [23] and study the production of kinks in the linear to zigzag transition in ion crystals by deriving a time dependent Ginzburg-Landau theory for the transition. With this result, we determine the scaling of the number of defects with the rate of change of the transverse frequency which follows from the KZM and compare it with numerical simulations. In a ring trap when the interparticle distance is uniform, one recovers the standard KZM. By contrast, in a linear Paul trap the density of ions is inhomogeneous, higher at the center and more rarefied at the edges. When decreasing in time the transverse trap frequency, the value of the mechanical instability is first crossed in the center of the chain, and one can identify a front crossing the transition from the linear to the zigzag phase which moves from the center towards the edges, and whose velocity depends on the rate of the quench. If the velocity of this front is smaller than the sound velocity, at which a perturbation propagates along the chain, the dynamics becomes adiabatic. This result holds even in the thermodynamic limit. Otherwise, the nucleation of defects is partially suppressed with respect to the homogeneous case and a novel scaling of the density of kinks with the rate of the transition emerges in the inhomogeneous extension of KZM [21,22]. The whole body of experimental results in [19] has been aimed at the verification of the paradigmatic KZM in a homogeneous scenario (HKZM). By contrast, the inhomogeneous KZM (IKZM) [21,22] lacks to date experimental verification, and we argue in the following that ion chains in a Paul trap constitute an ideal system for such goal.

Ginzburg-Landau equation in presence of laser cooling for the order parameter
In this section we start from the theory presented in [24,9], and derive a Ginzburg-Landau equation for the order parameter of the linear-zigzag structural transition. The order parameter corresponds to the position offset of the ions from the trap axis. For a large number of ions, using a local density approximation we can approximate the crystal as a continuum, so that the order parameter is a field. The theory is extended to the case in which the crystal motion is laser cooled, and the coupling to an external reservoir is described using the theoretical model developed in [25]. Defect formation by quenching the control parameter, here the transverse trap frequency, across the structural instability is studied by using the corresponding Langevin equation for the order parameter. This theoretical model allows us to estimate the density of defects at the end of the quench.

Ginzburg-Landau Lagrangian for the structural phase transition
The system we consider consists of N ions of mass m, charge Q and coordinates r n = (x n , y n , z n ) which are confined in a quasi one dimensional trap with tight harmonic confinement at frequency ν t in the plane yz. In this section we assume a ring trap with large radius R, so that we may impose periodic boundary conditions along the x-axis. The Lagrangian describing the dynamics of the ions is where the kinetic and potential energies take the form, respectively, At sufficiently low temperatures and sufficiently large transverse confinement the ions crystallize around the stable equilibrium points of the potential V , which are aligned along the x-axis, r 0n = (na, 0, 0), with a the interparticle distance in the x-direction such that a = 2πR/N . The stability of the linear chain along the x axis requires a transverse trap frequency exceeding a threshold value ν (c) t , which scales with the characteristic frequency ω 0 = Q 2 /ma 3 . At ν (c) t the configuration has a structural instability, such that for ν t < ν (c) t the ions are organized in a planar structure with equilibrium positions r 0n = (na, (−1) n (b/2) cos θ, (−1) n (b/2) sin θ) with θ ∈ [0; 2π] the angle between the crystal plane and the plane xy. The structure has the form of a zigzag line, which joins the charges along the plane, with b the transverse size of the zigzag. An appropriate thermodynamic limit can be defined, letting the number of ions N → ∞ while keeping the value ν (c) t and the characteristic frequency ω 0 fixed. This corresponds to keeping fixed the linear interparticle distance a or equivalently the charge density [24]. One finds that the structural instability is a second-order phase transition, with control field ν t (or alternatively, a) and order parameter b [9]. In particular, the critical value of the transverse frequency in the thermodynamic limit is given by ν (c) t = ω 0 7ζ(3)/2 = (2.051 . . .)ω 0 , ζ being the Riemann-zeta function. From the Landau theory of the structural phase transition in Ref. [9], one can develop a continuum model for the field describing the order parameter applying standard assumptions, thereby obtaining the Ginzburg-Landau equation. We sketch below the assumptions made, in order to clarify the range of validity of the model which forms the basis of this work.
In the linear chain, in the harmonic limit, transverse and axial modes are decoupled. For periodic boundary conditions the modes have a well defined quasimomentum k which takes values within the Brillouin zone. For N ions, the transverse, normal modes at wave vector k are the real and imaginary part of the mode where σ n = y n , z n is the transverse displacement of the ion at the equilibrium position r 0n and k ∈ [0, π/a]. Using decomposition (4), as shown in Ref. [9], where all the details of the derivation are provided, the transverse potential, obtained by expanding Eq. (3) around the equilibrium positions, takes the form where the label indicates the order of the expansion. In particular, and where we omitted to write the potential term for the axial modes. The terms V (3) and V (4) correspond to third-and fourth-order expansion in the fluctuations around the equilibrium position, and contain the coupling between the axial and radial modes. We now make the assumption that the ions are pinned in the axial direction and can only oscillate in the transverse direction and discard the coupling to the axial modes. In this regime, it was found in Ref. [9] that V (3) = 0 and the fourth order term reads: In the stability regime one finds that β(k) is minimum at k 0 = π/a, corresponding to the so-called zigzag mode. The critical value of the transverse trap frequency ν (c) t is found from the relation β(k 0 ) = 0. When ν t < ν (c) t the linear chain is unstable, and the ions order in a zigzag structure across the trap axis.
Hence, sufficiently close to the critical value ν (c) t , the mode of the linear chain which first becomes unstable and determines the equilibrium positions of the zigzag structure (soft mode) is the transverse mode of the linear chain with the shortest wave length (zigzag mode), at wave vector k 0 = π/a. When ν t is sufficiently close to the critical value ν (c) t , an effective potential can be derived for the transverse normal modes Ψ σ k with wave vectork = k 0 − δk, such that aδk ≪ 1. The effective potential is composed by a quadratic component, given in Eq. (6), with coefficient β(k), Eq. (7) such that [10] and which can take negative values depending on the value of ν t . Within the same assumptions, we can approximate: where A = (93ζ(5)/32)ω 2 0 /a 2 and we obtain for the fourth order term: Note that the prime index in the sum over the quasimomenta refers to the fact that the sum is restricted to the quasi momentak j = k 0 ± δk, with δka ≪ 1. The potential form V in Eq. (5) is valid at second order in the small parameter δka ≪ 1, when the amplitude of the transverse oscillations is much smaller than the interparticle distance a.
The limit aδk ≪ 1 corresponds to a long-wavelength expansion made with respect to the zigzag mode: the phase difference of the considered transverse modes from the zigzag modes varies slowly from ion to ion, and a continuum description can be introduced. In this limit x n → x, where x is a continuous variable, and (−1) n σ n → ψ σ (x), where ψ σ (x) is now the position-dependent order parameter, which is related to ψ σ k through the Fourier transform while the factor 1/a in the right-hand side gives the density of states. For a large number of ions the discrete sum overk vectors in the potential energy becomes an integral according to the rule k → d(δk)N a/(2π). The second order potential term becomes: (14) and using the integral representation of the Dirac function: we obtain: A similar calculations leads to: Finally, we obtain the Lagrangian L = dxL(x) where L(x) is the Lagrangian density and reads We now discuss individually the parameters entering Eq. (18). The parameter determines whether the ground state of the system is a linear chain or a zigzag, depending on its sign. The parameter h = ω 0 a √ log 2 is a velocity, and determines the speed with which a transverse perturbation propagates along the chain. Finally, the parameter A is positive and determines the value of the order parameter when δ < 0 as will be shown below. The Lagrangian density we derived has the form of a Ginzburg-Landau equation. It is valid for the modes of the linear chain close to instability and extends the theory presented in [9] for the specific case in which the coupling between axial and transverse modes can be discarded. The minimal energy solution of Eq. (18) is constant in space and fulfills the relation: Equation (20) always admits the solution ψ σ = 0, corresponding to all the ions laying on the x axis. This solution is clearly stable only for δ > 0, i.e. in the linear chain phase. For δ < 0 there is an infinite family of solutions satisfying ̺ = ± −δ/2A, with ̺ = ψ y2 + ψ z2 [9]. These solutions correspond to zigzag structures with the same amplitude but laid on different planes. From the Ginzburg-Landau equation one can derive the critical exponent of the correlation length close to the critical point and at T = 0. In particular, under the assumption of a small spatial and static deformation of the order parameter at a certain position, the field autocorrelation function at distance x decays exponentially as ∼ exp(−|x|/ξ) [14], where ξ = h/δ 1/2 is the correlation length, such that and it diverges as ξ ∼ (ν t − ν (c) t ) −1/2 with the corresponding critical exponent is 1/2. The same exponent governs the behaviour of the correlation length in the zigzag phase close to the critical point.

Equation of motion of the slow modes in presence of laser cooling
We now consider the physical situation, in which the crystal motion is laser cooled. For the moment, we consider the discrete distribution of ion charges forming a linear chain (so that the effective potential for the transverse modes is essentially quadratic). For Doppler cooling, an effective equation for the crystal modes can be derived, which has the form of a Fokker-Planck equation for the energy distribution of the normal modes [25]. This can be put in terms of a Langevin equation of the form which is here reported for the transverse modes. The damping rate η σ,k depends on the corresponding mode frequency and it can take different values depending on the propagation direction of the cooling lasers. The scalar ε σ k (t) represents the corresponding Langevin force, such that its moments fulfill the relations where κ B is the Boltzmann constant and T is the temperature which determines the thermal state due to cooling. Equation (22) relies on the assumption that the dynamics of the electronic degrees of freedom is much faster than that of the motional degrees of freedom so they can be adiabatically eliminated from the equations for the degrees of freedom of the crystal normal modes, which is generally true when all ions are driven by the cooling laser [26].
In the following we focus on the situation in which the trap frequency along the z-axis is much larger than that along the y-axis. In this limit we drop the σ-label, and write an equation for the transverse modes along y and close to the instability point, at the transition from a linear chain to a zigzag structure in the xy plane. Here, the parameters η k and ε k are slowly varying and can be assumed to be constant. Under these assumptions, the Euler-Lagrange equation for ψ = ψ y (x) obtained from Eq. (18) reads Note that h can be obtained from the group velocity s = h 2 k √ δ+h 2 k 2 at the critical point δ = 0. Damping introduces a characteristic relaxation time scale τ , with which the system reaches equilibrium. Sufficiently close to the transition point, η 2 ≫ δ and see [14] and the appendix. Therefore the relaxation time τ diverges when the transverse frequency ν t approaches the critical value ν The same scaling is valid also in the zigzag phase δ < 0. The quantities ξ and τ and their scaling with δ are crucial for determining the scaling of the defect production during a quench of the transverse frequency as explained in the next section.

Quenching the trap frequency through the critical value
Within the Ginzburg-Landau description we now assume that the transverse trap frequency ν t undergoes a change in time in the interval [−τ Q , τ Q ], such that and with δ 0 > 0 and δ 0 ≪ ν (c) t . The transverse trap frequency value is swept through the mechanical instability of the linear-zigzag chain, such that δ(−τ Q ) = δ 0 > 0 and δ(τ Q ) = −δ 0 < 0. Correspondingly, the equation for the field now reads where for a given time t the ground state and thermodynamic properties at the phase transition are well defined. The time-dependence of δ(t) gives now a non-equilibrium problem. In principle, this should be solved considering the differential equation with time-dependent parameters. Nevertheless one might reach an understanding of the essential features by separating it into two domains. In the first domain, one can use the equilibrium solution provided that one may assume the adiabatic approximation, namely, |δ(t)/δ(t)| ≫ τ (t), where τ (t) is the relaxation time for the equilibrium situation at the given value δ = δ(t). The second domain is known as the impulse region, where as a result of the critical slowing down the order paramer is assumed to be frozen, ceasing to react to the external quench [17]. This stage is separated from the adiabatic regime by the freeze-out time scalet, i.e., a time scale in which the adiabatic condition ceases to be valid. This time scale can be estimated by setting |δ(t)/δ(t)| = τ (t), which clearly sets a lower limit tot. On this time-scale the order parameter is correlated over domains of characteristic lengthξ = ξ(δ(t)) given by the correlation length at the freeze-out time. The density of defects d (number of defects over the total number of ions) can then be estimated by the relation d ∼ 1/ξ. We now determine the density of defects for the chosen time variation of the parameter δ, Eq. (28). In this case, an estimate for the freeze-out timet is given by the instant of timet at which the rate of change of δ equals the relaxation time: Equation (30) allows for simple solutions in two specific limits, which will be analyzed in this paper. In the limit, in which δ(t) ≪ η or equivalently η 3 ≫ δ 0 /τ Q , the damping overcomes the oscillations associated with the frequency at the freeze-out time scale, δ(t). Note that this condition is imposed at the freeze-out-time before crossing the transition. It is then when the properties of the broken symmetry phase are determined. We denote this regime by "overdamped limit", following the definition introduced in Ref. [28]. In this case, using Eqs. (28) and (26) one findst = √ τ 0 τ Q with τ 0 = η/δ 0 , which sets δ(t) = ηδ 0 /τ Q . The density of defects in the overdamped limit takes the form where we used Eq. (21), assuming that the temperature is sufficiently low to neglect finite temperature effects. Another limit, which we will consider, is the one in which δ(t) ≫ η, namely, the Ginzburg-Landau equation is the one of an underdamped oscillator at the freeze-out time. We denote this regime by "underdamped" limit according to Ref. [28]. In this case,t = (τ 2 0 τ Q ) 1/3 with τ 0 = 1/ √ δ 0 , δ(t) = (δ 0 /τ Q ) 2/3 , and the density of defects is given by The algebraic scaling in Eq. (32) coincides with the one found in Ref. [28]. The different power scaling of δ 0 , with respect to Ref. [28], is due to the fact that in the present expressions we have written explicitly the proportionality factors (which in [28] are summarized in the parameter ξ 0 ). We finally note that, the fact that the , which in our case corresponds to a temperature of the order of 20 phonons. Note that the ratio of the number of defects over the total number of ions N is related to the spatial density of defects by a factor l/N . density of defects exhibits a power-law behaviour as a function of τ Q is obviously due to the specific choice of the quench as a function of time and of the considered regime (determined by the ratio η/δ(t)). The result can be understood as the domain size, determining the density of defects is given by the speed at which a perturbation propagates along the chain, aω 0 , multiplied by the relaxation time scale.
The scenario just discussed describes the dynamics of the structural phase transition for ions placed on a ring with equilibrium positions r 0n , where the interparticle spacing is homogeneous. We now verify numerically the HKZM prediction by integrating numerically the Euler-Lagrange equations obtained by minimizing the Lagrangian in Eq. (1) for a finite number of ions, and whose motion is damped in presence of a Langevin force. Assuming a ring confinement and pinning of one ion, the relevant degrees of freedom in the critical region are the transverse coordinates y n . At t = 0 the ions are in a linear chain configuration with y n = 0. The system is then driven by a linear quench of the form given in Eq. (28), and the density of defects d is computed at some asymptotic time, after which d remains practically constant. A typical evolution instance with 4 defects is shown in Fig. 2. These defects resemble the non-massive kinks of the Frenkel-Kontorova model with a transversal degree of freedom that can be described by an effective φ 4 theory for the translational displacement [29]. The classical and quantum behavior of these kinks was studied in [13].
In the numerics, for each set of parameters, the density of defects is averaged over many different realizations such as the one in Fig. 2. We then study the dependence of d on τ Q . A least-squares fit to the list of data is used to compute the scaling exponent. Figure 3 shows the scaling of the density of defects d as a function of the frequency quench time τ Q in the overdamped regime, in agreement with the homogeneous KZM. Deviations from KZM occur for fast quenches, for which the density of defects saturates as a result of the interactions between different kinks. In addition, finite size effects limits also the validity of the KZM scaling for very slow quenches that can lead to an adiabatic dynamics. Nucleation of defects is expected to be suppressed whenever the correlation length at the freeze-out time scale exceeds the size of the system,ξ > N a, both in the underdamped and overdamped regimes.

Ion chain inside a linear Paul trap: inhomogeneous effects
The standard scaling of topological defects for homogeneous phase transitions should be revised whenever the quench is local or there is a spatial dependence of the relative frequency δ [22,21]. In this section we discuss the dynamics of the structural phase transition for ion chains with open boundaries, for which both transverse and axial trapping potentials are harmonic. The ion chain in a linear Paul trap is characterized by a linear charge distribution which is inhomogeneous, with ion density increasing towards the center of the trap [30]. This leads to a stronger Coulomb repulsion for the ions near the center of the chain. Correspondingly, the radial short-wavelength modes have amplitude which is larger at the center, and vanishes at the edges [24]. The structural instability is hence first visible at the center of the trap, such that at the critical value of the transverse frequency the central ions are displaced from the trap axis forming a zigzag chain, as illustrated in Fig. 4. At lower values of the transverse frequency the number of ions which are displaced from the center increases towards the edges of the chain, until all the chain is in a planar configuration [7]. For a sufficiently long chain, one could associate with this behaviour a spatially-dependent critical frequency, determining the transition to the zigzag, such that it is largest at the center and smaller at the edges. Correspondingly, if a quench of the transverse frequency is applied, the transition point is crossed at different instant of times along the chain (from the center to the edges). The velocity, at which this front propagates, is the so-called front velocity [22,23]. In this case, the ratio between the front velocity and the sound velocity determines nucleation of defects. In the following, we compute the number of kinks formed in an ion crystal in a linear Paul trap in such scenario. The mechanism resembles the formation of solitons in a cigar-shaped Bose-Einstein condensate recently discussed by Zurek [21].
More precisely, we now consider the case, in which the ions are also trapped in the x-direction by a harmonic potential of frequency ν, such that the aspect ratio ν/ν t is sufficiently small to allow for low-dimensional crystalline structures. We will now consider the case in which the number N of particles is finite, but it is sufficiently large to allow for a local density approximation. In this limit and away from the chain edges the linear density n(x) is well approximated by the function [30] n(x) = 3 4 with L the half-length of the chain and x the distance from the center. Within this treatment, the interparticle spacing a(x) is a slowly-varying function of the position, such that a(x) = 1/n(x). In the thermodynamic limit, in which a(0) is fixed as the number of particles goes to infinity, N → ∞, one recovers the statistical mechanics and dynamical properties of the ion chain in a infinite-radius ring trap [24,9]. For N finite, the transition from a linear to a zigzag chain can be estimated with the value ν (c) t ≈ 3N ν/(4 √ log N ), which was found by taking only nearest-neighbours coupling and where the corrections scale with powers of 1/ log N [24]. Apart from a factor of order unity, it corresponds to the relation with a(0) = 1/n(0). In order to determine a Ginzburg-Landau equation for this case, we first assume that the spatial variation is very slow, such that sufficiently far away from the edges, the length scale over which the interparticle distance changes is much larger than the coarse graining length δx ≫ a(x) with which we study the dynamics: In this limit, we can make a slowly-varying ansatz for the short-wavelength eigenmodes of the linear chain, such that for a given eigenmode we extend the treatment for the homogeneous case starting from Eq. (4) and write σ n = α n e ikna , with α n slowlyvarying amplitude [24]. Within a local-density approximation, the ions contained in a region of size δx and centered at position x, become unstable at the position-dependent critical transverse frequency given by which has the same form as the one for the homogeneous case but with the appropriate interparticle distance a(x). Moreover, extending to the non-homogeneous case the Lagrangian in Eq. (18), we find the Lagrangian L ′ = dxL ′ (x), with the Lagrangian density In Eq. (37), ρ(x) = mn(x) is the linear mass density, and the spatial dependence of the coefficients δ(x), h(x), and A(x) is found by using the position-dependent interparticle distance a(x) in the corresponding formulas for the homogeneous case. When deriving the Lagrangian density of Eq. (37) we have neglected the coupling between axial and transverse degrees of freedom, assuming very small fluctuations about the critical point.
In order to derive the Euler-Lagrange equation for the field, we assume that the linear density, and hence the interparticle distance a(x), does not depend on the value of the transverse trap frequency and thus remains constant when quenching ν t through the critical point. Terms proportional to the spatial gradient da/dx are also neglected, assuming that sufficiently far away from the edges the inequality |ψ ′′ (x)| ≫ |ψ ′ (x)a ′ (x)/a(x)| holds (which is consistent with the approximation made for deriving the Lagrange density). Within this limit, the equation of motion for the field ψ = ψ(x, t) reads where now and we considered the time-dependence given in Eq. (28). We now consider defect formation when the value of the transverse frequency is quenched through the critical point. This will happen in the central region of the chain first, giving rise to a propagating front along the axis, whose coordinates (x F , t F ) satisfy δ(x F , t F ) = 0. The front velocity v F , at which the instability propagates, can be found by taking the ratio between the characteristic length of the control parameter, (∂ x δ(x, t)/δ(x, t)) −1 , over the characteristic time scale at which it changes, An explicit dependence on the physical parameter can be found using the spatiallydependent critical frequency in the expression with X = x/L. For the front velocity one obtains Whenever the transition is homogeneous, v F becomes infinite and the standard KZM applies, allowing for the nucleation of defects in the whole system. Otherwise, as we show in the following, nucleation of kinks will only take place in an restricted fraction of the chain, with length 2X * L, where the front velocity v F is larger than the characteristic velocityv x with which a perturbation propagates along the chain at the freeze-out time. We will find that the scaling of kinks density with the quenching rate will be in this case different from the scaling found in the homogeneous case.
In order to determine the density of defects, we now evaluate the velocityv x with which a perturbation propagates along the chain at the freeze-out time. An upperbound forv x can be found by considering the ratio of the correlation lengthξ x and the relaxation timeτ where the x-subindex underlines the spatial dependence due to the inhomogeneous nature of the system. This also corresponds to the speed of sound at the freezeout point at an energy which is one over the relaxation time. In order to compute it, first note that the relative frequency can be written with reference to t F as t (x) 2 ]/δ 0 . One can find the instantt, relative to t F , at which the dynamics stops being adiabatic by equating the time scale δ/δ to the relaxation time τ x = η/δ(x, t).
In the overdamped regime,t = (ητ Q /δ 0 ) 1/2 , which sets the freezed-out correlation The key insight in the IKZM is that nucleation of defects is only expected whenever the front velocity is larger thanv x , namely, so that spatially separated regions of the chain are causally disconnected. The violation of this inequality opens the possibility of driving a truly adiabatic transition. The reason is that the choice of the ground state in the broken symmetry phase is not independent in different regions of the system whenever v F is small enough with respect tov x . As a result the inhomogeneous nature of the transition allows for an adiabatic crossing, where a zigzag chain results after the quench with no defects. This holds even in thermodynamic limit and dramatically differs from the HKZM where defects would always be expected for infinite systems.
In the overdamped regime, nucleation of kinks is still possible whenever with This will generally be fulfilled in a limited region of the system, 2X * L, where the homogeneous KZM applies. One can estimate the effective size of this region where defect nucleation is possible by setting v F /v x = 1, and assumingX * ≪ 1. Then, it follows thatX * o ≃ A o which leads to the scaling law of the density of kinks Though the absolute density of defects is reduced with respect to the homogeneous case (sinceX * < 1), the dependence on the quenching rate is enhanced. Nonetheless, such IKZM scaling breaks down for fast quenches due to the finite size of the chain, and for very slow quenches due to relevant defects losses as will be discussed below in more detail.
An analogous description applies to the underdamped case, where the relaxation time is independent of the dissipation and diverges as τ x = 1/ |δ(x, t)|, leading to the freeze-out timet = (τ Q /δ 0 ) 1/3 . In this time scale, the correlation length freezes, with respect to the quench time scale τ Q , at a valueξ x = aω 0 (τ Q /δ 0 ) 1/3 leading to a uniform sound velocityv x =ξ x /τ x = aω 0 . Again, the transition remains adiabatic as long as where we introduced the parameter Using the same argument employed in the overdamped regime, kinks are found in a region of sizeX * u ≃ A u , so that To test the IKZM, we study numerically in Figure 5 the scaling of the density of kinks as a function of the quenching rate both in the overdamped and underdamped regimes. The results show a good agreement with the predictions in Eqs. (47) and (50). Nonetheless, there is a saturation of the density of kinks at high quenching rates due to the interactions between the kinks.
We further note that the ratio between the number of defects N that nucleate according to the IKZM and HKZM obeys the relation with a well define scaling with respect to the quenching rate, and which applies both in the underdamped and overdamped regimes. In addition, the HKZM is known to overestimate the number of defects by a numerical factor f of order unity, f ∼ 4.5 in our case.

Discussion
In this section we discuss different mechanisms responsible for deviations from the IKZM scaling. Two new effects come into play with respect to the ring configuration: a) coupling between axial and transverse modes b) enhanced defect transport. We shall dwell on the implications of these two effects in KZM scaling.
Axial-transverse mode coupling -The (inhomogeneous) trapping potential allows the ions to shift in the axial direction as the structural phase transition takes place. As a result, during the course of the transition the axial density of ions increases near the center of the trap, an effect which hinders the study of the IKZM, as expected from the derivation of the GLE in the thermodynamic limit and further suggested by numerical simulations.
Dynamical losses of defects -In order to minimize the axial-transverse mode coupling it is desirable to drive the transition just in the center of the chain, so that the ions at the edges of the trap remain in the linear configuration. For the ground state, the amplitude of the transverse displacement of the ions increases monotonically from each of the edges of the chain towards the center of the trap. Hence, the effective Peierls-Nabarro potential [31,29] seen by a kink decays as one approaches the ends of the chain. As a consequence there is transport of defects, which provides a mechanism for their losses near the edges of the trap. Note that defect transport remains even if the longitudinal degrees of freedom of the ions are frozen on a lattice: the transverse motion suffices for its dynamics. A way of minimizing these defect losses is by making the inter-ion spacing homogeneous, which in turn makes the Peierls-Nabarro potential periodic along the chain. Different potentials have been proposed to achieve homogeneous inter-ion spacing [32]. For an ion chain in a linear Paul trap, the optimal axial trapping potential can be found by fitting the local Coulomb potential in an homogeneous chain, U ({x n }) = n =n ′ 1 |xn−x n ′ | = 1 a n =n ′ 1 |n−n ′ | where a is the desired value of the inter-ion spacing. Under such axial confinement, the transition becomes then more homogeneous. Nonetheless, there is a local correction to the transverse critical frequency, ν 2 r,ef f (n) = ν r (t) 2 − n =n ′ e 2 /m |xn−x n ′ | 3 which varies as a function of the position along the chain. Hence even when the ions are homogeneously spaced, the transition remains inhomogeneous. The underlying defect dynamics is also driven by the fact that pairs of defects with the same topological charge repel each other and attract otherwise. Scattering between kinks and anti-kinks (of positive and negative topological charge) can occur leading to their annihilation. Nonetheless defects stop seeing each other when they are separated by few ions (∼ 5), which motivates the study of thes scaling with moderate densities. Under the conditions in Fig. (5) defects are stable on a scale ∼ 10ν −1 which suffices for its imaging in the laboratory. The mechanisms for defect losses are particularly relevant in the underdamped regime, where higher scaling coefficients than those predicted by the KZM are observed and the fluctuations in the number of kinks increase. Generally, whenever defects losses are relevant the connection of the scaling of the density of defects with KZM can be questioned, since these process are non-universal and disregarded in the mechanism. As a result, the overdamped regime is single out as optimal for studies of the IKZM.
We close by noticing that the minimal length of the chain required to verify the IKZM scaling is constrained by two finite-size effects: a) The existence of the adiabatic dynamics whenever the correlation length at the freeze out time,ξ = ξ(t), equals the length of the part of the chain where defects are counted, and more fundamentally, where defects can nucleate 2X * . Defects nucleate as long as 2X * L >ξ. b) The breakdown of the IKZM scaling at fast rates due to a saturation of the average density of defects. If one wishes to check the scaling by varying the quenching time between τ i and τ f (ideally ranging over few orders of magnitude), it should be possible to achieve the corresponding densities of defects (d i and d f ) related as d i = (τ f /τ i ) α d f , where α is the IKZM scaling (1 in the overdamped regime, 4/3 in the underdamped regime), without jumping into the adiabatic dynamics. Long enough chains were created in various setups. In [2] a chain of up to 5 × 10 4 ions was created with an axial trapping frequency of 1 MHz. Ring configuration were created dynamically in [33]. In the same setup static structures of up to 50 ions were created at the temperature of 1 mK.

Conclusions
In this work we analyzed the formation of defects in a one dimensional Coulomb crystal during a frequency quench from the linear chain to the zigzag structure. We studied the cases of ions crystals in a linear Paul trap and in a ring-shaped trap and predicted the defects production rate. Our study shows that a Coulomb crystal is a particularly neat and controllable system where the homogeneous and inhomogeneous KZM can be tested. Despite the experimental work in [19] addressing the homogeneous KZM (HKZM), the inhomogeneous KZM (IKZM) [22,23] lacks to date experimental verification, and the ion chains in a Paul trap are put forward here as an ideal system for such goal.
where for simplicity we just study the dynamics of the ions in one direction y and we inverted Eq. (4) for the expansion of the particles coordinates in terms of normal modes. In the overdamped regime: η ≫ |β(k)| we can neglect the second order derivative in Eq. (22) whose solution becomes Substituting this solution in the definition of g t (τ ) Eq.(A.1) we get: Now for fixed time t, all the components with different momentum in the autocorrelation function g t (τ ) decay exponentially with the time separation τ . The largest time scale, corresponding to the minimum frequency δ = min k β(k), defines the relaxation time: where the subscript o stresses that this result has been obtained in the overdamped regime. Note also that defects are localised excitations and a correction depending on h can be obtained from the two-point correlation function.