Investigation of stabilization and survival of skyrmion vortices in the presence of magnetic field disorder in two-dimensional lattices: a case study for Janus dichalcogenides

Due to the lack of inversion symmetry, very large Dzyaloshinskii–Moriya interaction (DMI) has been reported for a series of Janus monolayers of manganese dichalcogenides within the framework of first-principles calculations (Liang et al 2020 Phys. Rev. B 101 184401). However, from the viewpoint of potential applications, the current ongoing research mainly focuses on the magnetism in pristine two-dimensional (2D) materials exhibiting non-zero DMI, and the effects of disorder in such systems remain an open problem since the influence of randomness may create some drastic effects on the magnetism of low dimensional systems. Here, we present Monte Carlo simulation results regarding the magnetic properties of a 2D manganese based Janus dichalcogenide material MnSTe in the presence of quenched random magnetic fields where the local field variables have been sampled from a Gaussian distribution. For the selected benchmark material, it has been found that the magnetic skyrmion vortexes emerging at (10 K, 3 T) may survive in the presence of weak and moderate quenched randomness, which is important from the viewpoint of technological applications. In both the pristine and random field cases, the stabilization of magnetic skyrmions are achieved by the major contribution of the ferromagnetic exchange energy to the total energy of the system, and the materials exhibiting large DMI/exchange ratios may exhibit resilient magnetic skyrmion vortexes in the presence of weak and moderate amount of randomness.


Introduction
Since the experimental realization of graphene as a twodimensional (2D) semi-metallic material in the form of a Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence.Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.single layer of carbon atoms arranged on a hexagonal lattice [1], 2D materials engineering has gained great interest from both experimental and academic perspectives.Since then, the search for 2D materials with intrinsic magnetic properties has been the focus of ongoing research on the exploration of magnetism in low-dimensional systems.In this regard, it is important to fabricate ferromagnetic (FM) materials with high values of Curie temperatures [2].It is also a vital task to enhance the magnetic phase transition temperature to overcome the thermal fluctuations, as it was previously predicted by Mermin-Wagner (MW) theorem [3,4] which states that the ferromagnetism and antiferromagnetism cannot survive in 2D materials described by isotropic Heisenberg model with SU(2) symmetry.In order to overcome the suppression of magnetic order in monolayer (ML) limit, one can count on the magnetic anisotropy [5], and on the experimental side based on the magneto-optical Kerr effect, MW theorem was only recently beaten by the discovery of feromagnetism in ML materials CrI 3 [6] and Cr 2 Ge 2 Te 6 [7].Since then, intrinsic long-range magnetic order in 2D has been reported for a wide variety of material class, including van der Vaals (vdW) crystals [8][9][10], FM semiconductors [11][12][13][14], antiFM (AFM) semiconductors [15][16][17][18], FM and AFM metals [19][20][21][22][23], and hetero-structures, as well as Janus MLs [24][25][26].
Most of the 2D magnets are centrosymmetric, which means that they exhibit inversion symmetry [27].However, the lack of inversion symmetry leads to some peculiar magnetic phenomena as a result of emerging spin-orbit coupling (SOC).For such systems, Dzyaloshinskii-Moriya type interaction (DMI) as an outcome of the broken inversion symmetry and SOC favors the formation of chiral magnetic structures.Moreover, the interplay between DMI, exchange coupling and magnetic field may lead to skyrmion states.Keeping the fact in mind that 2D materials may possess outstanding magnetic properties in comparison to their bulk counterparts, these materials would provide a new paradigm [28], and few studies have reported the emergence of non-zero DMI and stabilization of skyrmion vortices in the ML limit.In this regard, the formation of skyrmions and strong DMI have been reported for Janus MLs [29][30][31][32][33][34][35][36][37][38], multi-ferroics [39,40], and vdW materials [41][42][43][44][45][46].Despite the fact that both FM and AFM materials can host magnetic skyrmions, each magnetic ordering class has its own deficiencies.For instance, skyrmions in FM materials are prone to Magnus effect [47] whereas AFM skyrmions have not been experimentally discovered yet, with the exception of some thin-film materials such as Fe 2 O 3 − Pt [32].For the sake of technological applications, in addition to the enhanced Curie point, materials exhibiting small-sized skyrmions are highly demanded.Therefore, controlling and manipulating the vortex size remains as a vital and notorious achievement.In this context, Behera et al [48] showed for ML CrI 3 that the size of skyrmions gradually reduces with the increase in the magnetic field.Strain engineering also promises an alternative approach for experimentally tuning the diameter and density of skyrmions [30,33,40,48].Upon application of tensile strain, it has been reported in [30,40] that both FM exchange and magnetic anisotropy energies increase whereas the magnitude of DMI is found to be reduced for Janus ML materials which causes the domain walls (DWs) to evolve from chiral DW to FM states under tensile strain.Apart from these, interacting skyrmions have been identified in [49] and it was shown that two skyrmions with the same chirality repel each other.
On the theoretical side, the anisotropic Heisenberg model with DMI parameter is well suited for the investigation of skyrmion formation in ML materials.The 2D Heisenberg model without DMI has been investigated using the Monte Carlo method, and it was found that the anisotropic model exhibits a continuous phase transition obeying Ising universality with critical exponents identical to the 2D Ising model [50][51][52].It is also a crucial task to understand and clarify the influence of randomness on the stabilization of the skyrmion phase in ML magnets exhibiting the lack of inversion symmetry for the spintronic systems which have potential in practical applications.To the best of our knowledge, the magnetic properties of 2D FM lattices with random disorders have been scarcely studied in the literature [53][54][55][56][57][58][59][60][61].Moreover, the problem in the non-zero DMI case still remains as an open problem, which is the primary objective of the present work.Therefore, we aim to clarify the formation and survival of skyrmion states emerging on the 2D Janus ML lattices in the presence of quenched disorder associated with the applied magnetic field.The outline of the paper is as follows: In section 2 we present the details regarding our model.Section 3 is devoted to numerical results and physical discussions.Finally, section 4 contains our concluding remarks.

Model and formulation
In order to simulate the magnetic properties of Janus MLs, we consider the following atomistic spin Hamiltonian where S i is a classical spin vector with unit magnitude located on the nodes of a triangular lattice with six nearest-neighbor spins.In equation ( 1), the former two terms respectively represent the FM isotropic exchange coupling and the anisotropic symmetric exchange energies, the third term denotes the single-ion anisotropy energy, the next term is DMI energy, and finally the last term stands for the Zeeman energy where µ is the magnetic dipole moment of residing atom and B i is the local magnetic field acting on the lattice site i, in units of Tesla (T).Summations have been carried out over the nearestneighbor sites except for the third and the last terms in which all lattice points have been considered.Basically, bond and site disorders are the two types of disorder that can be introduced in magnetism, and one can consider the randomness due to the applied magnetic field as a type of site disorder.The effects of disorder due to the presence of randomness undoubtedly produce drastic effects on the magnetic properties of materials.In this context, in order to investigate the survival of skyrmion states in the presence of disorder, a random field distribution (Gaussian distribution) for the local field variable is introduced as where B is the center of the distribution and the parameter σ controls the degree of disorder.Note that in the limiting case σ → 0, we recover pristine material properties.All the interaction constants used in the simulations have been extracted from [29,30], and they are summarized in table 1.A close inspection of table 1 shows that in the presence of tensile strain (sample identities 2, 6, 8), DMI is reduced whereas other energy contributions are enhanced.Therefore, skyrmions tend to disappear even at low temperature region in the presence of strain.On the other hand, comparison of |d // /J| ratios for the sample IDs (1, 3, 4, 5, 7) suggests that the sample-3 with d // /J = 0.25 would be the best material to consider for the remainder of this work.Therefore, sample-3 is selected as the benchmark material for the simulations.
In order to investigate the evolution of the magnetic skyrmion phase of Janus MLs in the presence of disorder given by the probability distribution (2), we perform extensive Monte Carlo simulations (MCSs) to simulate the system defined by equation ( 1) using the Marsaglia method to update the direction of the classical spin vector [62,63].To reduce the finite-size effects, we impose periodic boundary conditions in each direction.For the calculation of the magnetic properties, we work on a triangular lattice with lateral dimension L = 160 and we consider 10 5 MCSs per lattice site for which the first 20% of the overall data have been discarded for thermalization.Note that an MCS corresponds to a complete sequential sweep of lattice sites.We take into account 20 independent realizations for sample averaging at each temperature to reduce the statistical errors.Initially, we start from a random initial orientation of spins corresponding to the high temperature region, and we gradually decrease the system temperature down to 10K which can be assumed to be the ground state of the system.We record the following properties of interest at each temperature: where N = L 2 is the total number of spins and α = x, y, z. • The magnetization of the system can be calculated using • Magnetic susceptibility is given by • The discretized skyrmion number χ Q and the total chirality χ L are respectively given by In equation ( 6), A ri is the local area of the surface spanned by three spins on every elementary triangle (r i , r a , r b ) and χ is defined as local chirality (see [64] for details).Generally, Bloch and Néel type skyrmions are respectively observed in bulk and ML materials with broken inversion symmetry [48].Therefore, we consider the formation of Néel type skyrmions in our investigations.

Pristine case revisited: magnetic profiles and skyrmion formation
Let us start our investigation by revisiting the magnetic properties of the pristine case.In order to clarify the magnetic behavior and to estimate the Curie temperatures for each sample tabulated in table 1, we plot the temperature dependent magnetization and magnetic susceptibility curves in figure 1 in the absence of DMI.According to our simulation results, each material exhibits a second-order phase transition between paramagnetic and FM states with Curie points ranging between the respective upper (404 K) and lower bounds (129 K).Note that the obtained T c values are slightly different from those reported in [29,30] due to the enhanced resolution of our simulation data, as we have considered more temperature points in our simulations which consequently improves the results given by those works.Next, in figure 2, we focus our attention on the formation of skyrmion vortices when DMI is turned on.As we noted before, we pay particular attention to sample-3 as a benchmark material in the remaining parts of this work due to the large d // /J ratio.In the absence of magnetic field (figure 2(a)), DMI is dominant against the Zeeman energy and it favors the formation of chiral magnetic structures.For B = 0, the real-space spin configurations of sample-3 exhibit labyrinth domains with worm-like structure where DWs are shown by thin green lines, which separate the domains of up and down magnetizations from each other.Upon application of 3 T magnetic field (figure 2(b)), we observe an ensemble of Néel type skyrmions (cf see the inset).For sufficiently larger magnetic fields such as 5 T (figure 2(c)), the majority of spins align with the magnetic field, and it becomes possible to delete skyrmions from the sample, and in this case we observe a few isolated skyrmions in the system.These outcomes are in good agreement with those reported in [29] for MnSTe (i.e.sample-3) and we conclude that our program successfully recovers the limiting pristine case, i.e. σ → 0 limit according to equation (2).Before proceeding further, we present the variation of the discretized skyrmion number χ Q defined by equation ( 6) as a function of temperature and applied magnetic field.The corresponding phase diagram is depicted in figure 2(d).Nonzero values of |χ Q | gives an approximate value for the number of skyrmions emerging on the lattice.We can see from figure 2(d) that the maximum number of skyrmions in the pristine MnSTe originate around B = 3 T magnetic field and T = 10 K temperature region, and that no skyrmion formation takes place above 150 K. Therefore, we will focus our attention on (10 K, 3 T) parameter set when investigating the influence of random magnetic fields governed by the probability distribution (2).

Skyrmion behavior in the presence of random magnetic fields
In figure 3, the thermal magnetization profile of sample-3 has been plotted as a function of disorder parameter σ with the Gaussian distribution center B = 3 T.In the pristine case, we observe a phase transition between paramagnetic and skyrmion phases.Besides, it is clear that the magnetization is depressed to zero with increasing randomness, and it is expected that the skyrmion vortices will vanish with increasing sigma.However, for σ < 20 no considerable variation takes place.This outcome shows that the magnetic skyrmion vortexes emerging at (10 K, 3 T) may survive in the presence of moderate quenched randomness.
Figure 4 shows a comparison of the pristine and random field cases.Here, we show the histogram distributions of the sampled local fields B i , corresponding real space spin configurations, as well as individual energy contributions calculated at 10 K coming from different interaction terms including exchange, DMI and Zeeman shares for the pristine case, as well as for a number of σ values.In all cases, the major contribution to the total energy of the system comes from the FM exchange whereas the exchange anisotropy and single-   ion anisotropy contributions are generally small and do not show any significant variation with increasing randomness.For the pristine system, DMI energy is dominant against the Zeeman term, and as a result of a competition mechanism between these two, the skyrmion phase is stabilized.However, the competition between DMI and Zeeman energies becomes more prominent with increasing σ.For the pristine case, enhanced Zeeman energy favors a spin-polarized state, whereas for large enough values of disorder (such as σ = 10), we see that ideally circular vortex geometry of the skyrmion phase evolves into distorted closed domains due to the growing frustration effects.Accordingly, due to the major contribution of FM exchange energy to the total energy of the system, the skyrmion phase may be stabilized in the presence of weak and moderate random field effects.These results can be compared with other materials enlisted in table 1.For example, it is worth to mention that DMI to exchange ratio for MnSeTe (sample-1) is given by |d // /J| ≈ 0.16, and the material exhibits isolated skyrmion vortices in the pristine case [29].According to our simulations, in contrast to MnSTe (sample-3), |d // /J| ratio of sample-1 is not sufficiently large to produce resilient magnetic skyrmion vortices for moderate quenched randomness.Therefore, we may claim that sample-3 and similar materials exhibiting large |d // /J| ratio may exhibit robust skyrmion patterns in the presence of weak and moderate amount of randomness.
For extremely large amounts of randomness (i.e.σ ⩾ 20, see figure 5), Zeeman energy gains considerable dominance against the DMI energy.However, due to the increased amount of randomness, some of the lattice sites are under the influence of extremely positive valued local magnetic field whereas for the remaining sites, the local magnetic field is acting in the opposite direction.The topology of skyrmions becomes evidently randomized and the skyrmion vortex geometry evolves into disordered stripe domains accompanied by zero net magnetization.

Conclusions
In conclusion, we investigated the random field effects on the formation of skyrmion vortices in two dimensional manganese based Janus dichalcogenides.For the pristine case with σ → 0, we successfully recover the results of the benchmark material MnSTe (sample-3) which is a promising material for spintronic applications with a transition temperature T c = 129 K, and the emergence of strong DMI.Benefiting from the phase diagrams, it was found that the maximum number of skyrmions emerges at (10 K, 3 T) for the pristine system, and magnetic skyrmion vortices are found to survive in the presence of weak and moderate quenched randomness.For both the pristine and random field cases, the major contribution to the total energy of the system comes from the FM exchange energy.Therefore, the interplay between DMI and the Zeeman energies leads to the stabilization of the skyrmion phase for weak and moderate distributions of quenched random fields.The stabilization of the skyrmion state in the presence of weak/moderate disorder is important for spintronic device applications such as information carriers, as the skyrmion vortex becomes immune to stray fields.Last but not least, we found that in the extremely large disorder regime, distorted skyrmions evolve into striped domains, and the skyrmion state vanishes.
We hope that the results presented in this work would stimulate further experimental and theoretical interests in the exploration of magnetic properties for 2D materials research.

Figure 1 .
Figure 1.Zero field magnetization and magnetic susceptibility curves as functions of the temperature for the simulated materials in the absence of DMI.

Figure 3 .
Figure 3. Variation of magnetization as a function of temperature in the presence of randomness.Each curve is accompanied by a particular value of σ standing for the degree of disorder.

Figure 4 .
Figure 4. Left panel: Distribution of sampled random magnetic field using equation (2).Middle panel: Real space spin configurations corresponding to different amount of randomness.Right panel: Ground state energy contributions coming from several components.Each row corresponds to a particular value of randomness parameter σ.

Table 1 .
Simulation parameters defining the materials.