Two-dimensional Mo decorated borophenes: high critical temperature, large magnetic anisotropy, and stacking-dependent magnetism

Two-dimensional magnetic materials with high critical temperature, large magnetic anisotropy energy and intrinsic magnetism hold great promise for advancements in spintronics. However, synergizing these attributes within a single material remains challenging. Through the application of swarm-intelligence-based structure searching along with first-principles calculations, we identify two Mo decorated borophene variants, denoted as MoB4 and MoB6, are such candidates with high thermal and dynamical stabilities. MoB4 and MoB6 are characterized as either ferromagnetic or antiferromagnetic metals. Notably, both MoB4 and MoB6 display sizable magnetic anisotropy energy—924 and 932 μeV per Mo atom, respectively—surpassing that of the widely studied CrI3 monolayer, which measures 685 μeV per Cr atom. Monte Carlo simulation suggests the Curie temperature of MoB4 sheet is 390 K, which is above room temperature. Our examination uncovers that bilayer Mo x B y formations exhibit layer-specific van der Waals interactions, contrasting with bilayer borophenes produced experimentally, which display robust interlayer chemical bonding. We determine that the stacking order profoundly influence both the magnetic anisotropy energy and critical temperatures of the material. Specifically, the magnetic anisotropy energy for both structures doubles in their bilayer configurations, with AB-stacked MoB4 and AC-stacked MoB6 demonstrating critical temperatures of 550 K and 320 K, respectively. The exceptional electronic and magnetic characteristics of the Mo x B y monolayers position them as favorable candidates for future spintronic devices.

In this work, extensive structure searches and density functional theory (DFT) simulations led us to the successful discovery of two stable TM borides-the Mo x B y sheets-distinguished by high critical temperatures and large MAE, making them candidates for spintronic applications.We identified two magnetic Mo x B y monolayers: MoB 4 , a ferromagnet, and MoB 6 , an antiferromagnet.Electronic band structure calculations reveal MoB 4 and MoB 6 are both metallic, and both manifest sizable MAE, approaching the order of 1 meV.Monte Carlo (MC) simulations underscore MoB 4 's Curie temperature at 390 K, outpacing many known 2D magnetic materials.Furthermore, bilayer configurations exhibit stacking-dependent band structures, magnetic states, and tunable critical temperatures, enhancing their versatility in device contexts.Our structure search also yielded three stable nonmagnetic Mo x B y configurations, which, although not the focus of this study, are discussed in detail in the supplementary material (figures S1-S6).

Computational methods
The structural search was performed based on the particle-swarm optimization (PSO) technique as implemented in the CALYPSO code [30].2D Mo x B y compounds with several Mo concentrations were systematically studied and the stable phases in MoB, Mo 2 B, MoB 2 , MoB 4 , MoB 6 , Mo 2 B 8 , and Mo 2 B 12 were systematically explored.The population size was set to 30% and 60% of structures in each generation were evolved into the next generation by PSO.The buckled structures were considered by setting the control parameter of buffering thickness to 0.5 Å.The subsequent structural relaxations were performed on basis of DFT as implemented in the Vienna ab initio Simulation Package [31].The electronic exchange-correlation functional was treated by the generalized gradient approximation proposed by Perdew, Burke and Ernzerhof (PBE) [32].To deal with the strongly correlated 4d orbitals of Mo, the Hubbard correction was used with the value U of 4.0 eV [33].To account for the van der Waals (vdW) interactions in the layered structures, we used the Grimme DFT-D3 approach [34].The energy cutoff of the plane waves was set to 500 eV.The structures were fully relaxed until the maximum force on each atom was less than 0.005 eV Å −1 .The energy convergence criterion in the self-consistent calculations was set to 10 −6 eV.For the MAE calculation, the spin-orbit coupling (SOC) was included.A Gamma-centered Monkhort-Pack k-point mesh with a resolution of 2π × 0.03 Å −1 was used for geometry optimization and self-consistent calculations.A vacuum slab of at least 20 Å in z direction was adopted to avoid artificial interactions between the neighboring layers.The phonon dispersion was computed by using the Phonopy code [35] within the density functional perturbation theory [36].In phonon calculations, a finer k-point grid of 2π × 0.02 Å −1 was employed.Ab initio molecular dynamics (AIMD) simulations were performed at 300 K, 1000 K and 1500 K to evaluate the thermally stability of 2D Mo x B y .
For the calculation of critical temperature, MC simulations [37] were performed based on the classical Heisenberg spin Hamiltonian: Here, S i and S j are the unit vector describing the orientation of the local spin moment at site i and j, respectively.k is the anisotropy constant and S z i represents components of S along z (out-of-plane).Note that the anisotropic exchange coupling (λ), which plays a role in the critical temperature of certain materials, was not considered in our Heisenberg model [38,39].For our predicted structures, the simulation geometry was of the size of 20 nm × 20 nm, which contains 2354 and 1613 Mo atoms for MoB 4 and MoB 6 , respectively.The spins were initialized along the easy axis and thermalized for 2 × 10 4 MC steps followed by 5 × 10 4 steps for time averaging at each temperature.The Hinzke-Nowak combinational algorithm [40] was used in the MC simulations for fast relaxation to thermal equilibrium, and the periodic boundary conditions were also included.Table 1.Lattice parameter Lp (in Å), magnetic ground state, space group, magnetic moment M (in µB) per Mo atom, cohesive energy E coh (in eV atom −1 ), formation energy E f (in eV atom −1 ) and the differences in the band edge energy (in eV) between the two spin components for the MoB4 and MoB6 monolayer using PBE + U method.To assess the bonding characteristics within these materials, we calculated the Electron Localization Function (ELF), which is depicted in figure S7.In both structures, ELF peaks are observed between B-B atoms, indicating robust covalent connections, while the ELF is diminished between B-Mo atoms, suggesting a more ionic bonding character.A Bader charge analysis corroborates this, showing that Mo atoms transfer approximately 0.38 and 0.24 electrons to the adjoining B atoms within MoB 4 and MoB 6 , respectively.

Magnetic ground states
Before a detailed investigation of its properties, we first determined their magnetic ground states by comparing the energies of various possible magnetic states at the PBE + U level (see figure S8 and table S3 in supplementary material).We have utilized a 2 × 2 supercell for MoB 4 and a 4 × 1 supercell for MoB 6 , which are sufficiently large to accommodate three antiferromagnetic (AFM) orderings within these structures.We determined the energies of ferromagnetic (FM) and AFM configurations for each structure by fully relaxing the lattices of different magnetic states, as shown in table S3.In the case of MoB 4 , the FM configuration exhibited the lowest energy, confirming its FM nature (figure 2(a)).Conversely, for MoB 6 , the stripy-AFM (AFM2) configuration manifested as energetically favorable over the FM and other AFM states (table S4), signifying an AFM ground state (figure 2(b)).The respective magnetic moments for MoB 4 and MoB 6 -2.19 µ B and 2.82 µ B -are markedly greater than those of other Mo-based compounds previously reported (0.14 and 0.048 µ B ) [43,44].

Stabilities
To assess the structural stabilities of MoB 4 and MoB 6 monolayers, we first calculated their cohesive energies E coh based on the equation: where x and y indicate the number of atoms, E (Mo x B y ) is the total energies of the MoB 4 or MoB 6 layer, E(Mo) and E(B) are the energy of the isolated Mo and the B atoms, respectively.We find that the E coh for MoB 4 and MoB 6 monolayers are −5.41 and −5.55 eV per atom, respectively.The calculated values, derived from the employment of the same method, are lower than the energies of other monolayer metal borides like MnB (−4.09 eV per atom) [45] and Mn 2 C ( −3.40 eV per atom) [46] as indicated by our calculations.We also calculated the formation energy E f for the two sheets.We define it with respect to borophene and solid molybdenum (see details in the supplementary material) and found the E f for MoB 4 and MoB 6 monolayers are −0.056 and −0.053 eV atom −1 , respectively.The negative values indicate that they can be thermodynamically stable Mo−B phases.We also selected the highly stable 3D MoB 2 crystal (space group: R-3M) as the reference point and calculated the energies above the convex hull for the MoB 4 and MoB 6 monolayers (figure S9).Our findings show that the values are 0.204 eV atom −1 for MoB 4 and 0.186 eV atom −1 for MoB 6 , which generally fall below the threshold used to classify the thermodynamic stability of materials [47].
According to the Born-Huang's elastic stability criteria [48], stable 2D square lattices should satisfy C 11 , C 22 , C 66 > 0 and C 11 + C 22 − C 12 > 0, where C ij are the elastic constants.The calculated values for MoB 4 and MoB 6 summarized in table S5 fully meet the criteria, confirming their mechanical stability.
To ascertain the dynamic stability of the MoB 4 and MoB 6 monolayers, we conducted phonon dispersion calculations.As depicted in figures 2(c) and (d), the absence of imaginary phonon modes within the first Brillouin zone signifies that these monolayers are dynamically stable.The thermal stability of the sheets was further explored through AIMD simulations at temperatures of 500 K and 1000 K.The simulations, depicted in figure S10, revealed no significant structural deterioration after 20 ps, affirming the robustness of these materials against high temperatures up to 1000 K. Overall, the pronounced stability of the MoB 4 and MoB 6 sheets suggests that they are good candidates for experimental realization.

Electronic properties
After establishing the pronounced stability of these monolayers, we proceeded to explore their electronic properties.For the FM MoB 4 sheet, it shows metallic feature (figure 3(a)).The projected density of states (PDOS) indicates that the Mo d-orbitals and B p-orbitals strongly hybridize, significantly contributing to the electronic states near the Fermi level (E F ) (figure 4(a)).The asymmetry of spin states in figure 3(a)  unambiguously confirms the presence of magnetic moments in the MoB 4 monolayer.With a total magnetic moment of 3.83 µ B per unit cell, the magnetic contribution chiefly arises from Mo atoms.Each Mo atom carries an on-site magnetic moment of 2.19 µ B , with boron atoms contributing only marginally (−0.068 µ B ).The SOC displays small impact on the electronic property of the material (figures S11-13).For the MoB 6 monolayer, it also exhibits metallic behavior in the band dispersions of both spin-up and spin-down channels (figure 3(c)).Due to its AFM nature, the bands for both spin channels are degenerate.PDOS analysis (figure 4(c)) shows that the Mo d-orbitals and B p-states are predominantly responsible for the states around the Fermi level.

Magnetic anisotropy energy
Magnetic anisotropy is a fundamental aspect critical for establishing long-range magnetic order in nanoscale materials, as magnetization remains vulnerable to spin-rotational fluctuations.Notably, significant magnetic anisotropy is often indicative of enhanced thermal stability of the magnetic orientation.MAE is quantified as the energy difference between magnetic orientations along various axes: where θ is the polar angle and φ is the azimuthal angle.We have charted the computed MAE for the MoB 4 and MoB 6 monolayers as the spins rotate within the x-y and x-z planes, as presented in figure 5. Since MoB 4 exhibits isotropy within the x-y plane, we have focused on illustrating its MAE in the  [46], MnP (166 µeV) [52] and bulk Co (65 µeV per Co atom) [53], and are two orders of magnitude higher than those of bulk Fe (1.4 µeV per atom) and Ni (2.7 µeV per atom) [54].The substantial MAE endows both the MoB 4 and MoB 6 monolayers with considerable potential for use in magnetoelectronic applications.

Critical temperature
To develop a practical spintronic device, the magnetic transition temperature of materials must be comparable to or above the room temperature.Now, we estimate the critical temperature for the predicted monolayers by performing the MC simulations [37].The high reliability of this method has been validated by testing the Curie temperature of CrI 3 monolayer (figure S14).For the MoB 4 monolayer, the nearest neighbor exchange interaction J 1 and the next nearest neighbor J 2 were considered and estimated by the following equations: Here, E 0 represents the free magnetic coupling energy, E(FM), E(AFM1) and E(AFM2) denote the total energy of the FM, AFM1and AFM2 states, respectively.
Similarly, the exchange interaction for the MoB 6 monolayer can be solved by the following equations:  Thus, the values of J 1 (J 2 ) for MoB 4 and MoB 6 are 44.9 meV (2.3 meV), and −106.2 meV (-12.76 meV), respectively.Subsequently, we derived the temperature-dependent magnetization, as illustrated in figure 5(d).The T c for MoB 4 is determined to be 390 K, surpassing room temperature and exceeding the T c of both experimentally realized 2D materials CrI 3 and Cr 2 Ge 2 Te 6 [1,2], and previously predicted MnB monolayer [45,55].The next-next nearest neighbor interaction J 3 would have a minimal impact on the critical temperature of the monolayers (see figure S15), it is therefore excluded from the subsequent discussion.The estimated Néel temperature (T N ) for MoB 6 registers at 80 K, which is comparatively modest.However, as we will discuss in the following section, stacking it into a bilayer configuration can considerably elevate the T N .

Mo x B y bilayers
The recent successful synthesis of bilayer borophenes [56,57] has triggered the investigations of boron materials beyond monolayer limit [58][59][60][61][62].However, the bilayer borophenes that have been experimentally produced exhibit covalent bonds between their atomic layers.This bonding situation is distinctly different from other common bilayer structures such as bilayer graphene or bilayer MoS 2 , which are characterized by vdW interactions between their layers.The strong interlayer covalent bonding in bilayer borophene restricts the mobility of the atoms, thus hampering property manipulation through stacking approaches and, consequently, limiting our capacity to modify the electronic and magnetic properties of the material.In contrast, we have successfully fashioned vdW Mo x B y bilayers and have demonstrated that their magnetic characteristics can be modulated by varying the stacking orders.
For the MoB 4 bilayer, the simplest stacking configuration is the AA stacking, where the top layer is directly above the bottom layer, displaced along the z direction.The AB stacking is achieved by sliding the top layer by the minimum atom-to-atom distance along the direction m 1 , as shown in figure 6(b).The interlayer distances calculated for each stacking are summarized in table 2 and illustrated in figure S16.The ELF results show a clear intraplanar distribution of the bonding states, with an absence of such between the neighboring layers, confirming that the interlayer interactions are of a vdW type.
As depicted in figure 6, diverse stacking orders markedly influence the magnetic ground states of bilayer arrangements.In the case of bilayer MoB 6 , different stacking orders induce various AFM alignments.For AA, AB, and AC stacking, AFM ordering is evident within each layer, with spins alternating in direction between the top and bottom layers for AA and AB stackings, while spins remain parallel for AC stacking.In contrast, bilayer MoB 4 consistently shows FM ground states within the plane, but with antiparallel spin alignment between neighboring layers.Notably, the preferred magnetic easy axis for both structures remains unaltered, oriented along the out-of-plane z direction (figure 6).
Subsequently, we investigated the stacking dependent MAE for the two materials.As shown in table 2, the MAE values for different stacking configurations are in the range of 625-1023 µeV atom −1 .Importantly, the values for AB and AC stacking MoB 6 reach up to the order of 1 meV, indicating the stacking orders have great influence on MAE.By performing MC simulations, we assessed the T N for various stacking configurations of the bilayers.The corresponding values of exchange interactions J 1 and J 2 are collated in table 2. From figure 7, it is evident that the bilayer forms exhibit elevated T N when compared to the monolayers.Specifically, for AB-stacked MoB 4 and AC-stacked MoB 6 , the T N values surpass room temperature, with temperatures of 550 K and 320 K, respectively.Lastly, in figure 6, band structure calculations demonstrate that both MoB 4 and MoB 6 bilayers maintain metallic state across all stacking orders, coherent with single layers.

Conclusion
We have conducted a comprehensive study on the structural, electronic, and magnetic properties of MoB 4 and MoB 6 monolayers based on ab initio calculations.Phonon spectra and AIMD simulations affirm the dynamical and thermal stability of these materials.The MoB 6 sheet is characterized by an AFM magnetic ground state, while the MoB 4 sheet exhibits a FM ground state.Both sheets exhibit a preferred out-of-plane easy magnetization axis and significant magnetic anisotropy, highlighting their potential in magnetic data storage applications.MC simulations corroborate that the Curie temperature of MoB 4 exceeds room temperature, making it particularly noteworthy.
For bilayer structures, we demonstrate that the MAE, critical temperature, and electronic structure are highly dependent on the stacking configuration.Furthermore, both the MAE and critical temperature of MoB 4 and MoB 6 sheets are markedly improved when in bilayer form.Our investigation of 2D Mo x B y sheets offers an encouraging framework for the exploration of complex magnetism in low-dimensional systems and sets the stage for advancements in practical spintronic applications.

Figure 1 .
Figure 1.Atomic structures of monolayer (a) MoB4 and (b) MoB6 in the top and side views.

Figure 2 .
Figure 2. (a) and (b) The magnetic ground states and (c) and (d) phonon dispersions of the MoB4 and MoB6 monolayers, respectively.

Figure 3 .
Figure 3. Band structures of the (a) MoB4 and (b) MoB6 monolayers based on HSE06 method.

Figure 5 .
Figure 5.The MAE of (a) MoB4, and (b) and (c) MoB6 monolayers by rotating the spin within the x-y and x-z planes, respectively.(d) The normalized magnetic moment of MoB4 and MoB6 monolayers as a function of temperature by MC simulations.

Figure 6 .
Figure 6.Crystal structures for different stacking configurations of bilayer (a) and (b) MoB4, and (c)-(e) MoB6.The yellow area in side view represents electronic localization.The isosurface is 0.6 eÅ −3 .HSE06 band structure of bilayer (f) and (g) MoB4 and (h)-(j) MoB6 with AA, AB and AC stacking orders.

Figure 7 .
Figure 7.The simulated normalized magnetic moment of bilayer (a) MoB4, and (b) MoB6 monolayers as a function of temperature by MC simulations.

Table 2 .
The interlayer distance d (in Å), magnetic ground state (GS), magnetic moments M (in µB per Mo atom), magnetic easy axis, MAE (in µeV per Mo atom), and exchange coupling parameters J1 and J2 (in meV) for the bilayer MoB4 and MoB6 with various stacking configurations.