Dynamics of skyrmion in disordered chiral magnet of thin film form

Magnetic skyrmion is a topological spin texture characterized by the mapping from the two dimensional real space to the unit sphere. It is realized in chiral magnets under an external magnetic field in the plane perpendicular to it. In thin film samples, which are most relevant to the applications, the thickness of the system parallel to the magnetic field is finite, and a skyrmion turns into a skyrmion string, which is often assumed to be a straight rod. There are phenomena related to the internal degrees of freedom along the string, e.g., the monopole and anti-monopole creation/annihilation, corresponding to the change in the skyrmion number. However, the role of this finite thickness in the topological stability and dynamics has not been explored yet. Here we study theoretically the current-driven dynamics of a skyrmion string under disorder potential by systematically changing the thickness of the sample to reveal the dynamical phase diagram in the plane of current density and thickness. We found the three regions, i.e., (i) pinned skyrmion string, (ii) moving depinned skyrmion string, and (iii) annihilation of skyrmion string, for thin and thick limits while (iii) is missing in the intermediate case. This indicates that there is the optimal range of thickness for the topological stability of skyrmion string enhanced compared with a two-dimensional skyrmion. This result provides a way to design and control skyrmions in thin films and interfaces of finite thickness.

Magnetic skyrmion 1-16 is a two-dimensional topological spin texture characterized by an integer called skyrmion number which counts how many times the spin directions wrap the unit sphere. This topological nature leads to the variety of novel properties of skyrmions such as their topological stability, high mobility with very small threshold current density, and topological Hall effect. The crystal form of skyrmions has been first observed in three dimensional bulk crystal of chiral magnet MnSi by the neutron scattering experiment 6 , where the skyrmionic phase is limited to the small region in the plane of temperature and magnetic field. In three-dimensions, the skyrmions form the rod-like string along the direction of the external magnetic field. Later it has been recognized that the skyrmion phase is more stable in two-dimensional systems such as thin film by the real-space observation in terms of the Lorentz electron microscopy in (Fe,Co)Si 8 . In two-dimensions, a skyrmion is a point-like particle. From the viewpoint of applications, the proposals of skyrmionic devices have been discussed (e.g., refs [17][18][19], and interfaces or thin films [20][21][22][23][24][25][26][27][28] offer an important laboratory since it naturally introduces the Dzyaloshinskii-Moriya (DM) interaction [29][30][31] in a controllable way. However, the thickness of the sample is usually larger than the atomic lattice constant even in thin films and interfaces, and the degrees of freedom along the external magnetic field are still important. The distortion of the skyrmion string along the magnetic field is an important issue especially in the current-driven motion with impurities. The singular spin configurations associated with the creation and annihilation of skyrmions are called hedgehog or antihedgehog, corresponding to the monopole and antimonopole of emergent magnetic field, respectively [32][33][34][35][36][37] . Compared with the skyrmion dynamics in two-dimensions, the dynamics of skyrmion strings are less studied up to now especially during the current-driven motion. The mechanism for the creation of the monopole and antimonopole is an important issue, which is not easy to access experimentally but can be studied by numerical simulations.
In this paper, we investigate the current-driven dynamics of skyrmion with disorder by numerically solving the Landau-Lifshitz-Gilbert (LLG) equation for spins. To examine the essential role of the string degree of freedom of the skyrmion along the external magnetic field, we focus on the system thickness or string length dependence on its dynamics. By the numerical simulation of a skyrmion string under current, we find that the impurity effect strongly depends on the length of the skyrmion string. In particular, the theoretical results show that the impurity pinning is strongly suppressed by the unique dynamics of the skyrmion string. For a longer skyrmion string, we find its annihilation by the monopole-antimonopole creation through the excitation of the collective

Results
Model and simulation. The  i i with x, ŷ and ẑ being the unit vectors in the x-, y-and z-directions, respectively. The normalized magnetic moment at site = r r r r ( , , ) x y z is expressed by n r . The lattice constant a is taken as the unit of length. The second term represents the antisymmetric DM interaction. We for the model Hamiltonian of the chiral magnets, so that Bloch skyrmion dynamics is discussed. For the interface driven DM systems with thin film form, we use , so that Néel skyrmion dynamics is discussed. In the case that = D 0 Z , at least in this minimal model Hamiltonian, the dynamics of the magnetic texture is the same as that in the system with =D x D X and =D y D Y . The single spin anisotropy K imp is introduced on the random sites ∈ Λ r i (Λ: set of the random sites), which represents the disorder effect. In reality, the impurity potential has finite spatial extent. To mimic such impurity potential profile, we introduce dense impurities where the impurity-impurity distance is much smaller than skyrmion size. For the skyrmion, therefore, the impurities act as rather a spatially slowly varying potential rather than the point-like one. In the system without the disorder, i.e., = K 0 imp , the competition between the exchange interaction J and DM interaction D results in the single-q helix state with = q D J / , for = h 0. With increasing h, the skyrmion crystal state and the ferromagnetic state appears successively 38,39 . In this paper, we put =D z D Z in most cases and we use the parameter set = = .
where the metastable skyrmion exists in the ferromagnetic background. For the impurities, we introduce = . K 0 2 imp with 10% concentration. (see also Methods). The magnetic skyrmions discussed here is topologically the same as the magnetic bubble where the long-range magnetic dipolar interaction 40-43 plays essential role. Usually, such skyrmionic state (or nomal magnetic bubble) has several hundreds nanometer to sub-micrometer scale size. In those cases, the twist of magnetic moment alignment between Bloch and Néel winding textures occurs along the thickness direction for the static skyrmionic states 40,41,[44][45][46] . For such larger size magnetic skyrmion(s), however, different degree of freedom becomes important, i.e., the flexibility of the magnetic winding texture causes substantial deformation of the skyrmionic state during its dynamics. Those flexibilities induces complex skyrmion dynamics 42,43 . When the DM interaction dominates the winding magnetic texture, the skyrmion size becomes small, and it is typically several tens nanometer scale or smaller. In those cases, the magnetic bubble characteristics discussed above are reduced. Here, we discuss the limiting case, i.e., the DM interaction driven magnetic skyrmion and its dynamcs. Even in this simplified systems, skyrmion string shows rich dynamical behaviors.
The LLG equation is given by: where α is the Gilbert damping constant (see also Methods). The last two terms in Eq. (2) represent the spin-transfer-torque effect due to the spin polarlized electric current: j s represents the spin current density given by the spin polarlized electric current and β is the coefficient of the non-adiabatic effect. We use J for the unit of h, D, and K imp , and 1/(γJ) for the unit of time t. The unit of the spin current density j s = |j s | is 2eγJ/(pa 2 ) for the spin polarized electric current density (p: polarization of magnet). (See Methods for computational details.) For the typical B20 skyrmion magnets 16 , e.g., MnSi, FeGe and so on, J ~ 10 −3 eV, 1/(γJ) ~ 0.7 ps for (g s : electron spin g-factor, μ B : Bohr magneton), and 2eγJ/(pa 2 ) ~ 1.0 × 10 13 A/m 2 in the case that polarization of magnet p = 0.2 and the lattice constant a = 5 Å. Figure 1 summarizes the dynamical phase diagram of the current-driven skyrmion string with disorder. In this study, we introduce the impurity sites with 10% concentration. Here, to examine the skyrmion-string length dependence on the current driven dynamics, we prepare system with thickness L Z = 1 ~ 100 (see Methods). We employ the open boundary condition along z-direction while the periodic one for x-and y-directions.
To examine the dynamics of skyrmion string, we define the quantity, with γ (=x, y, z) component of the emergent electric field (e-field) 47-49 We also define the skyrmion number on the xy-plane at r z by www.nature.com/scientificreports www.nature.com/scientificreports/ In the present study, for the skyrmion string without monopole/antimonopole inside system, N sk (r z ) is independent of r z , so that we use the notation N sk (r z ) = N sk . For such skyrmion string without distortion, it is found that where is the velocity of the skyrmion string in horizontal direction. In the present case, N sk = −1.
Single-layer system (L Z = 1). First, we discuss the skyrmion dynamics in L Z = 1 system. Figure 2(a-c) show the time dependence of the quantities E x /j s (blue line) and E y /j s (red line), for j s = 0.006, 0.008, and 0.8, respectively. The results shown in Fig. 2(a) indicate the pinning dynamics of the skyrmion, i.e., the dumped oscillation around the (meta) stable point. This is similar to the relaxation dynamics of the object with vorticity in the harmonic potential on two-dimensions. We find such pinned state of a skyrmion for j s ≤ 0.006. For j s = 0.008 which is larger than the critical depinning current density for this single-layer system .
, the skyrmion starts to move but immediately evaporates. As seen in Fig. 2(b), E x /j s and E y /j s have peaky structure at t ~ 580 by the drastic change in magnetic texture due to the skyrmion annihilation. Figure 2(d-f) show the snapshots of the magnetic texture at t = 0, 578, and 586, respectively. In the plot Fig. 2 This represents that the impurity effect is reduced enough by the large j s . For 0.8 ≤ j s , the motional narrowing effect reduces the effect of impurities, and the current driven skyrmion survives.
Multilayer systems (L Z ≥ 2). In L Z = 2 system, we find the substantial suppression of the impurity effect, i.e., the reduction of the skyrmion depinning critical current density j c , and the range of the current density where the skyrmion annihilates is reduced (see Fig. 1). Figure 3 shows the skyrmion dynamics in two-layer system. For j s = 0.006 (=j s,c,1 ), we see that the skyrmion is depinned and moves although the skyrmion is pinned in the single-layer system as shown in Fig. 2(a). For j s < 0.004, the skyrmion is pinned as in the case of the single-layer system. In the two-layer system, we find the moving skyrmion depinned state between its pinned and annihilation regions, while this behavior is hardly seen in the single-layer system with the strong disorder as studied in the present paper.
With increasing number of stacking layer (or thickness of the system) L Z , the depinning critical current density j s,c is further reduced and skyrmion annihilation region disappears at L Z ~ 10. As seen in the Fig. 1, at L Z ~ 40, we do not find skyrmion annihilation region of j s .
For L Z ≥ 60, we find that the skyrmion annihilation occurs again for 0.06 ≤ j s ≤ 0.1. The dynamics without (with) skyrmion annihilation for L Z = 100 are shown in Figs 4 and 5. Figure 5 summarizes an example of the skyrmion annihilation dynamics. In this case, the skyrmion string moves for a while, and it is broken at a point, in other words, monopole-antimonopole pair is created. The monoplole and antimonopole run through the string, and finally the (metastable) skyrmion string disappears. This skyrmion annihilation behavior is different from that in thin systems. Figure 6 shows the magnetic texture of the horizontal cross section of the system at r z = r m-am where the pair creation of monopole-antimonopole of the skyrmion string occurs. In the time evolution of the magnetic texture, Fig. 6(a-d), we find the oscillation of the skyrmion size. Along the current driven dynamics, the www.nature.com/scientificreports www.nature.com/scientificreports/ breathing mode of the skyrmion is excited due to the collision with the impurities. The characteristic collective mode propagates along the skyrmion string and shows a resonant behavior. Finally, the breathing in skyrmion size is enhanced and the part of the string is squeezed enough, so that the annihilation of the skyrmion, in other words, the monopole-antimonopole pair creation occurs (see Figs 5 and 6). The breathing dynamics of the skyrmion on the cross section of the skyrmion string at hight r z is represented by the time dependence of the average of z component of n r , i.e., z z z r r r x y X Y ,( , , ) x y z where L X × L Y is the horizontal cross-sectional area of the system. Figure 7(a) shows n z (t, r z = r m-am ) for the results shown in Fig. 6: At around  t 1650, n z (t, r z = r m-am ) drops and reaches smallest value during this dynamics. After that, in the time duration 1650 ≤ t ≤ 1750, n z (t, r z = r m-am ) grows rapidly and saturates → 1. This represents the overshoot behavior of the skyrmion breathing dynamics, i.e., the skyrmion is expanding largely first, and later it is squeezing enough to annihilate the skyrmion. On the skyrmion dynamics, the spin wave excitations 50,51 on the skyrmion string play important role. Figure 7 To show the magnetic texture clearly, impurity sites are not indicated. The brightness of the color represents the z-component of the magnetic moment n z , i.e., black is corresponding to n z = −1 but we do not put color for n z > 0.5, i.e., the outside of the skyrmion string, to make the 3-dimensional magnetic texture clear. www.nature.com/scientificreports www.nature.com/scientificreports/ for the time duration T (=1600 is used for Fig. 7(b)) before the monopole-antimonopole pair creation. This represents the collective mode dynamics along skyrmion string towards the monopole-antimonopole creation. The gray curve in Fig. 7(b) indicates the spin-wave dispersion relation of this system, x y z z for (0, 0, q z ) but an easy-axis anisotropy K = K imp × 0.1 is applied for every sites. Because of the DM interaction D Z (=ẑ D ) between the sites along ẑ direction, the spin-wave excitation is not symmetric in q z axis. It is indicated that the breathing collective mode is below the continuum of the spin waves for the positive q z while it is buried in the negative direction. This results in the large directional propagation of the breathing mode. For 0.1 < j s , we do not see the skyrmion annihilation, i.e., the spin-transfer-torque effect by j is large enough compare to the impurity effect. Therefore, the moving skyrmion string becomes rather straight compare to that for weak j s < 0.06 (see also Fig. 4(a,b)).
The impurity sites are given by the random number, so that the impurity configuration is layer dependent. This layer dependence of the impurity configuration affects the impurity effect. In the system with layer independent impurity configuration, skyrmion dynamics with D Z = 0 is the same as that in single-layer system.
As discussed above, the annihilation of the skyrmion string is achieved by the monopole and antimonopole motion along the string. We now consider the propagation of monopole and antimonopole along the skyrmion string. Figure 8 shows the monopole and antimonopole motion in the system without impurity K imp and current j s (see also Methods). Figure 8(a) shows the initial state, i.e., the (metastable) skyrmion string. At the top surface,  The gray curve is the dispersion relation of the spin-wave for (0, 0, q z ) with easy axis anisotropy K imp × 0.1 for every sites (see text). To see the weak intensity clearly, contour plots of log |S(q, ω)| 2 are presented.
www.nature.com/scientificreports www.nature.com/scientificreports/ we apply strong h = 1.0 for t ≥ 0 to annihilate the skyrmion only on the top surface. By this procedure, the monopole is created, which begins to move downward in the system, since the skyrmion is the metastable state in this parameter set. Figure 8(b) shows the snapshot of the magnetic texture at t = 400. Here, we call monopole (antimonopole) for the magnetic texture with positive (negative) divn r inside the system. In the same way, we apply strong h = 1.0 on the bottom surface and examine the antimonopole dynamics (see Fig. 8(c)). The monopole (antimonopole) runs through the skyrmion string, and finally the string totally disappears. Figure 8(d) shows the monopole (antimonopole) position as a function of time (see also Methods). To compare the monopole velocity and the antimonopole one, the red broken line is also presented in Fig. 8(d).
We examine the α dependence of the monopole and antimonopole motion (see Fig. 9). With decreasing α, the velocity of the monopole/antimonopole is enhanced. This large directional dependence of the (anti)monopole motion will lead to the large nonreciprocal responses of the system, which is left for future studies. We also examine the case with D Z = 0 where the directional dependence disappears as expected.

Discussion
It has been experimentally shown 52,53 that the skyrmions in three-dimensional bulk form samples are driven by electric current with current density ~10 6 A/m 2 being orders of magnitude smaller than that for magnetic domain wall motion. This was assumed to be an advantageous for technological application in the earlier studies 16,54 . However, it is reported 55 that the microscale-fabricated MnSi sample shows very large threshold current density  www.nature.com/scientificreports www.nature.com/scientificreports/ 10 8 A/m 2 which is two orders magnitude larger than that in bulk samples. The experimental studies for the manipulation of the skyrmion have been extensively developed after the discovery of the room-temperature skyrmion induced by interface DM interaction in artificial thin-film heterostructures [20][21][22][23][24][25][26][27][28] . Here, the spintronics friendly centrosymmetric materials are used, where D X and D Y are induced by the interface 56 while D Z = 0. In those systems, the threshold current density to drive skyrmion(s) is ~10 10 -10 11 A/m 2 and is 4 or 5 orders of magnitude larger than that of conventional bulk skyrmion magnets. In view of present numerical simulation study for the current driven skyrmion, the threshold current density is not only determined by impurity strength and concentration but also by thickness of the system, i.e., with decreasing thickness of the system, the threshold current density is several orders of magnitude enhanced. In addition to this, such two-dimensionality is rather enhanced in the thin system with D Z = 0 where the twist of magnetic moments in z direction by DM interaction is missing. For longer skyrmion string, although the threshold current density is reduced, the monopole-antimonopole creation instability appears. In other words, for the stability of the skyrmion, there exists optimum sample-thickness which is typically ~100 nm for lattice constant a~5Å.
In our numerical study, for the metastable skyrmion under the strong magnetic field h, the velocity of monopole/antimonopole (along skyrmion string) is estimated to be ~10 2 m/sec for typical parameters (see the lines below Eq. (2)) of bulk skyrmion materials. The velocity is not much faster than the velocity of current driven skyrmion(s) and is also getting slow for smaller h. Our numerical study also suggests that in bulk samples, the impurity effect on current driven skyrmion is suppressed and skyrmion annihilation instability with monopole/ antimonopole creation appears. The monopole/antimonopole motion causes emergent e-field. The experiments (e.g. ref. 55 ) to measure the characteristics of e-field responsible for the skyrmion string dynamics will also be useful to observe the monopole/antimonopole dynamics. Such experiments are desired.

Conclusion
We numerically studied the current driven dynamics of a magnetic skyrmion string in the film form system in the presence of pinning effect. It was shown that the pinning effect is strongly enhanced for shorter skyrmion string which is corresponding to the skyrmion dynamics in thinner samples. For longer skyrmion string, it was also found that the skyrmion string annihilation instability due to the monopole-antimonopole pair creation. For the instability, it was found that the collective breathing mode on the skyrmion string plays a crucial role. The theoretical results suggest that there is an optimum sample thickness for the stability of the current driven skyrmion string.

Methods
For the LLG simulation of the single skyrmion string dynamics, we use 60 × 30 × L Z (L Z = 1 ~ 100) finite size system with periodic (open) boundary condition in horizontal (vertical) direction. For the initial state, we first prepare the metastable at K imp = 0, and obtain the relaxed state. To obtain the phase diagram Fig. 1 The numerical condition α = β = 0.01 is used for the numerical simulations except for the results summarized in Figs 8 and 9. In the condition α = β, the Hall angle of current driven skyrmion is zero for K imp = 0, so that we can examine purely impurity driven scattering effect.
For the impurity sites ∈ Λ r i (Λ: set of the random sites), we use the random-number generator developed by M. Matsumoto and T. Nishimura (http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html).
The center-of-mass of the skyrmion at r z is calculated by , and plot the trajectory in Figs 2 and 3. The discontinuity of N sk (r z ) as function of t and r z is used to identify the broken point of the skyrmion string r m-am and monopole/antimonopole position in z axis.