Traveling skyrmions in chiral antiferromagnets

Skyrmions in antiferromagnetic (AFM) materials with the Dzyaloshinskii-Moriya (DM) interaction are expected to exist for essentially the same reasons as in DM ferromagnets (FM). It is shown that skyrmions in antiferromagnets with the DM interaction can be traveling as solitary waves with velocities up to a maximum value and their configuration is calculated in detail. The energy and the linear momentum of an AFM skyrmion are calculated and these lead to a proper definition of its mass. We give the energy-momentum dispersion of traveling skyrmions and explore their particle-like character. The skyrmion number, known to be linked to the dynamics of topological solitons in FM, is, here, unrelated to the dynamical behavior. As a result, the solitonic behavior of skyrmions in AFM is in stark contrast to the dynamical behavior of their FM counterparts.


I. INTRODUCTION
Topological magnetic solitons have been studied extensively for ferromagnets (FM) and weak ferromagnets (WFM). In both cases a nonvanishing magnetization develops in the ground state, albeit by a different physical mechanism, which allows a detailed experimental investigation by standard techniques [1,2]. In contrast, direct experimental evidence for pure antiferromagnetic (AFM) solitons is rare and has only been reported in recent years [3,4]. Theoretical arguments suggest that such solitons should exist for essentially the same reasons as in ordinary FM. It is therefore expected that vortices and chiral skyrmions that have been observed and studied in FM with the Dzyaloshinskii-Moriya (DM) interaction [5] have their counterparts in AFM [6]. On the other hand, the dynamics of solitons in AFM is expected to be substantially different [7].
The dynamics of the magnetic microstructure in AFM is governed by suitable extensions of the relativistic nonlinear σ model [8] instead of the LandauLifshitz equation in FM. The relevance of the σ model for the description of antiferromagnets became apparent through standard hydrodynamic approaches [9][10][11] but detailed applications to AFM solitons were carried out mostly in the Soviet literature reviewed in part in [2]. The type of dynamics of magnetic solitons supported by the σ model allows for traveling solitons and it is thus very different than the dynamics within the Landau-Lifshitz equation [12]. A topological soliton in a FM, such as a chiral skyrmion, is characterized by a topological number, called the skyrmion number in this context. A direct link between the skyrmion number and the dynamics of topological solitons in FM [12] was already apparent in the, so-called, Thiele equation for rigid vortex motion [13]. In contract to FM, the skyrmion number is not linked to the dynamics of AFM solitons; instead, a different topological number was shown to be linked to AFM soliton dynamics when an external magnetic field is applied [8].
The above remarks paint an intriguing picture for the dynamics of topological solitons in AFM. In this work, we will focus on AFM materials with the DM interaction, such as those studied in Ref. [14] and we will study the dynamics of chiral skyrmions. The existing results leave open the possibility for driving skyrmions in AFM as ordinary Newtonian particles, without experiencing a skew deflection (or Magnus force dynamics) seen for topological solitons in FM. In fact, traveling solitons can be readily found within the standard σ model by means of a Lorentz transformation, but the issue remains to be studied within an extension of the σ model for chiral magnets.
We show that topological solitons, such as skyrmions, present solitary wave behavior and they can propagate with a velocity up to a maximum value that depends on the DM parameter. We calculate the details of the traveling skyrmion configuration and we show that a rest mass can be naturally associated to a skyrmion. The particle-like character of AFM skyrmions is shown via their dispersion relation. Although our results are obtained within the conservative σ model the details can help guide any subsequent efforts to exploit the dynamics of skyrmions by applying external forces, such as spintorques.
The outline of the paper is as follows. Sec. II introduces the discrete model for spin dynamics in AFM and a continuous theory derived from this, that is, a nonlinear σ model. Sec. III presents numerical solutions for traveling AFM skyrmions within the σ model and the discrete model. Sec. IV contains our concluding remarks. In Appendix A we five details of the derivation of the σ model. In Appendix B we derive virial relations for traveling skyrmions.

II. THE NONLINEAR σ MODEL
A. The discrete model As a model for the magnet, we consider a square lattice of spins S i,j with a fixed length S 2 i,j = s 2 , where i, j are integer indices for the spin site. Magnetic materials with crystal structure of low symmetry present exchange interactions with both a symmetric and an antisymmetric part and the latter is usually called the Dzyaloshinskii-Moriya interaction [15,16]. We write a discrete Hamiltonian on the square lattice and we include symmetric exchange, a DM term, and an anisotropy term, The symmetric part of the exchange energy is where antiferromagnetic coupling has been assumed. For the DM interaction we will follow Ref. [14] and consider the material K 2 V 4 O 8 . We will use only a simplified version of the DM energy and not the full form allowed by the symmetry of the crystal. We set (3) whereê i , i = 1, 2, 3 denote the unit vectors in spin space. We consider an anisotropy term of the easy-axis type: where (S i,j ) 3 denotes the third component of a spin vector. The equation of motion for the spins is derived from the Hamiltonian and reads where F is the effective field. The first term on the right side of the equation conserves the energy and the second one is a damping term withα the dissipation constant. The explicit form of the effective field is

B. The continuum approximation
Further analysis will be greatly facilitated if we pass to a continuum model for the antiferromagnet. For the derivation of a continuum model we first need to define an appropriate order parameter with a continuum limit. In order to do this, we consider a tetramerization of the square lattice as shown in Fig. 1 and each tetramer is labelled by two indices α and β numbered consecutively in the horizontal and vertical directions respectively. At each tetramer, we denote the spin values by A α,β , B α,β , C α,β , D α,β thus defining four sublattices. A convenient set of fields is defined by the linear combinations, defined at each tetramer, The vector m gives the normalised magnetization at each tetramer and n is called the Néel vector. We consider a small parameter in terms of which the cartesian coordinates are where α 0 , β 0 are constants defining the central point of the lattice of tetramers. As → 0 the coordinates in Eq. (8) become continuous variables. In the same limit, we assume that the fields A α,β , B α,β , C α,β , D α,β and also those in Eq. (7) approach continuous limits. The derivation of the model is given in Appendix A. The continuous field n = n(x, y, τ ) with components (n 1 , n 2 , n 3 ) satisfies the constraint n 2 = 1, with the scaled time variable τ defined in Eq. (A2). In the conservative case (α = 0), it satisfies Eq. (A7), which is also written as where the dot denoted differentiation with respect to time τ , ∆ denotes the Laplacian in two dimensions, µν is the antisymmetric tensor with µ, ν = 1, 2, and the summation convention for repeated indices is adopted. The notation ∂ 1 , ∂ 2 denotes differentiation with respect to x, y respectively. Model (9) is an extension of the nonlinear σ model. It is Hamiltonian with energy (see also Refs. [6,17]) and the effective field in Eq. (9) is derived from f = −δV /δn. The fields m, k, l are auxiliary and they are given in terms of n by One should notice that the static sector of the σ model (9) for the Néel vector in an antiferromagnet is identical to the static sector of the Landau-Lifshitz equation for the magnetization vector of a ferromagnet with corresponding interactions (exchange, DM interaction, and anisotropy). We therefore expect that the static solitons (skyrmions, vortices, etc) obtained in an AFM precisely correspond to their counterparts in a FM. On the other hand, the dynamics of these solitons are expected to be dramatically different in AFM compared to FM because of the vast difference between the dynamical sectors of the σ model and the Landau-Lifshitz equation.

A. Numerical solutions
If we include the standard Gilbert damping, as it appears in the discrete Eq. (5), the continuum model (9) is extended as follows, where the damping constant in the discrete and in the continuous models are related by α = ( /2)α. We can derive a relaxation algorithm by assuming that the damping term dominates, α → ∞. This is equivalent to neglecting the second time derivative in Eq. (12) and setting α = 1 (or rescaling time) thus obtaining the relaxation algorithmṅ For any initial configuration, the above algorithm will lead to a local minimum of the energy in the limit t → ∞. Eq. (13) is identical to the relaxation algorithm used for magnetization configurations satisfying the Landau-Lifshitz equation. We have applied (13) in order to find static AFM skyrmion solutions. The result is, obviously, identical to the chiral skyrmion configurations or profiles found for FM [18,19], except that, here, the skyrmion configuration refers to the field n, and, according to Eq. (11), m = 0. We are interested in skyrmions traveling as solitary waves, that is, solutions of Eq. (9) of the form where υ is the velocity of propagation and we have chosen x as the direction of propagation. It is instructive to note that such solutions would be obtained in a straightforward way if the DM interaction were not present in Eq. (9). In that case, the model would be Lorentz invariant, i.e., for any static solution n 0 (x, y), a traveling solution would be obtained by applying the Lorentz transformation and the velocity of propagation can be chosen in the interval 0 ≤ υ < 1. When the DM interaction is present, one cannot obtain traveling solutions by simply invoking the Lorentz transformation, nevertheless, we will find numerically that traveling skyrmion solutions do exist for λ = 0. We insert the traveling wave form (14) in Eq. (9) and obtain Solutions of the latter equation can be found by using the relaxation algorithm (13) where instead of f we have to use f − υ 2 ∂ 2 1 n. We apply the numerical relaxation algorithm using as an initial condition a Belavin-Polyakov skyrmion. We typically use a lattice spacing ∆x = 0.1 and the parameter values We note that the choice κ = 1, used throughout this paper, is not restrictive, as changing κ amounts to simply rescaling the space variable. The relaxation algorithm converges to a skyrmion that is a solution of Eq. (16), for a range of velocities 0 ≤ υ < υ c . For the parameter values (17), we find a maximum velocity υ c 0.715. In Fig. 2, we give the configurations for the field n for a static skyrmion and for skyrmions traveling with various velocities. The traveling skyrmion configurations for n and the associated fields m, k, l, obtained by Eqs. (18) and (11), can be used to find the spins at each tetramer via Eq. (A8). As a check of consistency we have tested the dynamics of the spin configurations S i,j obtained from the configuration n of traveling skyrmions such as those in Fig. 2. We propagate in time under Eq. (5) a skyrmion configuration in the spin lattice and we verify that this propagates rigidly with a velocity 2 √ 2asJυ, where υ is the skyrmion velocity in the σ model and a the lattice spacing in the spin lattice. The factor 2 √ 2asJ is due to the definition of the scaled time (A2) in the σ model. This provides a verification of the consistency of the original equations (5) for the spins with the continuum approximation (9).
For small velocities, e.g., υ = 0.2, the configuration of n in Fig. 2 does not appear to differ a lot from that of the static skyrmion, given in the first entry of the figure. As the velocity increases, the skyrmion is elongated along the axis perpendicular to the direction of propagation. There is also a small elongation of the configuration along the axis of propagation. This is very different than the configuration of traveling solitons under the Lorentz transformation (15), where the soliton is actually just contracted along the direction of propagation. When the velocity approaches the critical velocity υ c the skyrmion core expands in space, apparently to become infinitely elongated in the limit υ → υ c .
The maximum velocity of propagation is υ c < 1, that is, it is lower for the model with DM interaction compared to the value υ c = 1 attained for the case of traveling solitons in the Lorentz invariant model (for λ = 0). We expect that υ c → 1 as λ → 0. Indeed, we find numerically that the maximum velocity υ c increases for lower values of the DM parameter. For example, for λ = 0.4, we obtain numerically that υ c 0.79.
A significant feature of traveling AFM skyrmions is seen in their magnetization vector given in Eq. (11). For traveling solutions of the form (14) we have  The detailed description of the configuration for a traveling skyrmion can serve as a guide for setting-up schemes to obtain these in experiments. Engineering the vector m in order to obtain configurations such as those in Fig. 3 could lead to such methods.

B. Particle-like character of skyrmions
One of the more interesting features of solitons is that they behave as particles. In order to explore the details of its particle-like character, we plot in Fig. 5, the energy (10) of the skyrmion as a function of velocity. For small velocities we may assume that the traveling skyrmion configuration is n(x, y, τ ; υ) ≈ n 0 (x − υτ, y). Then, the energy is where E 0 is the energy of the static skyrmion. We denote by E 0 ex the exchange energy of the static skyrmion and we have as n 0 is an axially symmetric configuration. We set and this will play the role of the rest mass of the skyrmion. We write the skyrmion energy for small velocities This has the form of the energy of a Newtonian particle with a mass M 0 . For the parameter values in Eq. (17) we find numerically M 0 = 14.29, E 0 = 9.73. A measure of the size of the skyrmion is given by the total number of reversals for the third component of the vector n, defined as In the simple case of a circular region of radius R where n 3 = −1, we would have R = N /(2π) and this formula can serve as a definition of the skyrmion radius.
The mass of the skyrmion as well as N depend on the skyrmion profile. Changing the DM parameter λ changes the skyrmion profile. In Fig. 4, we show the skyrmion rest mass M 0 and the Néel vector reversals N as functions of the DM parameter. For λ → 0 the skyrmion radius goes to zero and the skyrmion profile approaches that of the BP skyrmion [18]. Thus, in this limit, M 0 → 4π and N → 0. For λ → 2/π, the skyrmion radius goes to infinity and the skyrmion profile is described via a domain wall similar to the one-dimensional domain wall [19,20]. In this limit, M 0 is proportional to the skyrmion radius and N is proportional to the skyrmion area (radius squared) and they both diverge to infinity. The components of the linear momentum (P 1 , P 2 ) for the skyrmion are given by [8,21] P 1 = − ṅ · ∂ 1 n dxdy, P 2 = − ṅ · ∂ 2 n dxdy. (22) Only the first component P = P 1 is nonzero for skyrmions traveling in the x direction. For a traveling wave as in Eq. (14) we have where we have defined the skyrmion mass that depends on the velocity, M = M(υ). For small velocities the assumption n(x, y, τ ; υ) ≈ n 0 (x − υτ, y) leads to  From Eqs. (20) and (25), we obtain the energymomentum relation for small velocities Note that the latter approximate relation is consistent with the group velocity relation Eq. (27) is exact and it is verified by our numerical results in the entire range of linear momenta. For large momenta the numerical data are fitted very well by the energymomentum relation which is again consistent with the group velocity relation (27). In Fig. 6, we show the dispersion relation (energy vs momentum) for the traveling skyrmion. The combination of the two forms (26) and (28), which are also plotted in the figure, give an excellent approximation for almost the entire range of linear momenta. Virial relations for the traveling skyrmions are derived in Appendix B. Eq. (B10) for the energy can be written as The second term on the right hand side is negative and, in the case of a static skyrmion, it is equal to half the DM energy. We thus conclude that the DM interaction modifies the relation between energy and mass compared to standard Lorentz invariant models. Eq. (29) is taking an interesting form if we use Eq. (B9). We obtain E = Mυ 2 + ∂ 2 n · ∂ 2 n dxdy − λ ê 1 · (∂ 2 n × n)dxdy.
(30) Since the skyrmion gets strongly elongated in the y direction for large velocities, the two last terms on the right hand side become negligible compared to the first term in the limit υ → υ c . We find a very good approximation of the numerical data for large velocities using the simplified version of Eq. (30) This is in agreement with Eq. (28) when P ≈ Mυ c . Aiming to give an intuitive picture of the traveling skyrmion features, we plot in Fig. 7 the numerically obtained values of the mass M and the Néel vector reversals N as functions of the skyrmion velocity. For υ = 0 we have N = 8.59 ⇒ R = 1.17. Both quantities increase as the velocity increases and they seem to diverge to infinity as υ → υ c .

IV. CONCLUDING REMARKS
We have given a detailed description of traveling skyrmions in antiferromagnets with the Dzyaloshinskii-Moriya interaction. The study is based on a nonlinear σ model that is derived as the continuum approximation of the original discrete model for a lattice of spins with antiferromagnetic interactions. Traveling skyrmions are found for velocities 0 ≤ υ < υ c where υ c < 1 is a maximum velocity. We find numerically the configurations of traveling skyrmion solutions. We observe that they are elongated in the direction perpendicular to the direction of propagation and they apparently get infinitely elongated in the limit υ → υ c . We find that a net magnetization is developed with orientation perpendicular to the direction of propagation and we suggest that this could offer a measurable quantity in order to observe propagating AFM skyrmions.
We give the dispersion relation for traveling skyrmions and define their mass. We derive virial relations and obtain exact and approximate relations between the mass the energy and the linear momentum of skyrmions. These clarify their particle-like features and substantiate their Newtonian and relativistic character for low and for large momenta respectively.
We remind the reader that, in the case of an FM, the topological skyrmion with the standard skyrmion number Q = 1 shows Hall dynamics while a non-topological skyrmionium with Q = 0 shows Newtonian dynamics [22]. It has thus been clear that the dynamics of solitons in a FM depends crucially on their skyrmion number. In the presently studied AFM, we have seen that the dynamics of a skyrmion (with Q = 1) is Newtonian, at least for small velocities. It is thus dramatically different that the dynamics of a FM skyrmion. On the other hand, it appears, counterintuitively, to be similar to the dynamics of a FM skyrmionium. The similarities can be seen in the energy and momentum behavior shown in Fig. 5 as well as in the configuration of the traveling AFM skyrmion shown in Fig. 2 when these are compared with the corresponding figures for the FM skyrmionium [22].
The study of dynamics of solitons or other configurations in AFM within suitable σ models can give rise to a wealth of dynamical phenomena [23][24][25] that have not been there in FM studied within the Landau-Lifshitz equation. This opens wide perspectives for the study of AFM dynamics, provided measurable quantities of AFM magnetic order be identified.
Appendix A: Derivation of the continuum model The equations for the four spins at each tetramer in Fig. 1 with indices α, β are derived from Eq. (5) and they read We assume that the spin fields approach continuous limits when → 0 and we use the relations 2 A, and similar relations for the fields B, C and D. The notation ∂ 1 , ∂ 2 denotes differentiation with respect to x, y respectively. Taking appropriate combinations of the Eqs. (A1) we derive dynamical equations for the fields (7). For the set of equations to be consistent in the various orders of , we assume that n ∼ O(1) and m, k, l ∼ O( ). From the definitions in Eq. (7) we find that m · n = k · n = l · n = 0 in the limit → 0. We further rescale time according to In the equations for k, l the time derivative does not enter in the order O( ) and we obtain and The equation for n gives, to O( ), ṅ = 2 √ 2 m × n and this is solved for m to give where the dot denotes differentiation with respect to the rescaled time τ . Finally, the dynamical equation for m gives, in the order O( 2 ), where we have used Eqs. (A3), (A4). We introduce the rescaled parameters κ, λ defined from and we insert m from Eq. (A5) to obtain the equation for the Néel vector, The process of finding an actual AFM configuration proceeds as follows. One first solves Eq. (A7) and then the fields in Eqs. (A3), (A4), (A5) are calculated. Finally, relations (7) are inverted to give the spins at each tetramer A = s(m + n + k + l), B = s(m − n + k − l) C = s(m + n − k − l), D = s(m − n − k + l). (A8)
(B7) This is identical to the virial relation that can be obtained for static skyrmions through Derrick's scaling argument [27] (see also Refs. [22,28]). Thus, Eq. (B7) is satisfied by all static as well as traveling skyrmion solutions presented in this paper. A second virial relation is obtained if we take the combination σ 12 − σ 21 , that gives υ 2 ∂ 1 n · ∂ 2 n dxdy = λ e µ · (∂ µ n × n)dxdy. (B8) The term on the right hand side coincides with the formula for the bulk DM energy. We find numerically that both terms in the above relation are zero for every velocity υ. This is due to the parity symmetries that are apparent in all entries of Fig. 2 for the skyrmion configurations.
For the skyrmions in this paper (of Néel type) the second term on the right hand side is negative. A useful relation for the energy is obtained if we add the last of Eqs. (B6) to the energy and use Eq. (B7), E = (∂ 1 n) 2 dxdy + λ ê 2 · (∂ 1 n × n)dxdy. (B10) The first term on the right hand side is the mass of the traveling skyrmion defined in Eq. (24).
We have verified that all virial relations presented in this Appendix are satisfied by the numerically calculated traveling skyrmion solutions.

The Lorentz invariant model
We conclude this appendix by giving virial relations for the Lorentz invariant model, e.g., the one obtained if we omit the DM interaction in the model (9) (see also Ref. [21]). The following relations should be compared with those obtained above in this Appendix. Let n 0 (x, y) be a static solution of the Lorentz invariant model and n(x, y, τ ; υ) = n 0 (ξ, y), ξ = x where M 0 is the rest mass defined from Eq. (B11) for n(x, y, τ ) = n 0 (x, y). Eq. (B10), for λ = 0, gives The linear momentum is P = − ṅ · ∂ 1 n dxdy = υ (∂ 1 n) 2 dxdy = Mυ. (B14) Eqs. (B12), (B13) and (B14) give the well-known relativistic expression for the energy