Switching of chiral magnetic skyrmions by picosecond magnetic field pulses via transient topological states

Magnetic chiral skyrmions are vortex like spin structures that appear as stable or meta-stable states in magnetic materials due to the interplay between the symmetric and antisymmetric exchange interactions, applied magnetic field and/or uniaxial anisotropy. Their small size and internal stability make them prospective objects for data storage but for this, the controlled switching between skyrmion states of opposite polarity and topological charge is essential. Here we present a study of magnetic skyrmion switching by an applied magnetic field pulse based on a discrete model of classical spins and atomistic spin dynamics. We found a finite range of coupling parameters corresponding to the coexistence of two degenerate isolated skyrmions characterized by mutually inverted spin structures with opposite polarity and topological charge. We demonstrate how for a wide range of material parameters a short inclined magnetic field pulse can initiate the reliable switching between these states at GHz rates. Detailed analysis of the switching mechanism revealed the complex path of the system accompanied with the excitation of a chiral-achiral meron pair and the formation of an achiral skyrmion.

Magnetic chiral skyrmions are vortex like spin textures with particle like properties. They may appear as stable (hexagonal lattice of skyrmions) or metastable states (isolated skyrmions) in magnetic materials as a result of the interplay between Heisenberg exchange, Dzyaloshinskii-Moriya interaction (DMI) 1,2 , applied magnetic field and/ or uniaxial anisotropy. Their nontrivial topology gives rise to intriguing dynamic properties such as topological Hall effect 3 , skyrmion Hall effect 4 and unconventional electromagnetic behavior 5,6 . The theory of thermodynamic stability of magnetic skyrmions was developed by Bogdanov and coworkers [7][8][9][10][11][12] . The interest in skyrmion properties and potential applications strongly increased by several indirect 13 and direct [14][15][16] observations of skyrmions with different techniques in different materials. Magnetic skyrmions are attractive for use in spintronic devices because of their high mobility for low current densities 17 and internal stability 18 , as the relevant interactions, in general, do not depend on the size and shape of the sample. Recently, Fert and coauthors presented a conceptual idea of a spintronic device based on skyrmion motion driven by spin-polarized currents 19 , similar to the racetrack memory based on the domain wall motion 20 but much more energy efficient. Skyrmions in such a device are assumed to be stabilized on top of a ferromagnetic ground state of fixed magnetic polarization.
Here, we propose an alternative scheme which involves the manipulation of an isolated skyrmion (iSk) as a data bit localized in a finite-size domain, similar to an element of Magnetoresistive Random Access Memory (MRAM) 21 . We demonstrate the stability of such iSk in zero magnetic field and the possibility of switching between two degenerate skyrmion states characterized by opposite polarity and topological charge (Q) by an inclined magnetic pulse. We found that the switching between such states takes place via the excitation of a pair of chiral and achiral merons and the subsequent emergence of a transient intermediate achiral skyrmion. Our findings not only explain the microscopic details of the controlled switching process which can be achieved at GHz frequencies but also indicate the possibility of creating an MRAM type of device based on the manipulation of such skyrmion states without the necessity to apply a stabilizing magnetic field. Recently, there has been a surge

Results and Discussion
Metastable skyrmions at zero magnetic field. In our model, we consider a thin film of a chiral magnet, assuming a simple cubic lattice structure with lattice constant a, see Fig. 1a. Direct energy minimization of the model Hamiltonian, see Methods, allows one to identify the ground state of the system with respect to geometrical parameters, e.g. thickness of the magnetic film, and material parameters as coupling constants of Heisenberg exchange (J), DMI (D) and uniaxial anisotropy constant (K). Figure 1b shows the numerically calculated phase diagram of the ground state for an infinite magnetic film composed of three monolayers (thickness, L = 2a) at zero applied magnetic field, presented in terms of two reduced parameters D/J and K/J. The red line corresponds to the second order phase transition between a spin spiral (SS) and a saturated ferromagnetic (FM) state. The period of the SS state goes to infinity and effectively approaches the FM state when D/J or K/J approaches the phase transition line. In the case of a magnetic film of a finite thickness, the free surfaces where the magnetic atoms have a reduced number of neighbors provide an additional spatial modulation of the magnetization also through the whole film thickness 29,30 . Due to the complexity of such a three-dimensional state only numerical calculations allow an identification of the correct phase transitions. The details of the calculation of the phase diagram and a comparison with the analytical solution for the limiting cases of pure two-dimensional and bulk chiral magnets are given in Supplementary Materials S1. According to the micromagnetic continuum model, a solution for metastable iSks can be found for any D/J and K/J below the phase transition line within the FM ground state 8,9 . However, the characteristic sizes of such skyrmion solutions can be smaller than the lattice constant and, thereby, lose their physical meaning. In order to describe the stability of magnetic skyrmions properly, it is therefore crucial to use an adequate discrete model, which allows one to identify the collapse (blue) line in Fig. 1b. Therefore, iSks are stable only within a finite range of D/J and K/J, given by the shaded area in the phase diagram. Due to the absence of a reference magnetic field, the two iSks with mutually inverted spin structures with opposite sign of out-of-plane magnetization (polarity) and topological charge Q are degenerate, see inset in Fig. 1b. These metastable skyrmion states are separated by a finite energy barrier Δ E bar defined by a priori not known complex energy landscape, Fig. 1c, which strongly depends on the material parameters and geometry of the system. In Fig. 1d, the dependence of the size of an iSk is shown for fixed values of K/J marked by arrows in the phase diagram. The size of an iSk goes to infinity when the ratio D/J approaches the transition line, D t and becomes very small close to  29,30 . To investigate the generic features of switching without loss of generality we adjusted the material parameters such that the atomistic spin-dynamics simulations can be performed on a reasonably large domain of 100 × 100 × 3 spins.
Skyrmions switching in a square domain. We have used a time dependent Gaussian magnetic field pulse B p applied in a direction defined by the polar angle θ and the azimuthal angle ϕ, as shown in Fig. 1a, see also Methods. Figure 2 shows the time-dependent snapshots of the switching process for different pulse widths t w , damping parameter α and different boundary conditions: open (OBC) or periodic (PBC), see left panel in Fig. 2. Here and below, the origin of time is given relative to the time of maximum field pulse, t* = t − t p . The initial (left) image shows an iSk with the core magnetization pointing up (red color) while the host is a ferromagnetic state with magnetization pointing down (blue color). The final state shown at the right most image is an almost relaxed state at 150 ps after the pulse maximum. The initial and final states represent mutually inverted spin structures, see inset Fig. 1b.
Following the case of PBC, see the sequence in Fig. 2a, we explain the main mechanism of the skyrmion switching which is based on two simultaneous processes: i) an expansion of the skyrmion core and ii) a homogeneous rotation of magnetization of the surrounding ferromagnetic state. Due to the inclination of the applied field, the core of the skyrmion expands asymmetrically along the direction perpendicular to the projection of the magnetic field on the film surface, B xy p . Opposite to that direction one observes the formation of a region with magnetization opposite to the expanded core, see the blue area at t* = − 5 and 0 ps. Such an excited state has similarity to the vortex-antivortex pair, which can be observed in magnetic vortex core reversal dynamics 31 . Subsequently, the magnetization of the surroundings keeps on turning towards the field direction together with an expansion of the core. When the surrounding magnetization is completely flipped, only the small blue region remains with opposite polarity. Below, we show that this crescent shape object in the snapshots at t* = 10 and 15 ps corresponds to a non-axisymmetric achiral skyrmion which during the relaxation converges to an axisymmetric chiral skyrmion, see t* = 150 ps, with Q and polarity opposite to the initial state. In the case of OBC, Fig. 2b, the switching mechanism remains qualitatively the same as for the PBC, but the change in the polarity of the surrounding magnetization occurs in an inhomogeneous manner. First, one observes the appearance of an up-polarized magnetic region at the free boundary, see the red area at the bottom edge in Fig. 2b, t* = 0 ps. It expands and quickly propagates through the whole domain, which results in the switching of the surrounding magnetization, see t* = 10 and 15 ps. Finally, the spin structure relaxes to an iSk and an additional domain may attach to the edge, see the final image in Fig. 2b. Nevertheless, within a short time of about 100 ps, this domain disappears and only a single chiral skyrmion remains.
In case of realistic damping, α ≤ 0.3, a strong effect of spin-wave interference is observed, see Fig. 2c,d. The spin waves injected and reflected by the free edges of the structure interfere with each other and interact with the expanded core of the skyrmion. This may lead to the skyrmion collapse and/or nucleation of new skyrmions. The final state may appear as a mixed state, composed of multiple skyrmions and domain walls, see e.g. the final image in Fig. 2c. However, a successful one-to-one switching can be controlled, for example by adjusting the parameters of the magnetic field pulse. Figure 2d shows how, by adjusting the pulse width, the multiple skyrmion formation is suppressed and a one-to-one switching is observed.
Reverse switching of skyrmions can be controlled by flipping the polar angle of the applied pulse with respect to the layer plane: θ′ = π − θ. With that, one can achieve skyrmion switching rates of the order of a few GHz. For instance, in our simulations the time interval between two field pulses required for sequential switching has typically a value of 200 ps, which corresponds to a switching frequency of 5 GHz. In the Supplementary movies, we illustrate the back and forth skyrmion switching with an interval of 200 ps driven by a sequence of pulses with alternating polar angles θ = 45° and θ′ = 135°.
Such one-to-one skyrmion switching is by no means trivial: the two iSk states are metastable and the system switches between them across a complex energy landscape with various local minima and two global minima corresponding to the up and down ferromagnetic states. In order to understand better the general switching mechanism, we have studied the time dependence of the energy balance and the topological charge during the switching.  Fig. 3b shows two steps at about t* = 34 ps and 38 ps, where Q rapidly changes from + 1 to 0 and then from 0 to − 1. The snapshots in Fig. 3e-g illustrate the changes in the spin texture before and after these steps. To explain and clarify this step-like behavior, we schematically show the intermediate topological states through which the system passes during the switching in Fig. 4. Figure 4a corresponds to the initial state and Fig. 4b represents a pair of chiral and achiral merons (or half skyrmions) carrying an individual topological charge Q = + 1/2, while the total Q of the pair remains Q = + 1. The schematic picture of this meron pair is structurally equivalent to the state in Fig. 3d (see also Fig. 2a at t* = 0 ps). Because of the interaction with the inclined magnetic field, in the dynamical process the core of the chiral meron is much larger than the achiral one and both have distorted shapes. The state in Fig. 4c represents an ideal achiral skyrmion. The snapshot of an achiral skyrmion in our simulation is shown in Fig. 3e (see also Fig. 2a at t* = 15 ps). The emergence of an achiral skyrmion via excitation of a meron pair is the key stage of the skyrmion switching mechanism. It reflects the tendency of the system to conserve Q, which in turn emphasizes the dominance of the Heisenberg exchange interactions. The initial chiral skyrmion, meron pair and achiral skyrmion belong to the same homotopy class, all three states possess the same Q. It is worth to mention that the meron pair, Fig. 4b, is identical to the so-called bimeron state earlier discussed in ref. 25. Indeed, in case of in-plane anisotropy such meron pair (or bimeron) may appear as stable equilibrium state, while the radially symmetric skyrmions, Fig. 4a,h, become unstable and in such systems they may appear only as an excited state.
Due to the DMI, which provides the largest energy gain for those localized states with a single chirality, achiral skyrmion is energetically unfavorable and within a short time it switches to a chiral one. Indeed, for an ideal radially symmetric achiral skyrmion as in Fig. 4c, the energy contribution of the DMI equals precisely zero, while for mutually inverted chiral skyrmions as in Fig. 4a,h, the energy gain from DMI is the same because of the conserved chirality of the spin structures.
An achiral skyrmion with Q = 1 and a chiral skyrmion with Q = − 1 belong to different homotopy classes and the transition between them may occur only via formation of a singular point (SP), where magnetization locally vanishes and Q becomes undefined. A local violation in the chirality of the spin structure, see dashed square in Fig. 4d, precedes the formation of a SP shown in Fig. 4e. The position of the SP can be associated with the center of a finite size volume where the average magnetization totally vanishes, see solid square in Fig. 4e. Then, SP pushed out from the skyrmion into the surrounding ferromagnetic phase where where it smoothed out and disappear. This, in turn, results in the formation of a so-called half-switched skyrmion state with Q = 0, see Fig. 4f and the corresponding simulation snapshot in Fig. 3f. The system remains in such intermediate state only for a short time (a few picoseconds, see Q = 0 in Fig. 3b). The change in chirality on the other side of the half-switched skyrmion is also accompanied with the appearance of a SP, see Fig. 4g, and finally results in a transition to the chiral skyrmion, see Fig. 4h and snapshot in Fig. 3g. Such a transition results in the second jump of Q, from Q = 0 to − 1. In Fig. 3g, we have shown the chiral skyrmion just after the switching, which converges to an ideal axisymmetric skyrmion after a full relaxation.
The sequence of the states sketched above explains all the details of the evolution of the energy and magnetization in Fig. 3a-c. For instance, the excitation of a meron pair with a large number of spins pointing in the plane rather than in the easy out-of-plane direction results in a pronounced increase of the anisotropy energy E ani . Consequently, the chiral meron core expands which extends the inhomogeneous area and leads to a substantial lowering of the DMI energy E DMI as well as to an increase in the Heisenberg exchange energy E ex . Subsequently, the formation of an achiral skyrmion and its localization leads to a lowering in E ani and E ex but causes an increase in the E DMI . The following transition between achiral and chiral skyrmions involves only a small number of spins around the SP and does not affect the total energy significantly. Note, using micromagnetic simulations we checked that the dipole-dipole interaction does not play a role in the switching process, while for adequate description of dynamical states with SP the atomistic model, used in present work, is more reasonable.  Fig. 5a-i, the snapshots display the sequence of states taken at different moments in time (see corresponding vertical dashed lines in figures j-l) through which the system passes during the switching. Note that the topological charge of skyrmions in OBC is an ill-defined quantity. Thereby, Q for the initial state is found to be + 0.93 and − 0.93 after the switching, see the time dependence of Q presented in Fig. 5j. Similar to the switching mechanism discussed in the previous section, the simultaneous effects of skyrmion expansion and change in the polarity of the surrounding ferromagnetic state result in the formation of an vortex-antivortex pair with distorted shapes, see snapshot c. Spin waves injected from the free boundary start to propagate and interact with the excited skyrmion core. Comparing the switching processes in the square and disk domains, one may conclude that the effect of spin-wave interference in the disk shape domain is less pronounced, compare the snapshots in Figs 2 and 5. After some simulation time at about t* = 5 ps, a localized achiral skyrmion appears, while the magnetization of the surrounding state is almost flipped, see snapshot f. A sharp fall in Q, which appears between the time step marked with f and g represents the transition from an achiral skyrmion to a half-switched skyrmion. Contrary to the PBC, here for the case of OBC, the half-switched skyrmion may appear with nonzero topological charge, see the time variation of Q in Fig. 5j at about t* = 9 ps. The nonzero value of Q is due to the presence of domain attached to the edge, Switching diagrams. The precession of the spins around the applied magnetic field leads to additional excitations that affect the temporal behavior of the magnetization and energies. Such excitations become more pronounced in the case of OBC and show longer attenuation times for small damping α. The influence of damping α and the pulse width t w on the induced dynamics is illustrated in the switching diagram presented in Fig. 6. It shows four noticeable regions corresponding to the excitation of (i) the breathing mode, (ii) the mixed state, (iii) the one-to-one skyrmion switching and (iv) the skyrmion collapse. For short pulses,  t 6ps w , a skyrmion experiences only the breathing mode, i.e. the excitation causes a skyrmion core expansion and relaxation back to the initial state. In this regime, the response of the system to the magnetic pulse is not strong enough to excite the meron pair state and switch the polarity of the surrounding magnetization. On the other hand, for pulses that are too long the switched skyrmion appears when the field is still too strong, which causes the skyrmion to collapse, see the white region in Fig. 6. For the intermediate pulse width, we observe either successful one-to-one skyrmion switching or the nucleation of a mixed state consisting of iSks and domain walls, similar to the state shown in Fig. 2c at t* = 150 ps. The mixed state becomes more apparent for small damping, α < 0.1, due to the strong influence of the interference between spin waves injected by the free edges. However, for  α . 0 1, the generation of the spin waves and their interference are significantly suppressed and this results in an enlargement of the successful switching range. The effect of spin wave interference is also suppressed in the case of a circular shape domain, even for low damping. We identified the optimal polar angle θ of the magnetic field pulse for successful switching to be in range between 35° and 55°. Outside this range, switching is not successful, while the variation in the azimuthal angle ϕ does not affect significantly the switching mechanism. In the present simulations we used B 0 = 1 T. However, by adjusting the material parameters, in particular the anisotropy, stable switching can be obtained for lower fields B 0 ~ 100 mT but longer pulses t w ~ 100 ps.
The skyrmion switching turned to be robust within a wide range of material parameters. We have investigated the skyrmion switching in circular domain with varying exchange coupling constant, pulse width and intensity. The range for successful switching is presented in the switching diagrams, Fig. 7. The major energy contribution in the switching energetics comes from the exchange interaction J, which is usually about one order of magnitude higher than all other interactions such as D and K. In Fig. 7a, a range of parameters corresponding to the successful switching is identified varying pulse intensity B 0 and the absolute values of J for fixed relative values of K/J and D/J. Here, the damping parameter α = 0.1, Gaussian width t w = 8 ps and polar angle θ = 45° have been used. The successful switching regime (red area) in terms of B 0 is sufficiently wide, of about 1.2 T, while the lower and upper critical fields increase linearly with exchange coupling.
Thereafter, we calculated the switching diagram in the parameter space of B 0 and t w for fixed J = 5 meV, Fig. 7b. A wide (red) area within the range between 3 and 6 T and t w s between 3 to 14 ps corresponds to successful one-to-one skyrmion switching. Note, the rough boundaries of the switching area in Fig. 7b reflect the complex energy landscape of the system with a large number of metastable states while inside red region we always find stable one-to-one skyrmion switching. In conclusion, we presented a complete phase diagram for a thin magnetic film, using a discrete model, which includes the DMI and anisotropy energies, that shows stable isolated skyrmion solutions at zero applied field. In particular, we found a collapse line for an isolated skyrmion, which defines its finite parameter-range of existence. The degeneracy of the ground state allows the existence of two skyrmion solutions with mutually inverted spin structure, opposite polarity and topological charge. We demonstrated that the switching between them can be driven by a single inclined magnetic field pulse below 100 ps. The general mechanism of the chiral skyrmions switching follows a sequence of transient topological states: a chiral-achiral meron pair, an achiral skyrmion and a half-switched skyrmion. The newly proposed skyrmion switching mechanism can be achieved in a wide range of material and pulse parameters and allows a repetitive skyrmion toggling, on the GHz scale, which makes these findings of interest for potential applications in MRAM-like devices.

Methods
Model Hamiltonian. The total Hamiltonian of a thin film of chiral magnet 30 is given by where n i = m i /μ s is a unit vector of the magnetic moment at lattice site i; J is the exchange coupling constant; D ij is the Dzyaloshinskii-Moriya vector defined as D ij = Dr ij , with D a scalar constant and r ij a unit vector pointing from site i to site j, see Fig. 1a; K is the out-of-plane uniaxial anisotropy constant, and the last term, b = μ s B, describes the coupling of the magnetic moments to an external applied field B. Here we assume μ s = 2μ B . We restrict ourselves to nearest-neighbor interaction in order to work with a reduced number of variables and conserve generality of the results.

Atomistic Spin Dynamics.
To describe the skyrmion switching dynamics we use atomistic spin dynamic simulations based on the solution of the Landau-Lifschitz-Gilbert equation 32 : where B i is an effective magnetic field defined by  = −∂ ∂ B n / i i , γ is the gyromagnetic ratio and α is a dimensionless damping coefficient. Eq. (2) is solved using currently the most efficient algorithm proposed in ref. 33 realized in the juSpinX code 34 . We used a time-dependent magnetic field pulse defined by a Gaussian function, applied into a direction ê B inclined relative to the surface normal by a polar angle θ and an azimuthal angle ϕ, as shown in Fig. 1a. B 0 , t w and t p are the amplitude, Gaussian width and position of the maximum of the pulse, respectively, see inset in Fig. 1a. The pulse is assumed to penetrate uniformly through the whole domain. The time  The vertices i, j, k of each elementary triangle are numbered in a counter-clockwise sense relative to the surface normal vector n pointing in positive direction of the z-axis. The latter means that the numbering should satisfies the following condition n · (r ij × r ik ) > 0, where r ij is a connection vector directed from lattice site i to j.
The parameter a l = A l /4π can be thought as local topological charge, which takes values in the range of − 0.5 < a l < + 0.5. Note, according to Berg and Lüscher 35 , there is a set of exceptional spin configurations for which Q is not defined but still measurable as A l in Eq. (5) is defined for all possible spin configurations. Figure 8b illustrates our approach for defining the size of an iSk for a discrete model. The tangent line for the skyrmion profile fitted by a linear function at the point of inflection, i.e. with largest inclination angle β, see the straight doted line running through two nearest points. Diameter of the skyrmion is defined as twice the distance between the center of a skyrmion and the intersection point of the tangent with the x-axis, i.e. R iSk .