Dynamics of topological defects in ion Coulomb crystals

We study experimentally and theoretically the properties of structural defects (kink solitons) in two-dimensional ion Coulomb crystals. We show how different types of kink solitons with different physical properties can be realized, and transformed from one type into another by varying the aspect ratio of the trap confinement. Further, we discuss how impurities in ion Coulomb crystals, such as mass defects, can modify the dynamics of kink creation and their stability. For both pure and impure crystals, the experimentally observed kink dynamics are analyzed in detail and explained theoretically by numerical simulations and calculations of the Peierls-Nabarro potential. Finally, we show that static electric fields provide a handle to vary the influence of mass defects on kinks in a controlled way and allow for deterministic manipulation and creation of kinks.


Introduction
Recent experimental realizations of structural defects (kinks) in ion Coulomb crystals [1,2,3,4] have raised an increasing interest in nonlinear physics in laser cooled ions, and have been proposed as a tool for storing quantum information [5].
Topological defects are created when a symmetry breaking phase transition is crossed non-adiabatically, such as the second-order phase transition from the linear to zigzag configuration in ion Coulomb crystals [6,7]. Kinks were first observed during crystallization due to laser cooling [8]. Since then, several groups have succeeded in stable trapping of defects, either created through crystallization [1], or by axial [3] or radial [2,4] quenches of the confining potential. The storing of quasi-3D kinks and study of the interaction of different classes of defects involving mass defects are also under investigation [1,9,10] We have previously used the creation of topological defects during rapid radial quenches to study the dynamics of phase transitions and to test the predictions of Kibble-Zurek theory [2,3,11,12,13]. In this work, we investigate the stability and dynamics of kinks, including the influence of mass defects present in the ion crystal and the controlled manipulation of kinks using electric fields. The dynamics occur on microsecond timescales and are difficult to resolve experimentally. We interpret our observations at experimentally accessible timescales using simulations that provide information about the rapid dynamics and using numerical calculations of the Peierls-Nabarro (PN) potential [14], which describes the energy necessary for the kink to move within the crystal.
The paper is organized as follows. Section 2 introduces our experimental system and describes the kink structures. Section 3 treats the kinks as discrete solitons and introduces our method for calculating the PN potential. Section 4 discusses motion of kinks in the PN potential. Section 5 investigates kinks in the presence of mass defects; in our system, heavy ions can be created by the chemical reactions with background gas. In section 6 we discuss ways in which we can create and manipulate kinks using an applied dc electric field.

Creation of localized and extended topological defects
In the experiment we use an ion Coulomb crystal of around 30 172 Yb + ions confined in a segmented linear Paul trap which is designed for trapping large crystals with minimized micromotion [15,2]. Linear chains of ions are trapped with axial frequency ν z = 24.6 ± 0.5 kHz and radial frequencies ν x , ν y ≈ 500 kHz. When the radial confinement is weakened below a critical value such that ν x,y /ν z < 0.73 × N 0.86 where N is the number of ions [16], the linear chain buckles into a zigzag configuration. The kinks are created by rapidly quenching the radial confinement across this transition, which causes defects to occur probabilistically [2]. A radial asymmetry in the potential of ν y /ν x ≈ 1.3 is introduced by the end of the quench in order to confine the zigzag phase to a well-defined plane. Reducing this ratio leads to quasi-3d kink structures [1], however in this work we focus on two-dimensional configurations.
A kink can take on different configurations depending on the ratio of the transverse and axial confinement, ν x /ν z . Figure 1 shows typical experimental images of these kink structures. For high values of ν x /ν z (i.e. for small transverse inter-ion spacing) the kink has a localized structure shown in figure 1(a) which involves only ∼ 3 to 5 ions. We refer to this type of kink as "odd", using the terminology of [5]. For low values of ν x /ν z (large transverse inter-ion spacing) the kink has an extended structure as shown in 1(b) that involves ∼ 7 to 10 ions, which we refer to as an "extended" kink.
To create odd kinks the radial confinement ν x is relaxed from approximately 500 kHz to 204 kHz, and for extended kinks ν x is reduced from approximately 500 kHz to 140 kHz. After the quench, the spatially resolved ion configuration is imaged at a 45 degree angle in the transverse direction onto an electron multiplying CCD camera using a typical exposure time of 10 to 50 ms. Resulting images are evaluated for the presence and location of kinks. The quench parameter τ Q = t r /2, where t r is the length of the radial quench, compared with the timescale of the secular motion ν z ≈ 24.6 kHz determines the probability of kink creation (see [2]).
The kinks can propagate along the chain if they obtain kinetic energy ( 5 mK) from, for example, the ramp of the radial potential (this can be on the order of tens of mK), collisions with background gas particles, or a decrease in the laser cooling rate. Kinks that have moved sufficiently close to the end of the zigzag region are lost at the edges. Ion motion happens on the timescale of the secular frequency, making kink motion impossible to resolve with exposure times of tens of milliseconds. However, as indicated by a superimposed kink and pure zigzag image, as in figure 1(c) and (d), the trace of a kink can be observed moving out of the chain during the camera exposure. Experimentally, we find that the probability of observing an extended kink is higher than the probability of observing an odd kink, which suggests that odd kinks may exhibit higher losses.
In order to investigate the rapid dynamics of these losses we have performed molecular dynamics simulations, from which we observe that defects are first created as odd kinks in the early part of the quench, and transform into extended structures as the ratio of radial to axial confinement (ν x /ν z ) is lowered further. We observe that odd kinks are easily lost from the edge of the chain when they have enough energy to travel from site to site; in contrast, it is possible for extended kinks that are moving along the chain to remain in the crystal by oscillating back and forth and eventually stabilizing in the centre. Since the simulations show that when driving phase transitions via radial quenches, extended kinks are always created from odd kinks, the extended kinks provide a convenient stable probe of the statistics of defect creation (in finite systems with small kink probabilities d < 1). We have used this characteristic of the extended kinks to facilitate the testing of Kibble-Zurek scaling in this regime in past work [2].

The Peierls-Nabarro potential
The PN potential is the effective potential of the kink and is defined as the energy of the system as a function of the location of the kink centre, X. Before the PN potential is calculated it is necessary to understand qualitatively how a kink moves and to formulate a definition of the kink centre. To facilitate this discussion, for a zigzag we define the inter-ion spacing in the transverse direction as a and in the axial direction as b.
Odd kink. Figure 2(a) (right) illustrates the simulated propagation of an odd kink. For this type of kink, a < b and the kink propagates via the transverse motion of the two ions that are at the interface between the two zigzag domains. This motion is observed experimentally in figure 1(c). Let (z j , x j ) and (z j+1 , x j+1 ) be the coordinates of the ions at the interface between two domains. We define the centre of an odd kink as the location where a line connecting points (z j ,(−1) j x j ) and (z j+1 ,(−1) j+1 x j+1 ) intersects with the z-axis. This point is marked on figure 2(a) by a solid vertical line.
Extended kink. An extended kink has a different structure and mode of propagation. For a zigzag chain with an extended kink a > b, the kink contains an excess of one ion in either the bottom or the top row of the zigzag. The kink is identified as the region of the excess charge density and it propagates via the axial motion of ions, as illustrated in figure 2(b) (right). This motion is visible in the experiment in the exposure of figure 1(d). We take the kink centre as a geometric centre of the shift in axial coordinates that a kink introduces into the chain (see appendix A for a precise definition).
As the radial confinement is decreased during a quench, the transverse separation a of the zigzag increases, to create first odd (a < b) and then extended (a > b) kink structures. In the same manner, a quench of increasing radial confinement transforms an extended kink back into an odd kink. When a ≈ b the propagation of the kink involves both transverse and axial motions. Figure 2(c) illustrates this mode of propagation for a kink that occurs during the transformation between the odd and the extended type and can appear as either an odd or extended kink at different moments during the movement. This kink structure was also experimentally observed in [4]. The kink centre for this intermediate kink is determined using the same method as for the extended kink.
In order to find the PN potential for a kink in the zigzag ion crystal, we calculate the adiabatic trajectory connecting a configuration with a kink at one end of the crystal to a configuration with the kink at the other end. To find a minimum energy crystal configuration with a kink at position X we perform energy minimization constrained on the kink being at position X using the method of Lagrange multipliers (see appendix A). Effectively, we start with a kink configuration at a stationary point (a local or global minimum of the PN potential) and then move the kink against the energy gradient using constrained energy minimization. In this way the whole trajectory can be traced out in order to determine the PN potential.
Calculated PN potential curves are shown in figure 2 for the system of 30 172 Yb + ions with an axial confinement of 24.6 kHz and (weak) transverse confinement (ν x ) at selected values ranging from 140 to 236 kHz. First we consider the PN potential curves for the odd kinks shown in figure 2(a). It is expected that the PN potential should be periodic in lattice spacing, since this periodicity is a fundamental feature of discrete soliton systems. For an infinite and homogeneous system, the periodicity arises from the invariance of a kink under translations by an integer number of lattice spacings. However, for a finite sized system, this invariance is only approximately true and the form of the PN potential is modified by the boundaries of the system. This can be seen in figure 2(a): the PN potential is periodic in lattice spacing, but tends to decrease as the kink moves away from the centre of the chain. The height of the PN potential barriers, which indicate the minimum energy that is required to move a kink from one lattice site to the next, increases as ν x is decreased. In order for a kink to be trapped, its kinetic energy must be lower than the height of the barrier. During the early stages of the quench, the created odd kinks can easily be lost, since the PN barrier is smaller than the kinetic energy in the system and the PN gradient pushes the kink outward. As the quench progresses, the PN barrier increases until it is eventually able to trap kinks, which are then experimentally observable.
The PN potential curves for extended kinks are shown in figure 2(b). For lower values of ν x , this potential is also periodic, but with a period about two times larger than that for odd kinks, reflecting the more extended nature of the kink structure. The PN barriers decrease with decreasing ν x , since as the distance between the two rows of ions increases, the discreteness of the lattice becomes less prominent: the kink experiences an averaged out and thus more continuous charge distribution. In this case, the effect of the system boundaries is that the PN potential curves upward at the edges of the zigzag region, where the charge density is higher. Thus, an extended kink tends to be pushed away from the boundaries and trapped at the centre. When an extended kink is created from an odd kink, it will oscillate back and forth in this PN potential and remain trapped in the centre of the chain. This stabilization mechanism also leads to the annihilation of pairs of extended kinks in the centre of the crystal, so that the number of observable extended kinks saturates at one defect in an inhomogeneous (pure) crystal of sufficiently small size, as also noted in [1].
As the odd kink transforms into an extended kink, the PN potential must transform from the periodic potential that decreases at the edges into one with a global minimum. This happens temporally when a ≈ b and the defect becomes an intermediate kink, where the propagation of the kink involves the motion of ions in both the transverse and axial directions. With this mechanism of movement, the kink propagation is such that the kink moves in steps of two lattice sites, instead of one as in the odd kink case. This leads to the change in the periodicity of the PN potential in the intermediate stage.
In figure 2(c) we see that, although this transition from odd to intermediate to extended kinks happens temporally during the quench, there is also a spatial dependence of the transformation: at the crystal centre, intermediate and then extended kinks can form, while the kink continues to behave as an odd kink near the edges.
These numerically calculated potentials provide an explanation for the observed higher stability of extended kinks compared to odd kinks, as well as the oscillating behaviour seen in the simulations. The global minimum of the extended kink PN potential makes it more difficult for an extended kink to escape the crystal, while in the odd case, a kink that is already moving within the chain is easily lost as it approaches the edge of the chain. Therefore, for the same amount of kinetic energy in the system, the average probability for an odd kink to be lost is higher than that for an extended kink.

Kink motion
The PN potential predicts some qualitative features of the kink dynamics, such as trapping and losses of kinks at the boundaries. By itself the PN potential does not provide quantitative dynamical information such as the frequency with which a kink oscillates around a local minimum of the PN potential. To obtain such dynamic information it is necessary to derive the effective equation of motion for the kink centre. This can be done, for example, using the collective variable formalism (see Appendix B) [17]. The equations of motion for kinks typically contain an inertial term, a non-linear friction term that arises from the coupling between kink motion and the phonon bath, and the PN force term, which is obtained by differentiating the PN potential. Knowing the equations of motion for kinks is important for quantitative analysis of interactions between kink modes and phonon modes of the chain, as well as the analytic predictions of quantities such as the frequencies of kink oscillations. Here, we take a simple approach and use numerical simulations to measure the oscillation frequency of an extended kink. We also use numerical simulations to look at the response of the kink to thermal driving of an ion chain.
We consider an extended kink of 30 Yb + ions in a trap with ν x = 140 kHz and ν z = 24.6 kHz. The PN potential for this kink is shown in figure 2(b). Molecular dynamics simulations show that these kinks oscillate about X = 0, which is the point of global minimum of the PN potential. To observe this oscillation in its pure form, we initialize the positions of ions in a configuration taken from the adiabatic trajectory, set the initial ion velocities to zero, and numerically evolve the equations of motion without damping or stochastic forces. Once ion trajectories have been found, the position of the kink centre at each time instance is calculated using equation A.3. For initial conditions X(0) = 25 µm andẊ(0) = 0 the trajectory of a kink is shown in figure 3(a) (T = 0 mK line). The oscillation is very regular with no noticeable damping, which indicates weak coupling of the kink to the phonon modes of the crystal. The frequency of the oscillation is 11.9 kHz. Next, we consider the excitation of the motional kink mode by thermal noise. The thermal dynamics of the chain is modelled using the Langevin equation, where the effect of the thermal bath is modelled by the addition of frictional and stochastic forces whose statistical properties satisfy the fluctuation dissipation theorem (see [2] for more details concerning the simulation approach). The kink is initialized at X = 0 and its evolution is calculated numerically for the two temperatures T = 3 mK and T = 5 mK. A 400 µs sample interval of the thermalized evolution of kink centre is shown in figure 3. The kink oscillates with a frequency of ≈ 12 kHz, and the amplitude is larger for the crystal at higher temperature as expected. These simulations allow us to estimate the timescales and amount of heating necessary to move kinks within our ion Coulomb crystal.

Mass Defects
Next we consider how an impurity in the Coulomb crystal influences the kink solitons. A mass defect in the ponderomotive potential of a Paul trap experiences a different radial confinement than the rest of the ions, since the secular frequency ν x ∼ 1/m where m is the mass [19]. Because of this, the mass defect changes the configuration of the crystal and deforms the PN potential, consequently affecting the dynamics of the kinks. This modification of the PN potential depends on the location and the mass of the impurity species and the type of kink. We primarily consider extended kinks with mass defects heavier than the surrounding 172 Yb + ions in this analysis. Off-centre kinks were rarely observed in pure crystals. With molecules present, the probability to have an off-centre kink increases. This indicates trapping of kinks by molecules. bottom: Definition of the location of a kink as left, right, or centre. The centre of the kink must be outside of the middle three sites to be considered off-centre. (b) Scaling of measured kink density with the quench parameter τ Q , with and without one molecule present in the chain. Losses at slow ramps are reduced when a molecule is present. inset : Scaling of kink density for two molecules in the chain.
In our experiment, heavy mass defects are implanted by the creation of molecules during collisions of the 172 Yb + ions in the crystal with background gases. In a fluorescence image of the chain, a molecular mass defect appears as a dark ion. The 172 Yb + ions can react with water to produce YbO + , Yb(OH) + , and other molecules [20,21,22]. By performing parametric excitation of two to four ions including a mass defect, we have measured primarily YbO + or Yb(OH) + ; the addition of hydrogen was not resolvable by this method. These heavy masses distort the crystal as well as the kink configurations. The amount of distortion can be used as an alternative method to estimate the mass of the molecule by comparing the experimental images with crystal configurations inferred from molecular dynamics simulations. This analysis indicates the presence of even heavier molecules with one to three oxygen atoms involved. In the simulations, we investigate mass defects with representative masses m = 188 (YbO + ) and m = 220 (YbO + 3 ) (all masses are given in atomic mass units), consistent with observations in our system. Figure 7(a) shows three examples of structural defects including a molecule. The first image is of a crystal containing a light molecule (e.g. YbO + ), while the greater distortion evident in the other two images indicates the presence of a heavy molecule (e.g. YbO + 3 ). In the experiment, the position of mass defects within the crystal is not controlled and varies over the crystal.
Experimentally, we see that the spatial distribution of kink positions within the crystal is altered in the presence of a mass defect. This indicates that a mass defect can stabilize a kink at a different location than the centre of the crystal. Figure 4(a) shows a histogram of kink location with and without molecules present for, in each case, about 2000 observed extended kinks. In the case with molecules, the appearance of defects in the left and right sides of the crystal becomes much more probable. This indicates that kinks which would otherwise either escape the chain or settle in the centre of the PN potential after oscillating are instead stopped at the molecule.
We also see evidence that the probability of kinks to be created (i.e. the kink density d) and trapped during the radial quench is modified by the presence of molecules. Figure  4(b) shows a comparison of measured kink density as a function of the quench parameter τ Q for crystals with mass defects compared with the case of pure crystals.

Peierls-Nabarro potential with impurities
The modified spatial kink distribution can be understood by evaluating the PN potential in the presence of mass defects. Figure 5 shows the PN potential for kinks in a crystal with a molecule at fixed positions. These PN potentials have a local minimum at the mass defect. Figure 5(a) and (b) show that for extended kinks, the local minimum approaches a global minimum as ν x is lowered, whose axial location in the chain is modified from the centre to the location of the molecule. This is the case we have investigated experimentally. Since the molecule is less bound in the radial potential than the rest of the ions, it is easier for the molecule to escape into the transverse plane, thus lowering the charge density and the potential energy. This allows stable trapping of kinks at the location of the mass defect. Moreover, in the case of more than one mass defect, the creation of local minima at the positions of the mass defects allows trapping of more than one extended kink in such mixed species crystals, despite their small size (see black line in figure 5(b). Experimentally, we observed the stable trapping of two extended kinks in the presence of two heavy molecules (see figure 7(b)). Finally, we also calculate the PN potential for "odd" kinks in the presence of molecules. In this case, shown in figure 5(c), the molecule also induces a deep local minimum that can trap the kink before it leaves the chain. For both odd and extended kinks, the mass defect can therefore stabilize a kink nearby to its position.

Kink creation with impurities
In order to investigate the influence of mass defects on kink creation, the dynamics of the phase transition must be considered. This includes an interplay of the propagating phase front velocity and the speed of sound [2,23,24,25,26]. In addition, the weaker radial confinement of the mass defect causes the phase transition to set in at an earlier time at the position of the molecule. These effects introduce a dependency of the defect probability on the quench rate. The inhomogeneous nature of the crystal due to the harmonic trap confinement also leads to an expected dependence on the position of the mass defect. For this reason it becomes necessary to consider kink creation as a function of both the quench rate and the mass defect position. In our numerical simulations we typically vary the quench parameter τ Q from 15 to 100 µs and deduce the kink density with mass defects placed at fixed positions. A simulation of the density of extended kinks with a molecule of mass 220 at various positions along the ion chain for different ramp times is shown in figure 6 for a crystal of 30 ions. The position numbers indicate the number of ions counted from the left, where the leftmost ion is at position zero. Because of symmetry only half of the positions need to be simulated. Two cases are shown: first, including all kinks that are created during the ramp (figure 6(a)), and second, kinks that are still present, and therefore measurable, 400 µs after the ramp when a stable kink density is reached ( figure 6(b)). For each case, the kink density without molecules is plotted for comparison (black line). First, we consider the case of the created kinks in figure 6(a). For fast quenches the molecule has no effect on kink creation, but for slow quenches, the kink density exhibits a strong position dependence. When the molecule is in the centre (e.g. positions 13 or 15), kink formation is suppressed for slow radial quenches. Because the molecule is less bound radially than the rest of the crystal, it breaks out of the linear chain earlier and initiates the phase transition sooner than in a pure crystal. This leads to a prolonged time period in which the phase transition spreads out, suppressing defect formation for sufficiently slow quenches [23,24]. At intermediate sites (e.g. positions 9 and 11), the kink density is enhanced. Here, in addition to the quench induced phase transition, which always begins at the crystal centre due to its inhomogeneity, a second phase transition is initiated at the molecule's position. This occurs because of the lower radial confinement of the mass defect and the higher charge density near to the centre pushing the molecule out of the chain. These two phase transitions take place independently, bringing about an enhanced probability of kink formation due to the possible conflict when the phase fronts meet. With the molecule at the edges of the chain (positions 1-7), the kink density is unaffected by the presence of the molecule. In the outer parts of the crystal the molecule does not experience a sufficient repulsion in this region of lower charge density to initiate the independent phase transition. The molecule position's influence on kink creation behaviour is reduced at faster ramp times and eventually has no effect (i.e. for τ Q < 30 µs). In this regime the system behaves as a homogeneous system (the phase front is faster than the speed of sound), and this characteristic dominates the kink creation dynamics.
The case of experimentally accessible kinks is shown in figure 6(b). The density of kinks remaining after 400 µs is mainly governed by the altered rate of kink creation shown in figure 6(a), with two exceptions. First, in positions 9 to 11, the losses from the crystal are strongly reduced for slower ramps. For these intermediate positions, the PN potential at the molecule provides a very strong confinement that traps kinks formed at the centre as they start to move out of the crystal, reducing losses. Second, at position 5, high losses occur. When the molecule is around position 5, it is easy for the kink to be lost due to the reduced barrier of the PN potential caused by the molecule, as seen in figure 5(a).

Mass defects with electric fields
Finally, we discuss the influence of an additionally applied electric dc field on the creation and dynamics of topological defects in a crystal containing mass defects. A heavier ion experiencing a weaker radial confinement than the 172 Yb + ions is displaced farther from the trap axis by an external field. The radial displacement relative to the lighter ions can therefore be controlled, causing the mass defect to behave as a molecule with increased effective mass. The continuous tunability of the effective mass of the molecule provides a rich playground of methods to influence kink creation and behaviour. In this section we discuss several techniques for extended kink manipulation via the application of electric fields.
First we consider the case of applying a time-independent electric field in the xdirection during a quench of the radial potential. Here, in the linear phase the molecule is already displaced from the linear chain of ions in the preferred direction of the electric field. The amount of displacement depends on the mass of the molecule and the . Image (i) shows a kink in the center of the crystal after the radial quench. Image (ii) shows the crystal after the adiabatic application of an electric field in the plane of the zigzag (downward in the image). In (iii) the laser cooling light is shifted a few MHz closer to resonance, heating the Coulomb crystal and providing the kink with motional energy. In this image, motion during the exposure is visible when the kink begins to oscillate in the PN potential, gets trapped by the mass defect, and returns to the centre with the molecule. The kink and molecule then stabilize (iv).
amplitude of the electric field. Since the molecule is already displaced from the linear chain, it initiates the phase transition as soon as the quench is applied. This leads to an amplification of the position-dependent suppression and enhancement of kink formation discussed in the previous section for the molecule without an external field. As the radial confinement is lowered, not only is the mass defect displaced further from the crystal, but the entire ion chain moves to a new equilibrium position. Compared to the case without additional fields these rapid movements increase the kinetic energy in the system by approximately one order of magnitude, depending on the quench time τ Q and the strength of the electric field.
Next we describe a method for using the PN potential to control the dynamics by moving a molecule together with a kink to the centre of the crystal. To do this, we first apply the radial quench without an additional field to produce a kink and subsequently apply the dc electric field. This method circumvents the problem of excitations due to a fast translation of the whole crystal. In simulations we have observed that after the ramp, the PN potential minimum at the position of the molecule causes a passing kink in the presence of laser cooling to be trapped at the molecule. The emerging kink configuration pushes the molecule farther outside of the ion crystal (as in figure  7(a)). When the displacement of the molecule reaches a critical point, the influence of the overall axial trapping potential becomes stronger than the interaction with the neighbouring ions, causing the molecule to move to the centre. Because the PN potential minimum is always located at the position of the molecule (see figure 5(a)) the kink is dragged by the mass defect to the centre of the crystal. The field amplitude that is needed to induce this effect depends on the molecule mass and position. In simulations we found that the behaviour sets in at electric field values from ≈ 10 V/m (m = 220) to ≈ 100 V/m (m = 188), in contrast to the case without a kink involved, where approximately twice as much field strength is required to move the molecule to the centre of the chain. We also observed a position dependence due to the inhomogeneity of the crystal: defects near the ends of the zigzag region require a higher electric field strength than molecules near the center in order to trap a kink moving with the same kinetic energy.
Although the dynamics occur at microsecond timescales, we can induce and observe this behaviour experimentally in steps. Figure 7(c) shows an example of the experimental sequence with a molecule of m ≈ 188. After a quench of length 2τ Q , a series of pictures is taken with exposure times of 20 ms, with 20 ms between exposures. After a kink has been created by the ramp (i), we increase the electric field linearly to ≈ 110 V/m over 60 ms (ii). In (iii) we heat the entire ion chain by reducing the detuning of the cooling laser by a few MHz. This causes the kink to oscillate in the confining PN potential, similarly to what is shown in figure 3, until it reaches the position of the molecule. As the kink gets trapped at the molecule it pushes the molecule further out of the crystal where it is free to move to the centre of the axial potential. When the molecule changes its position, the PN potential minimum and the kink trapped therein move along with it. The kink and molecule finally stabilize at the centre of the crystal (figure 7(c)(iv)). We have verified that for the parameters used, the molecule does not move in the crystal without the presence of a kink.
Last we look at two techniques that could be used for deterministic kink creation that we have explored using simulations, and have also been observed in the experiment. First, using a time-independent electric field during the radial quench, we have the possibility for two molecules in the appropriate positions to create kinks deterministically, because the two molecules always move in the same preferred direction when entering the zigzag phase. For example, with an even (odd) number of ions between the molecules in the linear chain, a kink must (must not) form. Simulations indicate that an external electric field of 10 to 20 V/m is sufficient to produce this effect.
The second method is independent of the molecule position and starts with a crystal in the zigzag phase with a single mass defect. No quench of the trapping potential is necessary to produce the kink. To create a kink we ramp up the electrical field. Again, as the field increases the molecule is increasingly displaced from the zigzag, until it is sufficiently far outside the chain that it becomes decoupled and moves toward the centre of the crystal in the axial trapping potential. When the heavy ion leaves its position in the chain it produces two kink defects. One of them behaves like a kink in a chain without molecules and starts to move in the confining PN potential. The other one is directly located at the molecule and together they move to the centre of the chain, because the axial trapping potential dominates the motion of the molecule. The production of the two defects can be understood more easily in terms of lattice sites. As the molecule leaves its original lattice site there is a vacancy created in the lattice. This unstable vacancy then produces the first kink. The latter defect is created, because there is no available lattice site for the molecule. It can be interpreted as an interstitial defect, which then converts into another kink. This can also be seen as the production of a kink and anti-kink. Figure 8 shows two experimental time series demonstrating this method applied to crystals with one and two mass defects.

Conclusion
We have investigated the confining Peierls-Nabarro potentials of localized and extended kink solitons in two-dimensional ion Coulomb crystals. Each type of kink can be transformed into the other type by varying the aspect ratio of the trap potentials. Further, we have shown that a heavy mass defect can trap and stabilize a kink at its location by deforming the PN potential and creating local minima. This allows trapping of more than one extended kink soliton in the Coulomb crystal and offers the possibility to study kink-kink interaction in the future. The presence of mass defects displays a complex dynamics during the phase transition, leading to enhanced and suppressed kink formation. In our experiment we have used molecules to study the influence of mass defects; however, for a more deterministic control of mass defects, another addressable ion species can be used.
Tuning this influence of mass defects with controlled electric dc fields offers a rich playground for the manipulation of kink solitons. In particular, we have shown first evidence of kink solitons being deterministically created by electric field ramps and the manipulation of an extended kink attached to a mass defect. Since the extended defects can be transformed back into localized kinks, these techniques could also be applied to other various types of kink solitons. The techniques that we have established here help to open up the field of study in the rich physics of kink solitons in ion Coulomb crystals. Besides the dynamics and interaction of multiple kinks, this may also include the study of internal vibrational modes of kinks [1,5] and ultimately, the quantum dynamics of kinks and quantum statistical mechanics in ion chains [7].

Appendix A. Calculation of the Peierls-Nabarro potential
In this section we explain in detail the procedure that was used to calculate the PN potentials. First we define the function used to calculate the kink location of the kink centre. Then the problem of finding the adiabatic trajectory and hence PN potential is formulated as an optimization problem that can be solved using the method of Lagrange multipliers.
Appendix A.1. Definition of kink centre Consider a harmonic trap such that a single charged particle of mass m would oscillate in it with frequencies ω x and ω z . The potential energy of the system of N particles trapped in such a trap is given by where m j is the mass of jth ion. Suppose the charges are distributed in a zigzag kink configuration; the z-axis being the longitudinal direction and x-axis the transverse direction. We denote the charge configuration by a vector Q = (z 1 , . . . z N , x 1 , . . . x N , ), where the ions are numbered such that z 1 ≤ z 2 . . . ≤ z N . The centre of the kink, X, is a function of Q i.e. X = g(Q). For odd kinks we choose with x j x j+1 > 0 i.e. the kink is located between the jth and (j + 1)th ion. For extended kinks we choose where z (0) j denotes the z coordinate of an ion in chain in a configuration with no kink. Equation (A.3) is effectively a weighted average. It measures the location along the z-axis of the centre of charge distribution distortion that is introduced by the presence of the kink.

Appendix A.2. Evaluation of adiabatic trajectory
The adiabatic kink trajectory is defined as a trajectory that takes a kink at position X 1 to position X 2 with minimum amount of work. Let f(X) denote the configuration of particles corresponding to a kink at position X along the adiabatic trajectory. By definition f(X) has to minimize the energy function whilst satisfying g(f(X)) = X. Thus f(X) is found by solving the following optimization problem where ∇ ≡ (∂/∂x 1 , . . . , ∂/∂x N , ∂/∂z 1 , . . . , ∂/∂z N , ∂/∂λ) . Equation (A.6) can be solved numerically by using iterative gradient descent methods. It is best to select an initial guess that is close to the true solution in order to not get trapped in a wrong local minimum during the implementation of gradient descent algorithm. Typically, we take a kink configuration in a known stable stationary configuration f(X 0 ) and use it as an initial guess for finding f(X 0 +δX), where δX is small. Once f(X 0 +δX) is found, we use it as an initial guess for finding f(X 0 + 2δX) and so on. This way the whole adiabatic trajectory of a kink can be traced out. The potential energy of a kink along the adiabatic trajectory is the Peierls-Nabarro potential U(X) i.e. U(X) = V (f(X)). (A.7) Note that if we know that U(X) monotonically decreasing between some two points X 1 and X 2 then the adiabatic trajectory can simply be evaluated by initializing a kink at point X 1 and calculate its dynamics under heavy damping. The kink will "slide downhill" along the adiabatic trajectory from the point X 1 to the point X 2 . Where as Newton's laws always move the kink along the energy gradient, the constrained optimization can be used to push the kink with or against the energy gradient.

Appendix B. Equation of motion of the kink
Here we explain how the equation of motion of the kink can be derived using the collective-variable formalism. Some of the earliest work developing this approach was in the studies of solitons in the Klein-Gordon quantum field theories [27,28]. The theory was developed further in the context of discrete kinks [29] and used in a number of works analyzing the dynamics of discrete kinks such as Sine-Gordon kinks [30,31,32,33]. The first step in the derivation of the equation of motion is to write the vector that denotes the state of the system as a sum of two parts: a contribution from the kink configuration taken from the adiabatic trajectory and an additional perturbation We may refer to the configuration f(X) as a bare kink configuration and the vector q = (q 1 , . . . , q N ) as kink dressing, which may come for instance from the excited phonon modes of the chain. In the collective variable formalism the variable, X, is promoted to a dynamical variable and a corresponding canonical momentum, Π, is introduced. Thus to derive the equations of motion one must perform a canonical transformation from the original 4N variables Q j , P j to the new (4N + 2) variables q j , p j , X and Π. Since the number of degrees of freedom is increased by two one has to introduce two constraints. A good choice of the constraints is C 1 = f(X) · q = 0, (B.1) where prime denotes differentiation with respect to the argument e.g. f ′ (X) = (∂f 1 (X)/∂X, . . . , ∂f 2N (X)/∂X). Constrained (B.1) and (B.2) are chosen in order to minimize q j in the vicinity of the kink. The canonical transformation was performed in [30] and it was also proven that the result can be obtained using a simpler projection operator approach [29]. Here, we quote the result n (δ ln − P ln ) q n + f ′′  4) should give the kink dynamics, which is identical to the dynamics obtained by solving Newton's laws of motion. The collective-variable equations can give a better insight into kink dynamics since they explicitly keep track of the motion of the kink. They can serve as a starting point for predicting dynamical features such as, for example, the PN oscillation frequency and diffusion coefficients of kinks.