Anisotropic Dzyaloshinskii-Moriya interaction protected by D2d crystal symmetry in two-dimensional ternary compounds

Magnetic skyrmions, topologically protected chiral spin swirling quasiparticles, have attracted great attention in fundamental physics and applications. Recently, the discovery of two-dimensional (2D) van der Waals (vdW) magnets has aroused great interest due to their appealing physical properties. Moreover, both experimental and theoretical works have revealed that isotropic Dzyaloshinskii Moriya interaction (DMI) can be achieved in 2D magnets or ferromagnet-based heterostructures. However, 2D magnets with anisotropic DMI haven't been reported yet. Here, via using first-principles calculations, we unveil that anisotropic DMI protected by D2d crystal symmetry can exist in 2D ternary compounds MCuX2. Interestingly, by using micromagnetic simulations, we demonstrate that ferromagnetic (FM) antiskyrmions, FM bimerons, antiferromagnetic (AFM) antiskyrmions and AFM bimerons can be realized in MCuX2 family. Our discovery opens up an avenue to creating antiskyrmions and bimerons with anisotropic DMI protected by D2d crystal symmetry in 2D magnets.


INTRODUCTION
Topological non-trivial magnetic structures such as chiral domain walls, 1 merons, 2 bimerons 3,4 and skyrmions 5,6 have attracted great research interests due to their rich physical properties and widespread application prospects in spintronic devices. Among these spin textures, magnetic skyrmions have been extensively studied due to their small size, low energy consumption and low driving current density. 7 ultrathin epitaxial Au/Co/W(110) with C2v crystal symmetry, 22 and antiskyrmions with anisotropic DMI are reported in acentric tetragonal Heusler compounds with D2d crystal symmetry 23,24 and non-centrosymmetric tetragonal structure with S4 crystal symmetry. 25 In parallel with the development of hot study of skyrmions in these traditional bulk and multilayer thin films, 2D magnets, e.g., Fe3GeTe2, 26 CrI3, 27 CrGeTe3, 28 MnSe2 29 and VSe2 30 with longrange magnetic orderings have been extensively reported in the last few years, which have been providing an ideal platform to study fundamental properties of magnetism such as magnetooptical and magnetoelectric effect for ultracompact spintronics in reduced dimensions.
Moreover, recent works have proposed that Néel-type skyrmions with isotropic DMI can be realized in 2D Janus magnets, e.g., MnXY, 31 CrXY 32 and multiferroics structures, e.g., CrN,33 BaTiO3/SrRuO3, 34 In2Se3/MnBi2Se2Te2. 35 However, it is worth noting that anisotropic DMI has not been reported yet in 2D magnets. Different from previous materials with isotropic DMI vector along with x and y directions, MCuX2 with special crystal symmetry has an anisotropic DMI vector. Skyrmion Hall effect will cause FM skyrmions with opposite topological charges to propagate in opposite direction, instead of moving parrel to injected current. Antiskyrmions Hall angle is strongly dependent on the direction of apply current related to the internal spin texture of antiskyrmion. When applying spin-polarized current to drive an antiskyrmion, the propagation direction of the antiskyrmion will follow the current direction without topological skyrmion Hall effect. [36][37][38][39][40][41] Therefore, it is possible to achieve the zero antiskyrmion Hall angle in a critical current direction.
In experiments, in order to discover 2D materials, a lot of efforts have devoted to finding materials with characteristics of weak interlayer bonds, which allow their exfoliation down to single layer by mechanical or liquid-phase approach 42    ( , ) and higher t2 ( , , ) levels due to the influence of crystal field. The t2 ↔ ↔ super-exchange interaction favors to the appearance of FM coupling, while the ↔ direct exchange and ↔ direct exchange prefer to AFM coupling. When the d orbit is more than half-filled with electrons, the AFM coupling mainly benefits from ↔ direct exchange. Thus, in the MCuX2 family with d orbit no less than half-filled, the competitive AFM coupling is stronger. Similar results are also depicted in the zinc-blende binary transition metal compounds. 50 Collecting all 3d TM atoms magnetic moments in MCuX2 monolayer [see Table   1]. We can find that Mn atoms have the highest magnetic moments and they monotonically decrease on both sides of Mn. It is obvious that the overall trend across 3d TM row obeys the Hund's rule. 51 Figure 2 shows the calculated NN DMI of MCuX2 structures based on the chirality-dependent total energy difference approach. 52 It is found that all systems have anisotropic DMI and dx and dy have opposite signs along x and y directions, which is consistent with the DMI analysis at the beginning. Besides, the DMI strength varies from 0 to 15 meV/atom. These DMIs are very large compared to many state-of-the-art FM/HM heterostructures and 2D Janus structures, e.g., Co/Pt (∼3.0 meV) 52  shown in Supplementary Figure 3. The main reason is that the magnetic moment of Ni atoms in the system is very small, which leads to a small contribution to DMI. Of course, we also check the different U value from 2 to 4eV, a very small magnetic moments are obtained.

DMI of monolayer MCuX2
Next, according to the Moriya's rules 10  plane perpendicular to the propagation direction <110>. However, we ignore the NNN DMI in our theoretical calculations, because we find that the DMI between NN and NNN atoms differs by about two orders of magnitude by using the four-state energy-mapping analysis, 55 e.g.
VCuSe2 and MnCuSe2 (NNN DMI is -0.082meV and -0.075meV). Although the NNN DMI is neglected, we still observe the magnetic structure of Bloch-type helicoid from the results of the micromagnetic simulation.

Chiral spin textures of monolayer MCuX2
Furthermore, we performe the atomistic micromagnetic simulation based on first-principles calculated materials parameters are shown in Table 1 by using the VAMPIRE software package. 56 To get the dynamics of magnetization, the Landau-Lifshitz-Gilbert (LLG) equation was used with the Langevin dynamics as follows: where is the normal unit vector of ith magnetic atom, + is the gyromagnetic ratio and ' In addition, we also calculate magnetic parameters J1, J2, K and D of monolayer VCuSe2 when the tensile strain increases from 1% to 5%, as shown in Table 2. We find that NN and NNN FM exchange coupling strengths decrease a lot while DMI changes slightly, resulting in the large ratio of D/J. Therefore, antiskyrmion with smaller diameter is achieved under tensile strain. [7] The phase diagram of VCuSe2 monolayer under different stress and temperature is shown in Figure 5.

DFT calculations
First-principles calculations are carried out based on density functional theory (DFT) implemented in Vienna ab-initio Simulation Package (VASP). 61 We adopt Perdew-Burke-Ernzerhof (PBE) functionals of the generalized gradient approximation (GGA). 62 as the exchange correlation potential, and use projector augmented plane wave (PAW) method 63,64 to deal with the interaction between nuclear electrons and valence electrons. We set the cutoff energy of 520 eV for the plane wave basis set, and 24×24×1 k-point with Γ-centered meshes for the Brillouin zone integration. Partially occupied d orbitals of transition metal atoms are treated by GGA+U 65 with U = 3 eV for the 3d orbitals of M and Cu elements. We set a vacuum layer with thickness of 25Å in the z direction to ensure that there is no interaction between the periodic images. The convergence criteria of the total energy in the ion relaxation process and the Hellmann-Feynman force between atoms were set to be 1×10 -7 eV and 0.001eV/Å, respectively. To describe our magnetic system, we adopt the following Hamiltonian model: where Si (Sj) is the normal spin vector of ith (jth) magnetic atom, the J1 and J2 represent exchange coupling constants between Nearest-Neighbor (NN) and Next-Nearest-Neighbor (NNN) atoms, respectively. K is magnetic anisotropy constant and Dij is the DMI vectors. The methods to calculate J, K, D is described in the experimental section.

Magnetic parameters
Exchange coupling constant: we construct a 2×2×1 supercell to study three different magnetic where the positive/negative value corresponds to FM/AFM coupling.
Magnetic anisotropy energy K: magnetic anisotropy energy is defined as the energy difference between in-plane magnetized [100] axis and out-of-plane magnetized [001] axis: NN Dzyaloshinskii-Moriya interaction (NN-DMI) D: we performed DMI calculations using the chirality-dependent total energy difference method. 52 First of all, a 4× 1 × 1 supercell is constructed to obtain the charge distribution of system's ground state by solving the Kohn -Sham equations in the absence spin orbit coupling (SOC). Then, SOC is included and we set spin spirals to determine the self-consistent total energies in the clockwise and anticlockwise rotation. Finally, the energy difference between clockwise and anticlockwise rotation is calculated to obtain the anisotropic D. The DMI can be obtain by following formula:

Phonon spectrum
In calculations, based on the PHONOPY code, 66