Controlling skyrmion helicity via engineered Dzyaloshinskii–Moriya interactions

Skyrmions, swirling spin textures appearing in chiral magnets, are magnetic structures in which the spins point in all the directions wrapping a sphere [1, 2]. These topological spin textures have ignited a growing interest in spintronics [3] due to their rich phenomenology as well as novel potential applications. Their nano-scale size, topologically-protected stability [2], and the very low electric current densities needed to displace them [4] are among the best qualities that make them attractive candidates for information carriers in high-density data-storage technologies [5]. They have been explored in bulk magnets like MnSi [4, 6, 7], in thin films of Fe_1−xCo_xSi [8–11], Mn_1−xFe_xGe [12], FeGe [13], La_0.5Ba_0.5MnO₃ [14], CuOSeO₃ [15], and also in a Fe monoatomic layer on Ir(1 1 1) [16]. Controlling the temperature and the external magnetic field applied to the samples, spontaneous skyrmion phases have been realized and detected by neutron scattering [6], Lorentz transmission electron microscopy (LTEM) [11], and spin-resolved scanning tunneling microscopy [16] experiments.

Current-driven skyrmion dynamics also displays intriguing topological transport properties [7, 17–19], phenomena that result from the so-called spin-transfer torques exerted by carrier spins on the magnetization [20]. Other mechanisms such as thermal gradients [21], inhomogeneity in the fields [5, 22, 23], and magnon currents [24] have also been proposed to induce the motion of skyrmions. In these studies it was assumed that the evolution of the skyrmion magnetization could be appropriately treated using a reduced set of collective coordinates. Under this approach the Landau–Lifshitz–Gilbert (LLG) equation is mapped to a particle-like equation of motion, also known as Thiele's equation [25].

The spin texture of magnetic skyrmions can originate from various mechanisms. Depending on the particular magnetic system, skyrmions can be stabilized by long-ranged dipolar, Dzyaloshinskii–Moriya, frustrated exchange, or four-spin exchange interactions [3]. In chiral magnets—the focus of this paper—spin orbit coupling and lack of inversion symmetry cause the appearance of the Dzyaloshinskii–Moriya interaction (DMI) [26]. This interaction enters in the continuum free energy as terms linear in spatial derivatives of the magnetization which are responsible for the stability of static skyrmions in two and three dimensions [27]. In contrast, static skyrmions of the Skyrme model of particle physics are stabilized by an interaction term of fourth order in spatial derivatives [1]. Furthermore, a class of exact time-periodic skyrmions was discovered for the Skyrme model which, for a specific frequency value, admits a nonlinear superposition law [28]. It should be noted that in magnetic systems, single skyrmions as well as the skyrmion crystal can become the lowest energy configuration by adjusting temperature and the applied magnetic field.

The DMI not only plays a fundamental role in the nucleation and stability of skyrmions, but also in the helical and conical phases [2]. For instance, properties like the size and helicity of skyrmions are set by the magnitude and sign of the DMI coupling, respectively. It has been shown that single skyrmions can be nucleated by injection of a spin-polarized current [5, 16, 29] or local heating [30] on magnetic thin-films. Magnetic skyrmions may emerge as stable vortex- or hedgehog-type spin configurations [2] depending on the type of DMI, which could be due to bulk interactions or induced at the interface of magnetic films in contact with heavy metals [5]. Vortex-type (or Bloch-type) skyrmions, those for which the magnetization swirls around their center and perpendicular to their radial direction, have been observed in non- centrosymmetric magnets with B20-type crystal structure such as Mn_1−xFe_xGe [12]. Meanwhile, hedgehog-type (or Néel-type) skyrmions, whose spins point either radially outward or inward, have been observed in magnetic ultra-thin films of Fe/Ir(1 1 1) [16] or in the recently reported polar magnetic semiconductor GaV₄S₈ [31].

Recent experiments have revealed the dependence of the DMI in the chiral magnet alloys Mn_1−xFe_xGe [12, 32] and Fe_1−xCo_xSi [33–35] on the chemical composition x. For instance, measuring the wave vector of the helical spin texture, it has been shown [34] that the strength and even the sign of the DMI can be controlled by varying the composition of Co in Fe_1−xCo_xSi. Similar results were reported in [12], where by changing the composition x in Mn_1−xFe_xGe the size and helicity of magnetic skyrmions were modified. These experiments, together with the latest theoretical studies [36–38], have been devoted to elucidate the intriguing relation of the magnitude and sign of the DMI to chemical composition, thus opening the possibility to manipulate the skyrmion helicity by engineering the crystal chirality and spin–orbit coupling. Motivated by these experiments [12, 32–35], here we study the dependence of the helicity degree of freedom on single magnetic skyrmion dynamics in chiral magnets with an inhomogeneous DMI coupling.

The structure of this paper is the following. In section 2, an effective description of the skyrmion dynamics based on the collective coordinates method is developed assuming a general inhomogeneous DMI coupling. The generalized Thiele's equation thus obtained, an autonomous dynamical system for the skyrmion center of mass position and helicity, is then further reduced to the case of a DMI coupling with arbitrary one-dimensional spatial dependence. Two particular cases are then analyzed taking advantage of the general results derived in this section. Sections 3 and 4 dwell upon half-planes with opposite-sign DMI coupling and linearly varying DMI coupling, respectively. The role of the helicity and the set of possible trajectories allowed by the dynamical systems describing each case are determined. Finally, in section 5, these two cases are compared, possible measurements that could be performed to confirm our predictions are briefly discussed, and our conclusions are presented.

The time evolution of any magnetic texture, such as a skyrmion, can be determined by solving the field equations for the magnetization field $\boldsymbol{M}$. However, in the systems we address in this paper, not all of the uncountably infinite degrees of freedom in the magnetization field are relevant. For a static, axially-symmetric skyrmion centered at the origin, the direction of the magnetization, $\mathbf{\Omega}\equiv \boldsymbol{M}/M$, can be written as

$$\begin{eqnarray}&&\begin{array}{*{20}{l}}\mathbf{\Omega}\left(\boldsymbol{r},\gamma \right)=\sin \Theta \left(r,\gamma \right)\left[\cos \left(\varphi +\gamma \right)\hat{\boldsymbol{x}}+\sin \left(\varphi +\gamma \right)\hat{\boldsymbol{y}}\right]\\ +\cos \Theta \left(r,\gamma \right)\hat{\boldsymbol{z}},\end{array}\end{eqnarray} \tag{ 1 }$$

where r and φ are the usual cylindrical coordinates. It then follows that to completely specify this magnetic texture, only its helicity γ and radial profile $\Theta \left(r,\gamma \right)$ are required. At low temperatures and in the presence of an external magnetic field $\boldsymbol{B}$, these defining quantities can be obtained by minimizing the corresponding magnetic energy functional $\mathcal{U}$ which, up to an overall constant, is given by

$$\begin{eqnarray}&&\mathcal{U}={{M}^{2}}{\int}^{}\text{d}\boldsymbol{r}[J{{(\nabla \mathbf{\Omega})}^{2}}+2D(\boldsymbol{r})\mathbf{\Omega}\cdot \nabla \times \mathbf{\Omega}-\bar{\boldsymbol{B}}\cdot \mathbf{\Omega}],\end{eqnarray} \tag{ 2 }$$

where J denotes the spin stiffness and $\bar{\boldsymbol{B}}=\boldsymbol{B}/M$. In all the results presented here $\boldsymbol{B}=B\hat{\boldsymbol{z}}$, i.e. orthogonal to the sample. This model assumes the magnitude of the magnetization field is uniform across the sample, it neglects anisotropy energies, and allows for a spatially-dependent DMI coupling $D(\boldsymbol{r})$. It has been successfully used to describe skyrmions in chiral magnet thin films with uniform DMI coupling, i.e. when $D(\boldsymbol{r})=D$ [11].

The physical scenario we wish to study begins with a skyrmion already nucleated in a region of the sample with a locally uniform DMI coupling. Although this magnetic texture may get distorted as it evolves, owing to its topologically-protected stability, it will remain being a skyrmion. Furthermore, for low energy distortions, which correspond to long-term dynamics, the axial symmetry is retained [39]. Therefore, we expect the skyrmion to translate across the sample and, since the helicity is inextricably related to the DMI, its helicity is also expected to evolve in time as the skyrmion explores regions with varying DMI. The previous analysis and assumptions justify an effective description of the skyrmion dynamics in terms of the position of its center $\boldsymbol{r}_{{0}}=X\hat{\boldsymbol{x}}+Y\hat{\boldsymbol{y}}$ and its helicity γ.

2.1. Collective coordinates and generalized Thiele's equation

Now that we have identified the relevant degrees of freedom for the effective long-term skyrmion dynamics, we need to determine their equations of motion. To that end we will use the method of collective coordinates [40] and closely following [41] we will construct the corresponding generalized Thiele's equation.

The method of collective coordinates is based on the assumption that the magnetic texture can be described by a discrete set of time-dependent coordinates $\boldsymbol{\xi}{{(t)}^{T}}=\left({{\xi}_{1}}(t),\ldots,{{\xi}_{N}}(t)\right)$, so that the time dependence of the magnetization comes solely from them, i.e. $\boldsymbol{M}\left(\boldsymbol{r},t\right)=\boldsymbol{M}\left(\boldsymbol{r},\boldsymbol{\xi}(t)\right)$. For our effective skyrmion dynamics, as discussed above, we have three collective coordinates: ${{\xi}_{1}}(t)=X(t)$, ${{\xi}_{2}}(t)=Y(t)$, and ${{\xi}_{3}}(t)=\gamma (t)$. When the magnitude of the magnetization can be regarded as constant, the equations of motion for the collective coordinates follow from the dynamics of the direction of the magnetization $\mathbf{\Omega}$, governed by the Landau–Lifshitz–Gilbert (LLG) equation

$$\begin{eqnarray}&&\frac{\text{d}\mathbf{\Omega}}{\text{d}t}=\bar{\gamma}{{\mathbf{H}}_{\text{eff}}}\times \mathbf{\Omega}+\alpha \mathbf{\Omega}\times \frac{\text{d}\mathbf{\Omega}}{\text{d}t}+\boldsymbol{\tau}_{{\text{STT}}}+\boldsymbol{\tau}_{{\text{SOT}}}.\end{eqnarray} \tag{ 3 }$$

Here $\bar{\gamma}=g|e|/\left(2{{m}_{e}}\right)$ is the electron gyromagnetic ratio, the effective magnetic field is ${{\mathbf{H}}_{\text{eff}}}=-\frac{1}{M}\frac{\delta \mathcal{U}}{\delta \mathbf{\Omega}}$, and α is the Gilbert damping constant. For completeness, the conventional spin-transfer torque $\boldsymbol{\tau}_{{\text{STT}}}$, due to the transference of spin angular momentum from spin-polarized currents to the magnetization, and a spin–orbit torque $\boldsymbol{\tau}_{{\text{SOT}}}$ [42], due to the intrinsic spin–orbit coupling of the chiral magnet, have also been included. Using that the time dependence of $\mathbf{\Omega}$ enters through the collective coordinates, i.e. $\mathbf{\Omega}\left(\boldsymbol{r},t\right)=\mathbf{\Omega}\left(\boldsymbol{r},\boldsymbol{\xi}(t)\right)$, the left-hand side of this equation can be rewritten as $\text{d}\mathbf{\Omega}/\text{d}t=\dot{{{{\xi}_{j}}}}\partial \mathbf{\Omega}/\partial {{\xi}_{j}}$. Acting on the LLG equation with $\mathbf{\Omega}\times$, taking the dot product with $\partial \mathbf{\Omega}/\partial {{\xi}_{i}}$, multiplying by the angular momentum density ${\mathsf{J~}}\,=M/\bar{\gamma}$, and finally integrating over the volume, we arrive at the generalized Thiele's equation

$$\begin{eqnarray}&&{{G}_{ij}}{{\dot{{\xi}}}_{j}}-{{\mathcal{D}}_{ij}}{{\dot{{\xi}}}_{j}}+{{F}_{i}}+F_{i}^{s}=0,\end{eqnarray} \tag{ 4 }$$

where

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{G}_{ij}} & ={\mathsf{J}}{\int}^{}\text{d}\boldsymbol{r}\;\mathbf{\Omega}\cdot \left(\frac{\partial \mathbf{\Omega}}{\partial {{\xi}_{i}}}\times \frac{\partial \mathbf{\Omega}}{\partial {{\xi}_{j}}}\right), \end{array}\end{eqnarray} \tag{ 5 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\mathcal{D}}_{ij}} & =\alpha {\mathsf{J}}{\int}^{}\text{d}\boldsymbol{r}\left(\frac{\partial \mathbf{\Omega}}{\partial {{\xi}_{i}}}\cdot \frac{\partial \mathbf{\Omega}}{\partial {{\xi}_{j}}}\right), \end{array}\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{F}_{i}} & =-{\int}^{}\text{d}\boldsymbol{r}\frac{\delta \mathcal{U}}{\delta \mathbf{\Omega}}\cdot \frac{\partial \mathbf{\Omega}}{\partial {{\xi}_{i}}}=-\frac{\partial \mathcal{U}}{\partial {{\xi}_{i}}}, \end{array}\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} F_{i}^{s} & =-{\mathsf{J}}{\int}^{}\text{d}\boldsymbol{r}\frac{\partial \mathbf{\Omega}}{\partial {{\xi}_{i}}}\cdot \left[\mathbf{\Omega}\times \left(\boldsymbol{\tau}_{{\text{STT}}}+\boldsymbol{\tau}_{{\text{SOT}}}\right)\right]. \end{array}\end{eqnarray} \tag{ 8 }$$

Above, G_ij and ${{\mathcal{D}}_{ij}}$ denote the gyrotropic and dissipative tensor, respectively. $\boldsymbol{F}$ is called the generalized force. The spin force, $\boldsymbol{F}^{{s}}$, is determined by the two-dimensional electric current density and the specific form of $\boldsymbol{\tau}_{{\text{SOT}}}$. In the results presented here the electric current density is assumed to vanish, so this force will be neglected. Its role in skyrmion dynamics in systems with spatially inhomogeneous DMI will be discussed elsewhere [43]. By solving the generalized Thiele's equation, we obtain the time dependence of the collective coordinates necessary to describe the effective skyrmion dynamics. In order to do so, we first need to compute G_ij, ${{\mathcal{D}}_{ij}}$, and $\boldsymbol{F}$ which we proceed to do below.

2.2. Gyrotropic and dissipative tensors

As can be seen from their definition, while the gyrotropic tensor is antisymmetric, the dissipative tensor is symmetric. Moreover, these two tensors depend only on the collective coordinate-dependent magnetic texture. The magnetic energy $\mathcal{U}$ plays absolutely no role in their calculation. From our previous discussion, we will use the magnetic texture introduced above equation (1), but with its center and helicity now promoted to collective coordinate status: $\mathbf{\Omega}\left(\boldsymbol{r},t\right)=\mathbf{\Omega}\left(\boldsymbol{r}-\boldsymbol{r}_{{0}}(t),\gamma (t)\right)$.

We found that the non-zero components of the gyrotropic tensor are ${{G}_{12}}=-{{G}_{21}}=4\pi d{\mathsf{J}}$, with d the thickness of the magnetic film. Noting that G_ij is proportional to the topological charge in the ${{\xi}_{i}}{{\xi}_{j}}$ space, this result follows directly from the topological structure inherited from the parent magnetic texture $\mathbf{\Omega}\left(\boldsymbol{r},\gamma \right)$. Our calculations show that the dissipative tensor obeys ${{\mathcal{D}}_{ij}}=4\pi \alpha d{\mathsf{J}}{{\delta}_{ij}}{{\eta}_{j}}$, with ${{\eta}_{1}}={{\eta}_{2}}$ by the axial symmetry of the parent texture, and where

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\eta}_{1}}\left(\gamma \right) & =\frac{1}{4}{\int}_{0}^{\infty}r\text{d}r\left[{{\left(\frac{\partial \Theta \left(r,\gamma \right)}{\partial r}\right)}^{2}}+\frac{{{\sin}^{2}} \Theta \left(r,\gamma \right)}{{{r}^{2}}}\right], \end{array}\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\eta}_{3}}\left(\gamma \right) & =\frac{1}{2}{\int}_{0}^{\infty}r\text{d}r\left[{{\left(\frac{\partial \Theta \left(r,\gamma \right)}{\partial \gamma}\right)}^{2}}+{{\sin}^{2}} \Theta \left(r,\gamma \right)\right]. \end{array}\end{eqnarray} \tag{ 10 }$$

Clearly ${{\eta}_{1}},{{\eta}_{2}},{{\eta}_{3}}>0$, so ${{\mathcal{D}}_{ij}}$ is positive definite. The radial profile $\Theta \left(r,\gamma \right)$, required to compute the ${{\eta}_{i}}$'s, was obtained by numerically solving the corresponding Euler–Lagrange equation. Details of the calculation of the skyrmion radial profile are outlined in the appendix.

In figure 1 the skyrmion textures and profiles for different helicities and external magnetic fields are displayed. The magnetization texture is illustrated in panels (a)–(d), for $\gamma =\pi /2$, $\gamma =3{{\gamma}^{\ast}}$, $\gamma =1.5{{\gamma}^{\ast}}$ and $\gamma =3\pi /2$, respectively. A minimum value for the helicity, which is denoted by ${{\gamma}^{\ast}}$, has been introduced in our continuum model to account for the discrete nature of the magnetic skyrmion. Because the size of the skyrmion can not be less than the lattice constant, an elementary size for the skyrmion is considered and accordingly, a minimum value for the helicity is determined. This value is estimated for Fe_0.5Co_0.5Si [33] to be ${{\gamma}^{\ast}}=0.23$. Meanwhile, in (e) and (f), plots of helicity-dependent solutions of the adimensionalized profile equation, equation (A.1), for a single skyrmion in the presence of an adimensionalized uniform external magnetic field $b=J\bar{B}/2{{D}^{2}}$ are shown. In (e), for a helicity $\gamma =\pi /2$, it is observed that as the magnetic field increases, the size of the skyrmion core shrinks. A similar behavior holds for the skyrmion size, plotted in (f) for a magnetic field b  =  0.75, as the helicity decreases. The size of the skyrmion core is a smooth function of the helicity, being maximum when $\gamma =\pi /2$ (a) or $\gamma =3\pi /2$ (d) and minimum for $\gamma \in \left[0,{{\gamma}^{\ast}}\right]{\cup}^{}\left[\pi -{{\gamma}^{\ast}},\pi +{{\gamma}^{\ast}}\right]{\cup}^{}\left[2\pi -{{\gamma}^{\ast}},2\pi \right]$.

Zoom In Zoom Out Reset image size

The curves shown in figure 2 present the result obtained for the helicity-dependent parts of the non-zero components of the dissipative tensor, ${{\eta}_{1}}$ and ${{\eta}_{3}}$. The plots were performed for several values of the external magnetic field b. It is observed that as the magnetic field is decreased, the variation of the components of the dissipative tensor is enhanced. From equations (9) and (10), it is clear that this is caused by the increase in size of the skyrmion profile due to its strong dependence on the external magnetic field (see figure 1(f) and also figures 1(a)–(d)).

Zoom In Zoom Out Reset image size

2.3. Generalized force and equations of motion

In order to compute the generalized force we first need to evaluate the magnetic energy for the collective coordinate-dependent skyrmion texture $\mathbf{\Omega}\left(\boldsymbol{r}-\boldsymbol{r}_{{0}},\gamma \right)$. The magnetic energy thus obtained, now a function of the collective coordinates, reads $\mathcal{U}\left(\boldsymbol{r}_{{0}},\gamma \right)={{\mathcal{U}}_{\text{DM}}}\left(\boldsymbol{r}_{{0}},\gamma \right)+{{\mathcal{U}}_{\text{EX}-\text{Z}}}\left(\gamma \right)$, where ${{\mathcal{U}}_{\text{DM}}}\left(\boldsymbol{r}_{{0}},\gamma \right)=2d{{M}^{2}}\sin \gamma {\int}^{}{{\text{d}}^{2}}\boldsymbol{r}D\left(\boldsymbol{r}_{{0}}+\boldsymbol{r}\right)\,f\left(r,\gamma \right)$, ${{\mathcal{U}}_{\text{EX}-\text{Z}}}\left(\gamma \right)=$ $dJ{{M}^{2}}g\left(\gamma \right)$, and

$$\begin{eqnarray}\begin{array}{*{35}{l}} f\left(r,\gamma \right) & = & \frac{\partial \Theta \left(r,\gamma \right)}{\partial r}+\frac{\sin 2 \Theta \left(r,\gamma \right)}{2r}, \end{array}\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} g\left(\gamma \right) & = & 2\pi {\int}_{0}^{\infty}\text{d}rr\left[{{\left(\frac{\partial \Theta \left(r,\gamma \right)}{\partial r}\right)}^{2}}+\frac{{{\sin}^{2}} \Theta \left(r,\gamma \right)}{{{r}^{2}}}-\frac{{\bar{B}}}{J}\cos \Theta \left(r,\gamma \right)\right]. \end{array}\end{eqnarray} \tag{ 12 }$$

It is worth noting that while the exchange and Zeeman energy terms become helicity-dependent, ${{\mathcal{U}}_{\,\text{EX}-\text{Z}}}\left(\gamma \right)$, the DMI energy also acquires a dependence on the skyrmion position, ${{\mathcal{U}}_{\text{DM}}}\left(\boldsymbol{r}_{{0}},\gamma \right)$. This was expected since only the DMI coupling was position-dependent.

As outlined in section 2.1, the components of the generalized force follow from straightforward partial differentiation of $\mathcal{U}\,\left(\boldsymbol{r}_{{0}},\gamma \right)$. For the sake of clarity, hereafter we will use the collective coordinate symbols as subindices for the generalized force components instead of numbers. Thus, ${{F}_{\left\{1,2\right\}}}=-\partial \mathcal{U}\,/\partial \left\{X,Y\right\}={{F}_{\left\{X,Y\right\}}}$, and ${{F}_{3}}=-\partial \mathcal{U}\,/\partial \gamma ={{F}_{\gamma}}$. For a DMI coupling with arbitrary spatial dependence on the thin film plane, $D\left(\boldsymbol{r}\right)$, the generalized force components are given below

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{F}_{\left\{X,Y\right\}}} & = & -2d{{M}^{2}}\sin \gamma {\int}^{}{{\text{d}}^{2}}\boldsymbol{r}\frac{\partial D\left(\boldsymbol{r}_{{0}}+\boldsymbol{r}\right)}{\partial \left\{X,Y\right\}}f\left(r,\gamma \right), \end{array}\end{eqnarray} \tag{ 13 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{F}_{\gamma}} & = & -d{{M}^{2}}\frac{\partial}{\partial \gamma}\left[2\sin \gamma {\int}^{}{{\text{d}}^{2}}\boldsymbol{r}D\left(\boldsymbol{r}_{{0}}+\boldsymbol{r}\right)\,f\left(r,\gamma \right)+Jg\left(\gamma \right)\right]. \end{array}\end{eqnarray} \tag{ 14 }$$

Using the explicit expressions determined in the previous section for the gyrotropic and dissipative tensors, we can write the generalized Thiele's equation for our model as the following autonomous, third-order dynamical system

$$\begin{eqnarray}\begin{array}{*{35}{l}} \dot{{Y}}-\alpha {{\eta}_{1}}\left(\gamma \right)\dot{{X}}+\frac{1}{4\pi d{\mathsf{J}}}{{F}_{X}}\left(X,Y,\gamma \right) & = & 0, \end{array}\end{eqnarray} \tag{ 15 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} \dot{{X}}+\alpha {{\eta}_{1}}\left(\gamma \right)\dot{{Y}}-\frac{1}{4\pi d{\mathsf{J}}}{{F}_{Y}}\left(X,Y,\gamma \right) & = & 0, \end{array}\end{eqnarray} \tag{ 16 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} \alpha {{\eta}_{3}}\left(\gamma \right)\dot{{\gamma}}-\frac{1}{4\pi d{\mathsf{J}}}{{F}_{\gamma}}\left(X,Y,\gamma \right) & = & 0. \end{array}\end{eqnarray} \tag{ 17 }$$

Although expected when internal degrees of freedom associated with the change of shape of the skyrmion are incorporated, our equations do not contain an inertia term. Unlike the results presented in [44], where by integrating out the internal degrees of freedom an inertia term arises, the focus of our work is to explicitly track the helicity dynamics. If equation (17) were used to eliminate γ from equations (15) and (16), an inertia term, among others, should emerge.

Now that we have constructed the generalized Thiele's equation which describes the effective skyrmion dynamics in thin films with an arbitrary DMI coupling, we will focus exclusively on the special case of an engineered DMI coupling with one-dimensional spatial dependence. Besides being one of the simplest, we are interested in this particular case because skyrmion-supporting samples with this feature have already been reported. Without loss of generality, we choose the x-axis as the direction along which the DMI coupling varies, i.e. $D\left(\boldsymbol{r}\right)=D(x)$. This choice makes the magnetic energy, and consequently the generalized force, no longer dependent on Y. Our model can be further simplified by eliminating $\dot{{Y}}$ using (16)—now with F_Y  =  0—giving way to the following adimensionalized second-order dynamical system

$$\begin{eqnarray}\begin{array}{*{35}{l}} \frac{\text{d}\tilde{X}}{\text{d}\tau} & =\frac{\alpha {{{\tilde{\eta}}}_{1}}(\gamma )}{1+{{\alpha}^{2}}\tilde{\eta}_{1}^{2}(\gamma )}{{{\tilde{F}}}_{X}}(\tilde{X},\gamma ), \end{array}\end{eqnarray} \tag{ 18 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} \frac{\text{d}\gamma}{\text{d}\tau} & =\frac{1}{\alpha {{{\tilde{\eta}}}_{3}}(\gamma )}{{{\tilde{F}}}_{\gamma}}(\tilde{X},\gamma ), \end{array}\end{eqnarray} \tag{ 19 }$$

where all spatial variables and time have been scaled by R  =  J/D and $T=\pi {\mathsf{J}}R/(2D{{M}^{2}})$, respectively. Here we have also introduced $\tilde{X}=X/R$, ${{\tilde{\eta}}_{1}}={{\eta}_{1}}$, ${{\tilde{\eta}}_{3}}={{\eta}_{3}}/{{R}^{2}}$, ${{\tilde{F}}_{X}}={{F}_{X}}/(8dD{{M}^{2}})$, and ${{\tilde{F}}_{\gamma}}={{F}_{\gamma}}/(8dRD{{M}^{2}})$. Therefore, only $\tilde{X}$ and γ govern the effective skyrmion dynamics while $\tilde{Y}=Y/R$ is now a slave variable whose time evolution is determined from

$$\begin{eqnarray}&&\begin{array}{*{35}{l}} \frac{\text{d}\tilde{Y}}{\text{d}\tau}=-\tan \delta (\gamma )\frac{\text{d}\tilde{X}}{\text{d}\tau}, \end{array}\end{eqnarray} \tag{ 20 }$$

where the Hall angle $\delta (\gamma )$ satisfies $\tan \delta (\gamma )=1/[\alpha {{\tilde{\eta}}_{1}}(\gamma )]$. Because of the typical small values of the Gilbert damping, the Hall angle would attain an almost γ-independent value close to $\pi /2$, thus determining a large ratio between the Y- and X-velocities. However, if the variation of ${{\tilde{\eta}}_{1}}(\gamma )$ is made large enough, by decreasing the external magnetic field, it could be possible to observe the helicity dependence predicted for $\delta (\gamma )$.

The next two sections are devoted to the effective skyrmion dynamics in thin films, predicted by the model developed so far, in two particular scenarios of engineered DMI coupling. In both cases we will proceed to determine ${{\tilde{F}}_{X}}$ and ${{\tilde{F}}_{\gamma}}$, then the general equations (18)–(20) will be used to analyze the resulting dynamics.

Thin films of Mn_1−xFe_xGe with an engineered, spatially-varying DMI strength have been recently synthesized [12]. LTEM measurements revealed the existence of neighboring domains hosting skyrmions with opposite helicity: $\gamma =\pi /2$ and $\gamma =3\pi /2$, respectively. Skyrmion dynamics in such systems can be modeled by considering a thin film split into two half-planes with opposite-sign DMI couplings. If the interface between the two domains is located at x  =  0, the spatial dependence of the DMI coupling can be written as $D(\boldsymbol{r})=D\;\text{sgn}(x)$. Using equations (13) and (14), the adimensionalized components of the generalized force are calculated to be

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{{\tilde{F}}}_{X}}(\tilde{X},\gamma ) & = & -\sin \gamma {{{\tilde{f}}}_{1}}(\tilde{X},\gamma ), \end{array}\end{eqnarray} \tag{ 21 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\tilde{{F}}}_{\gamma}}(\tilde{{X}},\gamma ) & = & -\frac{\partial}{\partial \gamma}\left[\sin \gamma {{\tilde{{\,f}}}_{3}}(\tilde{{X}},\gamma )+\frac{1}{8}\tilde{{g}}(\gamma )\right], \end{array}\end{eqnarray} \tag{ 22 }$$

where

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{{\tilde{f}}}_{1}}(\tilde{X},\gamma ) & = & {\int}_{0}^{\infty}\text{d}v\tilde{f}(\sqrt{{{{\tilde{X}}}^{2}}+{{v}^{2}}},\gamma ), \end{array}\end{eqnarray} \tag{ 23 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{{\tilde{f}}}_{3}}(\tilde{X},\gamma ) & = & {\int}_{0}^{\infty}\text{d}v{\int}_{0}^{{\tilde{X}}}\text{d}u\tilde{f}(\sqrt{{{u}^{2}}+{{v}^{2}}},\gamma ), \end{array}\end{eqnarray} \tag{ 24 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} \tilde{{g}}\left(\gamma \right) & = & 2\pi {\int}_{0}^{\infty}\text{d}\rho \rho \left[{{\left(\frac{\partial \theta \left(\rho,\gamma \right)}{\partial \rho}\right)}^{2}}+\frac{{{\sin}^{2}}\theta \left(\rho,\gamma \right)}{{{\rho}^{2}}}-2b\cos \theta \left(\rho,\gamma \right)\right]. \end{array}\end{eqnarray} \tag{ 25 }$$

The dimensionless functions introduced above, $\tilde{g}\left(\gamma \right)=g\left(\gamma \right)$ and $\tilde{f}\left(\rho,\gamma \right)=Rf\left(r,\gamma \right)$, were obtained by scaling the radial coordinate as $\rho =r/R$, and then defining $\theta \left(\rho,\gamma \right)= \Theta \left(R\rho,\gamma \right)$.

In samples with uniform DMI coupling D, the static, axially-symmetric skyrmions that extremize the magnetic free energy in equation (2) have helicity $\pi /2$ or $3\pi /2$ for D  >  0 and D  <  0, respectively. A skyrmion with either of those helicities does not extremize the magnetic energy in the sample with nonuniform DMI coupling considered in this section, unless it is located at an infinite distance away from the interface. Therefore, the expected long-term dynamics of a skyrmion in this type of sample is to flow toward $(\tilde{X}=\infty,\gamma =\pi /2)$ or $(\tilde{X}=-\infty,\gamma =3\pi /2)$ depending on its initial condition. Indeed, as depicted in the phase portrait in figure 3(a), these are the only two attractors—stable fixed points—of the corresponding dynamical system. The unstable fixed points are $(\tilde{X}=0,\gamma =0,2\pi )$, $(\tilde{X}=0,\gamma =\pi )$. The separatrices, i.e. trajectories that separate regions of the phase portrait with qualitatively different solutions, are also shown, in red. A total of eight qualitatively different possible solutions emerge from the structure of fixed points and separatrices of this dynamical system. These correspond to the eight regions labeled in the phase portrait. A diagram sketching all types of possible trajectories in the XY plane is included in figure 3(b).

Zoom In Zoom Out Reset image size

According to the phase portrait of solutions, a skyrmion that has been nucleated to the right of the interface ($\tilde{X}>0$) could either remain within this domain and eventually flow away from the interface toward the attractor $(\tilde{X}=\infty,\gamma =\pi /2)$ or cross over to the left domain and flow asymptotically to $(\tilde{X}=-\infty,\gamma =3\pi /2)$. Among the possible trajectories that remain within the right domain, those whose initial condition is close to the interface could have the skyrmion approach the interface first and the turn around to finally flow toward the attractor. Interestingly, the time evolution of the helicity is strictly monotonic (clockwise/counterclockwise as in region 3/4 (figure 3(a))) for those skyrmions that remain to the right of the interface, whereas the helicity of skyrmions that cross over to the left domain can increase/decrease to a maximum/minimum value before decreasing/increasing to $3\pi /2$ at the attractor. A similar analysis can be made for skyrmions nucleated to the left of the interface ($\tilde{X}<0$).

Two features observed for skyrmions nucleated at large distances from the interface are worth noting. First, the further a skyrmion is nucleated from the interface, the smaller becomes its probability to cross over. This can be understood by observing how little regions 5–8 (figure 3(a)), which host solutions that cross over, such as that shown in figure 4(a), extend away from the interface. Second, for those initial conditions that do not result in the skyrmion crossing over (figure 4(b)), the time evolution appears to take place in two steps: a fast change in the helicity toward its asymptotic value with almost no change in the skyrmion position, followed by a slow flow of the skyrmion away from the interface.

Zoom In Zoom Out Reset image size

The skyrmion dynamics under a linear gradient of the DMI coupling is considered in this section. The spatially-dependent DMI coupling consists of a linear variation along the x-coordinate as $D(\boldsymbol{r})=D+{{D}_{0}}\,x/l$ for 0  <  x  <  l, connecting two regions at x  <  0 and x  >  l where $D(\boldsymbol{r})=D$ and $D(\boldsymbol{r})=D+{{D}_{0}}$, respectively. The coupling D is assumed positive in order to fix the skyrmion profile with helicity $\gamma =\pi /2$ far from the transition region of $D(\boldsymbol{r})$. As we previously mentioned in section 2, the spatial change in the DMI behaves as external forces acting on the helicity and center of mass position of the magnetic skyrmion. Thus, the generalized forces along the coordinates $\tilde{X}$ and γ explicitly read,

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\tilde{{F}}}_{X}}[\tilde{{X}},\gamma ] & =-\sin \gamma \frac{{{D}_{0}}}{4D}{{\tilde{{{\mathsf{I}}}}}_{1}}[\tilde{{X}},\gamma ], \end{array}\end{eqnarray} \tag{ 26 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\tilde{{F}}}_{\gamma}}[\tilde{{X}},\gamma] & =-\frac{1}{4}\frac{\partial}{\partial \gamma}\left[\tilde{{\mathcal{E}}}_{\text{DM}}^{0}(\gamma )+\sin \gamma \frac{{{D}_{0}}}{D}(\tilde{{X}}{{\tilde{{{\mathsf{I}}}}}_{1}}[\tilde{{X}},\gamma] \right. \\ {} & \left. \quad \quad +\,{{\tilde{{{\mathsf{I}}}}}_{2}}[\tilde{{X}},\gamma]+{{\tilde{{{\mathsf{I}}}}}_{3}}[\tilde{{X}},\gamma])+2\pi g(\gamma )\right], \end{array}\end{eqnarray} \tag{ 27 }$$

where the dimensionless integrals ${{\tilde{{\mathsf{I}}}}_{1}}$, ${{\tilde{{\mathsf{I}}}}_{2}}$, ${{\tilde{{\mathsf{I}}}}_{3}}$ and $\tilde{\mathcal{E}}_{\text{DM}}^{0}$ entering in the above expressions are defined as follows

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\tilde{{\mathsf{I}}}}_{1}}[\tilde{X},\gamma ] & =\frac{2}{\lambda}{\int}_{0}^{\infty}\text{d}{{{\tilde{y}}}^{\prime}}{\int}_{-\tilde{X}}^{\lambda -\tilde{X}}\text{d}{{{\tilde{x}}}^{\prime}}\tilde{f}\left({{\rho}^{\prime}},\gamma \right), \end{array}\end{eqnarray} \tag{ 28 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\tilde{{\mathsf{I}}}}_{2}}[\tilde{X},\gamma ] & =\frac{2}{\lambda}{\int}_{0}^{\infty}\text{d}{{{\tilde{y}}}^{\prime}}{\int}_{-\tilde{X}}^{\lambda -\tilde{X}}\text{d}{{{\tilde{x}}}^{\prime}}{{{\tilde{x}}}^{\prime}}\tilde{f}\left({{\rho}^{\prime}},\gamma \right), \end{array}\end{eqnarray} \tag{ 29 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} {{\tilde{{\mathsf{I}}}}_{3}}[\tilde{X},\gamma ] & =2{\int}_{0}^{\infty}\text{d}{{{\tilde{y}}}^{\prime}}{\int}_{\lambda -\tilde{X}}^{\infty}\text{d}{{{\tilde{x}}}^{\prime}}\tilde{f}\left({{\rho}^{\prime}},\gamma \right), \end{array}\end{eqnarray} \tag{ 30 }$$

$$\begin{eqnarray}\begin{array}{*{35}{l}} \tilde{{\mathcal{E}}}_{\text{DM}}^{0}\left(\gamma \right) & =2\pi \sin \gamma {\int}_{0}^{\infty}\text{d}{{\rho}^{\prime}}{{\rho}^{\prime}}\tilde{{f}}\left({{\rho}^{\prime}},\gamma \right), \end{array}\end{eqnarray} \tag{ 31 }$$

where $\tilde{f}\left(\rho,\gamma \right)$ was already introduced in section 3 and $\lambda =l/R$ is the ratio between the scale of variation of the DMI coupling and the relevant magnetic length scale.

Previous to explore the case devoted to this section, we will analyze first the uniform DMI coupling case, i.e. in the limit when the DMI gradient D₀ goes to zero. This assumption implies that the force ${{\tilde{F}}_{X}}=0$ and thus no motion of the center of mass of the skyrmion is expected. On the other hand, the force acting on the helicity turns out to derive from the potential energy $\mathcal{U}\,[{{D}_{0}}=0]=(2\pi g+\tilde{{\mathcal{E}}}_{\text{DM}}^{0})\,/4$, which is symmetric with a local minimum at $\gamma =\pi /2$. Therefore, the time evolution of the skyrmion helicity is driven to the equilibrium configuration at $\gamma =\pi /2$, for any given initial condition. The same argument holds for a negative DMI coupling, having the equilibrium skyrmion profile a helicity $\gamma =3\pi /2$. A similar behavior is observed for the case treated in section 3, which occurs far from the transition region located at $\tilde{X}=0$.

In figure 5(a) the phase portrait representation for the skyrmion dynamics under the linear gradient of DMI coupling is shown, (see figure 5(b)). The flow of trajectories obeying the generalized Thiele's equation are displayed as a function of the helicity γ and the center of mass position $\tilde{X}$. The computation was carried out by solving the set of integrals equations (28)–(31), for the parameter $\lambda =5$ and a skyrmion profile solution obtained for b  =  0.5. Alongside this, the generalized forces were determined for the parameters D₀ /D  =  5 and $\alpha =0.5$, whose behavior is displayed in figures 5(c) and (d). As was explained in section 3, at every single point on each trajectory in figure 5(a) there is a vector whose components represent the value for $(\text{d}\tilde{X}/\text{d}\tau,\text{d}\gamma /\text{d}\tau )$. Indeed, for a set of initial conditions dictated by $\gamma \ne 3\pi /2$ we see that all the trajectories converge to the solution for which the helicity is $\gamma =\pi /2$, represented by the right-directed red line. This is because D, as well as D₀, have been assumed positive and therefore, the skyrmion configuration for which the magnetic energy is extremized is that with helicity $\pi /2$. In particular, the motion for a helicity $\gamma =\pi /2$ (bottom red line) occurs in the rightward direction, i.e. opposite to the gradient of the DMI coupling. The above statement (${{\tilde{F}}_{X}}>0$) comes from the fact that the DMI energy ${{\mathcal{U}}_{\text{DM}}}$ is strictly negative, which is a consequence of $f(r,\gamma )<0$, see equation (11). Furthermore, the skyrmion preserves its size as it evolves along the trajectories at $\gamma =\pi /2$ and $\gamma =3\pi /2$. Indeed, this is because the computation of the skyrmion profile (for more details see the appendix) is accomplished as a function of the helicity and external magnetic field, but does not take into account the local variation of the DMI. Therefore, for those trajectories where the helicity is constant, there will be no change in the size or shape of the skyrmion.

Zoom In Zoom Out Reset image size

Finally, it is worth commenting about the connection between the linear case treated here and that described in section 3. The half-planes with opposite-sign DMI case can be obtained by taking an appropriate limit to the generalized forces ${{\tilde{F}}_{X}}$ and ${{\tilde{F}}_{\gamma}}$ described in this section. In fact, assuming that D  <  0, choosing that ${{D}_{0}}=-2|D|$ and taking the limit when $\lambda \to 0$, we make sure that the linear transition occurs between opposite values in a very narrow region. Under these considerations, the integral expressions for ${{{\mathsf{I}}}_{1}}$ and ${{{\mathsf{I}}}_{2}}$ turn out to obey ${{\tilde{{\mathsf{I}}}}_{1}}[\tilde{X},\gamma ]=2{{\tilde{f}}_{1}}(\tilde{X},\gamma )$ and ${{\tilde{{{\mathsf{I}}}}}_{2}}[\tilde{{X}},\gamma]=-2\tilde{{X}}{{\tilde{{\,f}}}_{1}}(\tilde{{X}},\gamma )$. Therefore, the results for ${{\tilde{F}}_{X}}$ and ${{\tilde{F}}_{\gamma}}$ given by equations (26) and (27) reduce to those defined in section 3.

We have investigated single skyrmion dynamics in chiral magnets with locally varying Dzyaloshinskii–Moriya interactions. An effective description, based on the collective coordinates method, was established in terms of the skyrmion center of mass and helicity, which were found to satisfy a generalized Thiele's equation. In the construction of this equation, only the static skyrmion radial profile, $\Theta (r,\gamma )$ was needed to compute the gyrotropic and the dissipative tensors. The spatially dependent DMI coupling only entered in the calculation of the generalized forces. This separation of dependencies was exploited to write a general dynamical system governing the skyrmion dynamics that could be used for the two cases discussed in the previous sections. Therefore, the dynamical system to be solved in each case was fully determined once the DMI coupling was specified so the generalized forces could be calculated. It should be noted that although the radius of the skyrmion was not included as a collective coordinate, by using a helicity-dependent skyrmion radial profile to derive the generalized Thiele's equation, the size of the skyrmion inherited its dynamics from the helicity. This is a general feature of our model independent of the particular DMI coupling considered, hence it was observed in the two studied cases. Another feature shared by both cases is the helicity dependence of the Hall angle. Since $\delta (\gamma )$ derives from the dissipative tensor, it will be the same for any form of the DMI coupling. As with the skyrmion radius, the Hall angle also inherits its time dependence from the helicity. Furthermore, for sufficiently small values of the Gilbert damping, the Hall angle becomes almost independent of the helicity with a value close to $\pi /2$.

As evidenced by their respective phase portraits, the dynamics of the case from section 3 is richer than that from section 4. Due to the choice of parameters in the case of linearly varying DMI, only one attractor was present and two types of qualitatively different trajectories, differing by the sense of rotation of the helicity, were allowed. On the other hand, by virtue of the change of sign in the DMI coupling, the other case had two attractors as well as eight qualitatively different possible trajectories.

The simple helicity dynamics allowed by the two cases we have investigated, should not be seen as mere academic exercises. On the contrary, the proposed engineered DMI coupling scenarios must be considered as the simplest testbeds where the nontrivial time evolution of the skyrmion helicity can be studied in a controlled fashion. Indeed, the features observed in the two cases discussed could be detected by the experimental technique used in [45] where the low-energy dynamics of a skyrmion lattice was investigated measuring its microwave response.

In this work, we predicted how a single magnetic skyrmion texture adjusts its helicity as it propagates under a spatially modulated DMI. The results presented here establish the road to single magnetic skyrmion manipulation with full control over all its relevant degrees of freedom. Moreover, we expect future applications to take advantage of the relation between the skyrmion position and helicity using the structure of their phase portrait as a valuable guide. Finally, our findings constitute the first step toward a more complete understanding of the physics of the skyrmion helicity, which opens numerous opportunities for an effective control of skyrmions in future high-density storage and logic devices. Nucleation of skyrmions in chiral magnets with inhomogeneous DMI as well as the role of electric currents on their dynamics, will be addressed elsewhere.

It is a pleasure to thank AS Nuñez and RA Duine for fruitful and stimulating discussions. RET acknowledges support by Proyecto Fondecyt Postdoctorado number 3150372. SD acknowledges partial support from the International Fulbright Science and Technology Award.

The radial profile of a static, axially symmetric skyrmion with helicity γ, $\Theta \left(r,\gamma \right)$, was at the core of our calculations for the terms in the generalized Thiele's equation. Assuming the skyrmion was nucleated in a region with a locally uniform DMI coupling, the magnetic energy equation (2), with $D\left(\boldsymbol{r}\right)=D$, was evaluated in $\boldsymbol{M}=M\mathbf{\Omega}\left(\boldsymbol{r},\gamma \right)$ with $\mathbf{\Omega}\left(\boldsymbol{r},\gamma \right)$ as in equation (1). Thus, a functional in the skyrmion radial profile was obtained. Extremizing this functional led to the following Euler–Lagrange equation

$$\begin{eqnarray}&&\begin{array}{*{35}{l}} \frac{{{\text{d}}^{2}}\theta}{\text{d}{{\rho}^{2}}}+\frac{1}{\rho}\frac{\text{d}\theta}{\text{d}\rho}+\sin \gamma \frac{2{{\sin}^{2}}\theta}{\rho}-\frac{\sin 2\theta}{2{{\rho}^{2}}}-b\sin \theta =0, \end{array}\end{eqnarray} \tag{ A.1 }$$

with boundary conditions $\theta \left(0,\gamma \right)=\pi$ and $\theta \left(\infty,\gamma \right)=0$. The above ODE has been adimensionalized by scaling the radial coordinate as $\rho =r/R$, where R  =  J/D, as well as defining $\theta \left(\rho,\gamma \right)= \Theta \left(R\rho,\gamma \right)$ and $b=\frac{\bar{B}{{R}^{2}}}{2J}$. It is worth pointing out that the helicity enters in equation (A.1) as a parameter, in the same way as the adimensionalized external magnetic field b. Only after solving this equation numerically, using the shooting method, γ was promoted to collective coordinate status.

The skyrmion radius is defined as that for which $\theta =\frac{\pi}{2}$. To account for the underlying lattice structure, we do not allow skyrmions with a radius smaller than ${{\rho}^{\ast}}=\sqrt{2}a/R$, where a is the lattice constant. This restriction is enforced demanding that for $\gamma \in \left[0,{{\gamma}^{\ast}}\right]{\cup}^{}\left[\pi -{{\gamma}^{\ast}},\pi +{{\gamma}^{\ast}}\right]{\cup}^{}\left[2\pi -{{\gamma}^{\ast}},2\pi \right]$, $\theta \left(\rho,\gamma \right)=\theta \left(\rho,{{\gamma}^{\ast}}\right)$, where ${{\gamma}^{\ast}}$ is such that $\theta \left({{\rho}^{\ast}},{{\gamma}^{\ast}}\right)=\frac{\pi}{2}$.

Schematics of the skyrmion magnetization field obtained are depicted in figure 1 for several values of the helicity (panels (a)–(d)). The dependence of the skyrmion radial profile on the external magnetic field and the helicity are shown in figures 1(e) and (f), respectively.