Walker-like domain wall breakdown in layered antiferromagnets driven by staggered spin–orbit fields

Within linear continuum theory, no magnetic texture can propagate faster than the maximum group velocity of the spin waves. Here, by atomistic spin dynamics simulations and supported by analytical theory, we report that a strongly non-linear transient regime due to the appearance of additional magnetic textures results in the breaking of the Lorentz translational invariance. This dynamical regime is akin to domain wall Walker-breakdown in ferromagnets and involves the nucleation of an antiferromagnetic domain wall pair. While one of the nucleated domain walls is accelerated beyond the magnonic limit, the remaining pair remains static. Under large spin–orbit fields, a cascade of multiple generation and recombination of domain walls are obtained. This result may clarify recent experiments on current pulse induced shattering of large domain structures into small fragmented domains and the subsequent slow recreation of large-scale domains. Antiferromagnetic systems are becoming an appealing alternative for spintronic-based devices due to the more rapid magnetisation dynamics when compared to their ferromagnetic counterparts. Here, using spin dynamic simulations, the authors demonstrate that the motion of domain walls can achieve supermagnonic speeds in an antiferromagnetic system by means of generation of additional domain wall pairs.

W ithin a general physical context, a topological defect is characterised by a local/core region where order is highly altered whereas far from that core region the order parameter varies smoothly in space 1 . Topological defects can be found in a variety of fields such as in superfluid helium 2 named vortices, in periodic crystal structures, they are called dislocations [3][4][5] and in magnetism, magnetic domain walls (DWs) 6,7 , vortices 8,9 and skyrmions 10 . In elementary dislocations, the underlying reason why these defects are called topological is because they cannot be made to disappear by a continuous deformation of the order parameter in infinite extended films. The measure of the distortion is associated to the so-called Burgers vector whose magnitude expresses the strength of the dislocation which is quantised to integer multiples of the atomic spacing of the lattice vectors 11 . Continuum-based textures such as magnetic textures like DWs, vortices and skyrmions are often characterised by the so-called winding number (a topological charge) whose value stands in direct relation to the magnitude of the normalised Burgers vector. However, topology is not at odds with instability. A relevant example of the previous is the plastic deformation and its role in dislocation dynamics 12 . Under low applied stresses, there is a quasilinear dependence between the dislocation velocity and the applied stress described by linear elasticity theory 13 . However, dislocation mobility is not restricted to subsonic regime under high strain rates. Surpassing of the sonic barrier was first postulated theoretically [14][15][16] and then confirmed by atomistic simulations in bcc tungsten 17,18 and in iron screw dislocations 19 by means of introducing secondary kinks via mother-daughter nucleation process with opposite Burgers vector 20 .
Within the micromagnetic approximation, in a twodimensional infinite magnetic system, the order parameter that breaks the symmetry is reflected by an average quantity: the magnetisation m ! ðx; yÞ at each site in a given lattice. The connection between topology and magnetic media and its role in dynamics became apparent by the pioneering work of several authors [21][22][23] , when studying the dynamics of magnetic topological defects. The analytical approach consists of assuming the centre of the texture as a dynamical variable where the equation of motion can be derived using the Lagrangian formalism. Concretely, in ferromagnets, DWs are defined by the transition region that separates two homogeneously and opposite magnetised regions. In a similar fashion as for dislocations, a quantised topological invariant is defined for the magnetisation field 24 . By framing the magnetic texture on a closed path, it allows to define a map of the magnetisation configuration. Distinct topological textures are characterised by different winding number, w, which counts the number of times the magnetisation is wrapped onto itself 25 . The winding number density along the line at each point in time t, is w x; t ð Þ ¼ À∇ x ϕ x; t ð Þ (with the total winding number, W, being 1 2π R w x; t ð Þdx, which for a 180°domain wall takes values ± 1 2 ). Here ϕ x; t ð Þ is the in-plane angle of the spin at location x at time t of the spin-configuration along the 1D-line. In the case of Bloch walls separating magnetic domains in ferromagnets topological defects arise in the form of Néel lines 26 . Since the 1990s, the study of DW mobility in ferromagnets has been the flagship in Spintronics because of its potential in data storage 27,28 and logic devices 29 . However, a major difficulty to implement such technological proposals arise as a result of the intrinsic instability the DW structure presents when it overcomes a certain threshold velocity. In this Article we show that this interesting feature of DW is present irrespective of whether the magnetic material is a ferromagnet or antiferromagnet. In ferromagnets, this regime is known as the Walker limit 30 and it is characterised by a DW mobility that is the result of the combination of translational and oscillatory dynamics 31 leading to periodic alteration of the DW chirality. Unfortunately, the Walker breakdown (WB) occurs at low domain wall velocities and imposes a limit on a reliable control of the DW motion. An alternative to prevent the WB consists of using nanotubes or nanocylinders 32,33 . These particular geometries have unveiled the spin analogue of the Cherenkov effect which involves the emission of electromagnetic waves by charge particles moving at speeds faster than the speed of light in the given medium 34 .
In the case of antiferromagnets, several studies have provided significant insights on DW motion in antiferromagnets where high speeds have been predicted thanks to prevention of the WB 35,36 , yet there is an upper theoretical limit which cannot be surpassed. Interestingly, the manner in which the spin-space is conformed in antiferromagnets obeys the relativistic kinematics of special relativity where the analogue role played by the photons in conventional space is played by the magnons in a magnetic system. Consequently, by comparison, the ultimate limiting velocity of any texture is given by the maximum group velocity of the magnons, v g . In this work, we reveal that in analogue fashion as it occurs for dislocations under plastic deformation where the dislocations propagate at supersonic speeds, super-magnonic regime of motion can occur through a sequence of multiple nucleated textures mimicking the socalled WB.

Results and discussion
To elucidate the existence of such type of magnetisation dynamics, we combined one-dimensional atomistic spin dynamics simulations of a layered antiferromagnet, Mn 2 Au 37 , with one-dimensional analytical theory. In order to induce magnetisation dynamics in Mn 2 Au (see Fig. 1a), we make use of the predicted staggered field-like torque in such crystal structures 38 , where the effective magnetic field resulting from a staggered-induced spin-density, H ! so , possesses opposite signs at each sub-lattice and gives rise to a spin-orbit torque (see Fig. 1b). Some representative dynamics of the temporal evolution of DW position is obtained by extracting the spatiotemporal evolution of the winding number density shown in Fig. 1c, d. These spatiotemporal evolutions were obtained by firstly ramping the spin-orbit field (SO-field) from 0 to its maximum in a time interval of 10 ps and then it is kept constant till the end of the simulation Fig. 1e, f. The intrinsic inertial character of the translational motion of the DW is revealed by a non-linear trend that appears within the first 4-5 ps. The ballistic nature of the DW motion is revealed by the narrowing of the DW as it approaches v g according to the relativistic kinematics framework. The DW width at a given velocity is given by: rest. This contraction translates into an increase of the DW exchange energy. It is noteworthy to comment that in conventional special relativity theory, the contraction observed from non-inertial reference frame does not result in any incremental (or detrimental) interaction among the constituents of the object. First, the DW is accelerated to a speed that approaches to v g . After these transient dynamics, two distinct behaviours are obtained depending on the applied SO-field. For low values of SO-fields (see Fig. 1c), a saturation effect of the DW speed is observed, characterised by the slope of the maximum winding number density position (DW position), q, with respect to time, t with v = ∂ t q. A linear variation of the DW speed as a function of the SO-field is consistent with previous reports 36 . In such a case, the energy dissipation induced by the viscous damping inhibits the DW to move faster as reported in previous work 35 . In contrast, at an onset of 65 mT, a different dynamical behaviour is observed (see Fig. 1d) which can be defined as the critical stress at which the DW suffers a threshold deformation which leads to the breeding of a DW-pair with trivial composite winding number, w, preserving the overall winding number (see Fig. 1f) . This strongly nonlinear regime which involves the generation of additional particles has been reported under several conditions in ferromagnets, for instance, domain wall WB mediated by the vortex-antivortex generation 33 , vortex-core reversal in spintorque oscillators in nanodots 39 and nanocontacts 40 . However, its occurrence in antiferromagnets has never been reported nor proposed. In ferromagnets, the physical origin for the appearance of new magnetic textures, concretely in the case of vortex-core reversal has been associated to the emergence of a effective magnetic field whose origin is purely dynamical and can be derived from the kinetic part of the Lagrangian density, L ¼ L kin À V. In our system, as the DWs are located in the ferromagnetic basal atomic planes where the coupling between the adjacent planes is antiferromagnetic one can derive the temporal evolution of the kinetic field using the Lagrangian formalismh where M s is the sub-lattice magnetisation and is given by 4μ B /2V where μ B is the Bohr magneton and V is the unit cell volume. Although this emergent field has already been proposed to explain the reversal of the vortex core polarity in ferromagnets where the dynamics are gyrotropic, we use this formalism also to translational motion of DWs in layered antiferromagnets. In antiferromagnets due to the relativistic nature of the dynamics, the DW suffers a substantial deformation due to its translational velocity. The extension and the kinetic field profile with respect to the moving DW depends upon the dynamical regime the DW is at a given instant (see Fig. 2a-d). During the first 4-5 ps of ultrafast acceleration, the kinetic field is located behind the DW and smoothly transitions from the back to the centre and then to the forefront of the DW while mimicking the DW profile. The maximum component of the kinetic field is found to be along thê z-direction and can reach values up~35 T in absolute value (see Supplementary Note 1). Contrary to ferromagnets where the magnetic texture absorbs energy from the kinetic field, we attribute the nucleation of DW pairs to the torque exerted by the kinetic field onto the local magnetisation promoting a reverse magnetic domain along the −x with an approximate extension just before the DW-pair nucleation of the exchange length, l ex = ðA=K 2k Þ 1=2 , where A is the effective exchange stiffness per unit area and K 2∥ is the uniaxial anisotropy (see Table 1). The evolution of this reverse magnetic domain is governed by the competition between the exchange and anisotropy energies which results into breaking the reverse magnetic domain into a DW-pair with each of its constituents holding an opposite winding number (thus the sum of their windings is zero). We identify the time involved in the nucleation to be around 1 ps. Different types of magnetic distribution could in principle appear obeying the boundary conditions over the generated magnetic domain, however, we observed that the initially moving DW and the nucleated DW located just in front share the same winding number. The origin of this arrangement can be attributed to the minimisation of the total Zeeman energy now comprised of the SO-field and the kinetic field (see Supplementary Note 1). Therefore, while the torque provided by the kinetic field leads to a reverse magnetic domain in front of the moving DW with the appropriate extension to accommodate two DWs, the competition between the magnetic energies and the Zeeman energy gives rise to the generation of a DW-pair with a specific ordering.
The temporal evolution of the system after the DW-pair is present in the system is dramatically affected. Hereafter, we call the initially moving DW, DW 1 and the nucleated DWs, DW 2 and DW 3 with DW 2 always being the one closest to DW 1 and DW 3 located in front of DW 2 . As can be seen in Fig. 3a, immediately after the nucleation, DW 1 and DW 2 perform a few oscillations and then get stuck evolving in time but not in space. As DW 1 and DW 2 have the same winding number, they repel each other due to the ferromagnetic exchange interaction (see Fig. 3a). However, the region in between the pair (DW 1 and DW 2 ) has opposite orientation with respect to H so and therefore, in order to reduce the Zeeman energy provided by the SO-field (note that the kinetic field is zero as the DWs are static), the distance between DW 1 and DW 2 has to be minimised. Consequently, as a result of the competition between the two forces, a local minimum appears that corresponds to a certain distance between the DW-pair. We verified this interpretation by calculating the energy barrier related to the distance among DW 1 and DW 2 (see Fig. 3b). The energy minimum corresponds to a given separation between the DWs' centres, which depends on the magnitude of the SO-field (see Supplementary Note 2). From numerical simulations, we obtain a stable distance of 32 nm, whereas from analytical calculation the distance is 29 nm for H so = 65 mT. Recent estimations of the effect of the current to produce a staggered local spin accumulations with opposite polarities are 2 mT per 10 7 A/cm 2 .
Another observation is that DW 3 propagates forward by undergoing a velocity boost that surpasses the maximum group velocity of the magnons. This regime is transient and lasts for several tens of picoseconds which is sufficient for DW 3 to cover over three microns. This mobility regime can be explained as follows. While the DW-pair is nucleated, DW 1 and DW 2 repel each other thanks to the ferromagnetic exchange interaction between them leading to a potential motion in opposite directions for each DW. Despite this fact, they barely move (see inset Fig. 3a) and therefore, all this stored energy is transferred from DW 1 to DW 3 in a similar fashion as reported by Tveten et al. 41 ,  where spin waves were able to transfer linear momentum to a given domain wall. Hence, the DW 3 mobility at the moment of the nucleation consist of the velocity provided by the DW 1 , in addition to the boost coming from the momentum transfer due to the repulsion between DW 1 and DW 2 . Using a collective coordinates approach combined with the Rayleigh dissipation function and the Euler-Lagrange equation of motion (see Supplementary Note 3), we obtain the DW 3 velocity as a function of its distance with respect to the DW-pair to be where α is the Gilbert damping, Δ DW 1 and Δ DW 3 are the DWs' width for DW 1 and DW 3 at the nucleation extracted from the atomistic simulations, respectively, and x(t) is the distance between DW 1 and DW 3 . From numerical simulations, the extracted distance between DW 1 and DW 3 is 17 nm when the nucleation starts occurring and the extracted speed value for DW 3 is around 177 km/s (see Fig. 3c) which is circa 4 times the maximum group velocity that magnons can theoretically attain in Mn 2 Au (see Supplementary Note 4). The maximum speed according to Eq. (2) is circa 133 km/s which corresponds to a distance of 16.1 nm between the DW 1 and DW 3 . It is appealing to think that the special theory of relativity is being violated as no magnetic texture is by principle allowed to propagate faster than the maximum velocity of the magnons in the spin space. However, as soon as the DW-pair (DW 2 and DW 3 ) is present in the system additional interactions need to be taken into account which result into a Lagrangian that is not Lorentz invariant and therefore, the constraint on the maximum speed must be lifted up. Besides, when DW 3 enters into a super-magnonic regime of motion, there is a sudden emission of spin-waves propagating together but never surpassing DW 3 (see Fig. 2d). Analogous effect has been reported in mechanical systems such as in edge dislocation dynamics under shear stress where spontaneous emission of radiation has also been observed 42,43 . We attribute the origin of the emitted spin-waves to a mixture of the so-called Bremsstrahlung effect 44 (or breaking radiation) and the spin-Cherenkov effect 45 . The Bremsstrahlung effect arises due to the deceleration of a charge particle. Therefore, one could speculate that such staggering deceleration could lead to excitations of the spin medium. As shown in Fig. 3d, DW 3 not only surpass the maximum group velocity of the magnons but it also enters in a mobility regime where it moves faster than the phase velocity of Inset represents the domain walls, DW 1 and DW 2 initial oscillations before reaching a stable distance between them. The separation is due to the competition between the Zeeman energy, E Zee , provided by the spin-orbit field, H ! so , which tries to pull them closer while the exchange energy, E exc , tries to separate them as both DWs have the same winding number. b Stable distance between DW 1 and DW 2 of 29 nm calculated for an applied H ! so = 65 mT. c DW 1 , DW 2 and DW 3 velocity as a function of time. At the moment the nucleation starts, the maximum registered speed of DW 3 travelling at supermagnonic speeds is 177 km/s. When DW 1 and DW 2 are fully formed, DW 3 propagates at 133 km/s. A smooth recovery of DW 3 towards submagnonic speeds is observed and lasts for about 100 ps. d Phase, v ph , and group velocity, v g , of the magnons as a function of the wave-vector, k, for an applied H ! so = 65 mT extracted from the dispersion relation. Notice that once DW 3 enters into the supermagnonic regime of motion it also surpasses the phase velocity of the magnons which is manifested by the spin-Cherenkov radiation.
COMMUNICATIONS PHYSICS | https://doi.org/10.1038/s42005-020-00456-5 ARTICLE COMMUNICATIONS PHYSICS | (2020) 3:190 | https://doi.org/10.1038/s42005-020-00456-5 | www.nature.com/commsphys the magnons. It is difficult in our case to clearly separate the existence of Bremsstrahlung and spin Cherenkov radiation. We note in Figs. 1d and 3a that the number of oscillations of DW 1 and DW 2 equals the number of spin-wave ripples travelling with DW 3 and further, that the decay of these ripples appear to stretch over a rather long time. Therefore, one possibility for co-existence of both types of radiation could be that the spin-Cherenkovoriginated spin waves acts in an anti-damping fashion onto the Bremsstrahlung in conjunction with the observation that the radiative ripples do not propagate in opposite direction to DW 3 . It is pertinent to explore how the system reacts when it has at its disposal additional pumped energy. We start by injecting an electric current with the time profile illustrated in Fig. 4a. The electric current has a rising time of 5 ps up to a peak value of H so = 100 mT. The value of H so is kept constant for the following 50 ps before reducing it to zero, with a falling time of 5 ps. The SO-field is set to zero for 50 ps before starting this pattern again three times in succession. The closed blue circles in Fig. 4a represents the maximum π/2-windings, i.e. the number of DWs in each pulse. In each cycle, we can observe an avalanche of DW-pairs. Once the primal DW nucleates the first DW-pair, the DW that propagates from the DW-pair (DW 1 and DW 2 ) at supermagnonic speeds becomes a new breeder and gives rise to new DW-pair. This phenomenon repeats for a 13 times leading to the appearance of 26 additional DWs in the system preserving at all times the overall topological charge. It is worth noticing that the DWs once nucleated do not reorganise into a more stable configuration while the SO-field is maximum giving rise to DWlattice-like structure whose inter-DW distance depends upon the SO-field pulse pattern. The average DW density per pulse is 0.013 DW/nm and it takes~20 ps to generate the whole DW lattice. Here closed orange circles denotes the number of DWs when the SO-field starts decreasing. Once this occurs, the lattice decompresses which, in combination to the attractive exchange interaction among DWs with opposite winding number, leads to recombination of some DWs. We observed a similar value for closed orange circles with pulse number meaning that the same number of DWs recombine once the SO-field reaches 0 mT. The closed green circles represent the number of DWs extracted once the SO-field is about to rise up again. It can be seen from Fig. 4a that only spatially close situated DWs of opposite winding numbers annihilate increasing closed green circles in each pulse. This implies that the accumulated number of remnant DWs in the system after each pulse increases meaning that not all the generated DWs annihilate. After the fourth pulse, the number of DWs in the system is 8 with a spacing among the DWs that ranges from hundreds of nanometres to few micrometers. Initially the system contained two magnetic domains (separate by a single DW) and after four pulses the number of magnetic domains N m , has increased by N m = N dws + 1, where, N dws , accounts for the number of DWs in the system. More complex DW-crystal lattices can be obtained by varying the pulse duration and the sign of the SO-field (see Supplementary Note 5) showing the rich magnetisation dynamics in these type of systems. We wish at this point to also stress some observations (not shown here), that there appears to be no clear correlation between the strength of H so and the number of DWs nucleated as well as the point in time and space where they are nucleated. This opens up the possibility that this system under the current excitation circumstances could potentially exhibit chaotic behaviour. This will be the topic of future work however. We note at this stage, that recent work on current-induced resistance changes in the antiferromagnet CuMnAs was, in conjunction with imaging attributed to the fragmentation and recovery of the domain structure 46,47 . In particular they observed a gradual increase in the resistance with repeated pulsing and a slow relaxation of the resistance towards lower values than after the last pulse when turning off the excitation completely. However, provided that finite size-effects can be excluded and that such a system observes global conservation of the total winding number, it should be possible to return to the original state by applying a reverse field below the AFM Walker-breakdown field, in order to recover completely the initial state of resistance, in other words, a complete reset of the system. Further, it is noteworthy that the relaxed resistance value after a long waiting time was higher than the original starting resistance. This could mean that residual DW are present. If the interpretation that the observed resistance changes is due to domain fragmentation and recovery/recombination is correct we cannot help but to speculate about possible measurements of our system under study. Although the mechanism of domain fragmentation in our case is due to the AFM WB in a perfect crystal and the results in the reported experimental work likely has heat as the main cause of fragmentation, this means that that such resistance variations in analogous systems where pinning and heating effects can be ruled out could be used as an indirect detection method for the occurrence of the AFM WB process and its associated generation of supermagnonic textures. It must be pointed out that for this particular pulse ramp-time no signature of DW-pair nucleation is obtained well after the SOfield has reached H crit = 65 mT, meaning that the critical feld depends upon the ramping time conditions rather than the absolute value of the energy pumped into the system. Figure 4b shows the transition from a system with a DW-lattice with 13 DW-pairs towards a system with only one DW-pair due to the multiple recombination when the SO-field is zero. It is noted that there is a distribution of distances that separate consecutive DWs with opposite winding number. This gives rise to a manifold of recombination times as shown in Fig. 4c. Analytically, it is possible to obtain the recombination time as a function of the distance by assuming that the decompression introduced by the absence of SO-field results in a transition phase from a DWlattice to a DW-gas for a transient time. Therefore, the recombination time of a given DW-pair is only governed by intrinsic parameters such as, the exchange interaction, the damping and the distance between the DWs ruling out the effect of the presence of the rest of DW-pairs in the system. We can derive the recombination time from a general expression of the exchange interaction between two domain wall pairs with opposite winding number as where Δ is the DW-width of each DW, l represents the distance each DW needs to travel before the recombination occurs at t 0 , l 0 = 1 nm and t r represents the recombination time. We note that the predicted recombination time reproduce quantitatively the recombination time values extracted from the atomistic spin dynamic simulations which validates the hypothesis of 1D DWpairwise gas approximation. For distances in the range of few hundreds of nanometres the recombination time is in the range of few picoseconds as the exchange interaction is a short-range interaction. However, we can see from Fig. 4b that there are DWs remaining in the system whose separation is in the range of few microns. For a given distance of 1.8 μm, the expected recombination time lies in the range of~6 days. This suggests that in an ideal scenario where pinning effects, thermal fluctuation effects are not present, the stability of such a configuration is granted thanks to the absence of long-range interaction. We can observe in Fig. 4b that once the SO-field is set to zero, each of the recombination processes produces an excitation in the continuum spectrum. This highly nonlinear process yields to a deformation amplitude that oscillates in time with a very precise frequency. The emergence of this breather-like excitation is due to the excess of kinetic energy carried by each DW involved in the collision whose value is not large enough to escape the attractive potential provided by its anti-particle. By mapping the breather into a simple damped harmonic oscillator, we obtained a good quantitative agreement on the breather lifetime (see Supplementary Note 6). The breather decay occurs in 20-30 ps and it is governed by the Gilbert damping. Moreover, the characteristic frequency of the breather is determined to be 544 GHz frequency, which lies within the linewidth of the calculated spin-wave band gap at H so = 100 mT, i.e. the frequency at zero wave vector, k = 0. We are thus compelled to assign the breather frequency to that of the band-gap.
In summary, we have studied by atomistic spin dynamics and theory, the dynamical properties of staggered field-driven DWs in an AFM. A Walker-breakdown-like process has been identified whereby the DW nucleates an additional DW pair with a total trivial winding number. One of these nucleated DW forms a π DW together with the original one, whereas the other nucleated DW travels away from the breakdown site at speeds far exceeding the maximum group velocity of spin-waves in the medium, thus achieving supermagnonic speeds. At such high speeds we observe a radiative tail travelling together with the supermagnonic texture. According to the computed dispersion relation, for this system, it is only at supermagnonic speeds where the conditions for spin-Cherenkov radiation is met. Thus one contributing source to the observed radiation, we attribute to spin-Cherenkov radiation. At the same time, we would expect a breaking radiation to be present due to a rapid deceleration of DW 1 and DW 2 . We observe associated oscillations in position of the said DWs as they come to a halt, with the number of oscillations equalling that of the number of radiative ripples travelling with DW 3 , opening the possibility of co-existence of the two types of radiation with the spin-Cherenkov radiation acting in an antidamping manner onto the breaking radiation, which could explain the rather slow decay of the radiative ripples. At higher spin-orbit field magnitudes, both multiple DW nucleations and DW recombinations occur, in order to keep conservation of the global winding number. Provided that the separation between nucleated magnetisation textures is sufficiently large, residual DWs can remain in the system after the pulse-excitation is turned off. Provided that the measured resistance in CuMnAs is of magnetic origin, in particular due to magnetic domain fragmentation 46,47 , our results imply that no Joule heating or temperatures approaching the Néel temperature are required in order to observe such resistance characteristics.

Methods
We perform numerical simulations via atomistic spin dynamics simulations. For this, the full Mn 2 Au crystal structure is taken into consideration. Two formula units as that in Fig. 1a in the main text is stacked on top of each and then replicated COMMUNICATIONS PHYSICS | https://doi.org/10.1038/s42005-020-00456-5 ARTICLE COMMUNICATIONS PHYSICS | (2020) 3:190 | https://doi.org/10.1038/s42005-020-00456-5 | www.nature.com/commsphys along the x-direction 60,000 times. The system has open boundaries along x and z while periodic boundary conditions are imposed along y. The time evolution of a unit vector spin at site i, S i , is simulated by solving the Landau-Lifshitz-Gilbert equation: where γ is the gyromagnetic ratio of a free electron (2.21 × 10 5 m/A s), α is the Gilbert damping set here to 0.001 and H eff i is the effective field resulting from all of the interaction energies. The energies taken into account are the three exchange interactions (two antiferromagnetic and one ferromagnetic), magneto crystalline energy contributions and the spin-orbit field. The total energy, E, considered is The first term on the right-hand side is the exchange energy where J ij is the exchange coefficient along the considered bonds (see Supplementary Note 4). The second and third terms are the uniaxial hard and easy anisotropies of strengths K 2⊥ and K 2∥ , respectively, while the fourth and fifth terms collectively describes tetragonal anisotropy. For the in-plane part of the tetragonal anisotropy, u 1 = 110 ½ and u 2 = 1 10 ½ . Finally, μ 0 and μ s are the magnetic permeability in vacuum and the magnetic moment, respectively. We have used μ s = 4μ B , with μ B being the Bohr magneton. It is to be noted that the tetragonal anisotropy was included for sake of completeness as it is present in this material but its role in the high speed dynamics is negligible due to the weak magnitude of its anisotropy constants. The effective field is then subsequently evaluated at each point in time using Eq. (5) as H eff  Table 1.

Data availability
The data generated for the current study are available from the corresponding author upon reasonable request.