Calculated magnetic exchange interactions in the van der Waals layered magnet CrSBr

Intrinsic van der Waals materials layered magnets have attracted much attention, especially the air-stable semiconductor CrSBr. Herein, we carry out a comprehensive investigation of both bulk and monolayer CrSBr using the first-principles linear-response method. Through the calculation of the magnetic exchange interactions, it is confirmed that the ground state of bulk CrSBr is A-type antiferromagnetic, while there are five sizable large intralayer exchange interactions with small magnetic frustration, which results in a relatively high magnetic transition temperature of both bulk and monolayer CrSBr. Moreover, the significant electron doping effect and strain effect are demonstrated, with further increased Curie temperature for monolayer CrSBr, as well as an antiferromagnetic to ferromagnetic phase transition for bulk CrSBr. We also calculate the magnon spectra using linear spin-wave theory. These features of CrSBr can be helpful to clarify the microscopic magnetic mechanism and promote the application in spintronics.

In recent years, ternary chromium thiohalide compound CrSBr has been extensively studied [19][20][21][22][23][24][25][26][27][28][29][30][31][32]. CrSBr is a 2D magnetic material with a van der Waals (vdW) layered structure along the c axis [22]. Scanning tunneling spectroscopy and photoluminescence studies indicate that CrSBr is a semiconductor with an electronic gap of 1.5 eV [22]. Bulk CrSBr is reported to have triaxial anisotropy, easy to magnetize b axis, middle to magnetize a axis and hard to magnetize c axis [22]. CrSBr bulk material has a high antiferromagnetic ordered temperature T N = 132 K [22], and many theoretical works predict that monolayer CrSBr has higher ferromagnetic ordered temperature [5,[19][20][21]. Driven by the theoretical prediction, Lee et al measured the T C of monolayer CrSBr at 146 K by using the second harmonic generation technique [23]. Although many theoretical studies have investigated the magnetic interactions of undoped CrSBr and estimated the magnetic transition temperature based on the magnetic interactions [5, 19-21, 29, 33]. However, most of them use a cluster approximation and map the total energy differences to a model Hamiltonian, which prevents detailed study of the long-range properties of the exchange interactions. These models generally consider three nearest neighbor interactions, and some studies even consider only two nearest neighbor interactions. Interestingly, the recent magnon dispersions of CrSBr, measured by inelastic neutron scattering, suggest a Heisenberg exchange model with seven nearest in-plane exchanges [31]. Therefore, in order to accurately study the magnetic interactions of CrSBr, we adopt the first-principles linear-response (FPLR) method [34,35].
In this work, using density functional theory (DFT) calculations, we systematically study the electronic and magnetic properties of both bulk and monolayer CrSBr. Our calculations show that bulk (monolayer) CrSBr is an semiconductor with a band gap of 1.42 eV (1.47 eV), which is in good agreement with the experimental results [22]. Using FPLR method, we calculate the magnetic exchange constants. There are five sizable intralayer magnetic exchange terms with small frustration. Although the calculated interlayer interaction J z1 of bulk CrSBr is very weak, it is indeed antiferromagnetic coupling, which is consistent with the experimental results [22]. Based on the calculated exchange constants, we estimate the magnetic transition temperature of bulk (monolayer) CrSBr at 178 K (211 K). In addition, we study the effect of electron doping and strain on the exchange interactions in CrSBr, and find that both strategies can increase T C of monolayer CrSBr. There is an antiferromagnetic to ferromagnetic phase transition for the doped bulk CrSBr, which is confirmed by both FPLR calculations and direct total energy calculations. We also calculated the magnon spectra using linear spin-wave theory.

Method
The electronic band structure calculations have been carried out by using the full potential linearized augmented plane wave method as implemented in WIEN2K package [36]. For the exchange-correlation potential, the generalized gradient approximation (GGA) is used. The vdW interactions in 2D materials can exhibit large many-body effects [37,38]. To better take into account the interlayer vdW forces, there are many approaches in use, including vdW density functional (vdW-DF) [39], DFT-D method [40], and many-body dispersion method [41,42]. For bulk CrSBr, we adopt the vdW-DF in the form of optB88-vdW [43,44] for structure related calculations. GGA + U calculations are also performed for including the effect of Coulomb repulsion in Cr-3d orbital [45]. Here, we use the values of U = 4 eV and J = 1 eV , which have been widely used in previous theoretical works [5,20,28,29]. Using the second-order variational procedure, we include the spin-orbit coupling interaction [46]. Based on the experimental lattice constants a = 3.50 Å, b = 4.76 Å, and c = 7.96 Å [22], we optimize the internal atomic coordinate for bulk CrSBr. The crystal structure of monolayer CrSBr is fully optimized, while the vacuum space is set to be 15 Å to avoid interactions with other neighboring layers. The phonon spectrum is calculated by using the PHONOPY code [47]. The basic functions were expanded to R mt × K max = 7, where R mt is the smallest of the muffin-tin sphere radii and K max is the largest reciprocal lattice vector used in the plane-wave expansion. The 13 × 9 × 6, 13 × 9 × 3, and 13 × 9 × 1 k-point meshes are used for the primitive cell, 1 × 1 × 2 supercell, and slab calculations, respectively. The self-consistent calculations are considered to be converged when the difference in the total energy of the crystal does not exceed 0.01 mRy at consecutive steps. In order to obtain the exact value of total energy, the convergence criteria for the energy difference is change to 0.0001 mRy.
The exchange constants J's are the basis for helping us understand the magnetic properties. Here, we use FPLR method to calculate the exchange interactions, which is based on a combination of the magnetic force theorem [34] and the linear response method [35]. We assumed a rigid rotation of atomic spin at sites R + τ and R ′ + τ ′ of the lattice (here R are the lattice translations and τ are the atoms in the basis). Then, the exchange constant J can be given as a second variation of the total energy induced by the rotation of atomic spin at sites R + τ and R ′ + τ ′ [35], where σ is Pauli matrix and B is the effective local magnetic field. ϵ is the one-electron energy while ψ is the corresponding wave function. This method directly computes the lattice Fourier transform J(q) of the exchange interaction J(R l ), so it is easy to calculate the exact long-range exchange interactions. This technique has been successfully used to evaluate magnetic interactions in a variety of materials [34,35,[48][49][50][51][52][53][54][55].

Bulk CrSBr
The crystal structure of bulk CrSBr belongs to space group Pmmn. The lattice constants are a = 3.50 Å, b = 4.76 Å, and c = 7.96 Å [22]. There are two Cr atoms in each cell. As shown in figure 1, the Cr atom is surrounded by four S atoms and two Br atoms, forming a distorted octahedra. CrS 4 Br 2 octahedra are connected by SBr edge-sharing along the a axis, SS edge-sharing along the ab direction and S corner-sharing along the b axis to form the 2D lattice. Magnetic measurements on single crystals indicate that the magnetic structure of CrSBr is A-type antiferromagnetic, where Cr atoms couple antiferromagnetically along the c axis (see figure 1(a)) [22]. We first perform the GGA + U calculations based on ferromagnetic (FM) configuration. The calculated magnetic moment on the Cr atom is 2.88 µ B , consistent with the high spin state of S = 3/2. As shown in figure 1, we depict the main magnetic interactions. With FPLR method, we estimate and give the magnetic exchange constants with bond lengths less than 8 Å in table 1. Among them, J 1 ,J 2 , and J 3 are the main magnetic interactions, and they are all ferromagnetic. These three exchange interactions determine the ferromagnetic order in the layer. On the other hand, the interlayer first-nearest-neighbor interaction J z1 , although three orders of magnitude weaker than the intralayer interaction, is antiferromagnetic, which causes ferromagnetism to be replaced by antiferromagnetic as the ground state.
Based on the ground state magnetic structure determined above, we perform GGA + U calculations using the 1 × 1 × 2 supercell, and show the band structure of bulk CrSBr in figure 2(a). Our calculations show that CrSBr is an semiconductor with a band gap of 1.42 eV, which is in good agreement with the experimental results (1.5 ± 0.2 eV) [22]. The calculated magnetic moment on the Cr atom is 2.88 µ B , which is the same as the magnetic moment calculated by FM order. The band structures of monolayer CrSBr are similar to those of bulk CrSBr and will be discussed below. The total energy of the antiferromagnetic (AFM) state is about 34 µeV f.u. −1 lower than that of FM state by the direct total energy calculations, confirming the ground state from the calculated magnetic interactions.
Using the FPLR method, we also calculate the exchange interactions for A-type AFM structure. We find that the values of exchange constants with different magnetic configurations are almost the same (the difference is less than 0.005 meV). The fitting magnetic exchange constants in the experimental works are also presented in table 1 for comparison [31]. Our J 1 = −2.31 meV, J 2 = −3.51 meV, and J 3 = −1.40 meV are consistent with the fitting results (J 1 = −1.90 meV, J 2 = −3.38 meV, J 3 = −1.67 meV) of the neutron-scattering measurements [31]. Moreover, J 4 and J 5 are very weak, but J 6 and J 7 cannot be neglected, which is consistent with the experiment [31]. In particular, J 6 is antiferromagnetic, which results in small frustration. As shown in table 1, our J 6 (0.38 meV) is very close to the fitted J 6 = 0.37 meV, while the fitted J 7 (−0.29 meV) is about two times of our J 7 (−0.136 meV).  Based on the calculated magnetic exchange constants, we calculate the magnetic transition temperature using the mean-field approximation theory [56]. T N is estimated to be 178 K, which is somewhat larger than the experimental value (132 K) [22]. Since the mean field theory often overestimates the magnetic transition temperatures, the exchange constants we calculated are considered to agree with the experimental results.
It is worth noting that although the CrSBr system has a global inversion center, most Cr-Cr bonds do not have inversion symmetry, therefore Dzyaloshinskii-Moriya (DM) interactions exist. Using the FPLR approach, we also calculate the DM interactions. D 1 is estimated to be 0.07 meV parallel to the b axis, and D 3 is estimated to be 0.18 meV parallel to the a axis. D 2 should be zero because its bond has an inversion center. The DM interactions of CrSBr are so weak (less than 0.8 meV) and hard to identify from the measured spin wave spectra [31]. It is also worth mentioning that the calculated value of magnetic anisotropy energy is less than 0.1 meV, so we ignore it here.

Monolayer CrSBr
The monolayer CrSBr is phase-stable and exhibits ferromagnetic order below 146 K [23]. Based on our optimized structure (a = 3.545 Å, and b = 4.733 Å), the phonon dispersions of monolayer CrSBr along high symmetry lines are calculated by using the PHONOPY code. As shown in figure 3, there are no imaginary frequencies in phonon dispersions, suggesting that the structure of monolayer CrSBr is dynamically stable.   Similarly, the magnetic exchange interactions of monolayer CrSBr are also calculated using FPLR method, as displayed in table 2. We find that the exchange constants of monolayer CrSBr only change slightly compared with those of bulk CrSBr. J 1 ,J 2 , and J 3 are all ferromagnetic, and there is no frustration between them, which results in the high T C of monolayer CrSBr. The J 1 (−2.94 meV) in monolayer CrSBr is stronger than that (−2.31 meV) in bulk CrSBr. J 2 = −3.70 meV, which is approximately equal to −3.51 meV in bulk CrSBr. We have J 3 = −1.98 meV, which is about one and a half times of −1.40 meV in bulk CrSBr. Also, J 6 is still antiferromagnetic, and the values of J 6 and J 7 in monolayer CrSBr are slightly smaller than those in bulk CrSBr.
Using the mean-field approximation theory [56], we also estimate the magnetic transition temperature of monolayer CrSBr, which is 211 K. The calculated magnetic transition temperature of monolayer CrSBr is slightly larger than that (178 K) of bulk CrSBr, which is consistent with experimental results [23].

Electron doping
In order to test the doping dependence, we perform a series of electron doping calculations using virtual crystal approximation. The number of doped electrons per unit cell varies from 0.1 to 0.7 . The exchange constants under electron doping are calculated by FPLR method. The main exchange constants of bulk and monolayer CrSBr as a function of the electron doping level are shown in figure 4.
For bulk CrSBr, J 1 and J 3 increase significantly with the increase of doping level, while J 2 increase slightly first and then decrease slightly. When the electron doping level of bulk CrSBr reaches 0.7 e/cell, J 1 becomes larger than J 2 . It is worth mentioning that when the electron doping level is above 0.1 e/cell, the interlayer interaction J z1 of bulk CrSBr changes from antiferromagnetic to ferromagnetic. As a result, the magnetic ground state changes from antiferromagnetic to ferromagnetic. To clarify the magnetic ground state of doped bulk CrSBr, we also compare the total energies of the FM and AFM states. When the electron doping level is 0.1 e/cell, the total energy of the FM state is about 15 µeV f.u. −1 lower than that of AFM state, confirming the phase transition from the calculated magnetic interactions.
Similarly, for monolayer CrSBr, J 1 and J 3 increase significantly with the increase of doping, while J 2 increase slightly first and then decrease slightly. J 1 becomes larger than J 2 , when the electron doping level reaches 0.5 e/cell for monolayer CrSBr. It can be expected that the magnetic transition temperature of monolayer CrSBr will increase with the increase of doped electrons. While increasing the electron doping level from 0 to 0.7 e/cell, T C of monolayer CrSBr increases from 211 K to 237 K.

Strain effect
The introduction of additional charge may cause lattice distortion [57], and theoretical studies of monolayer CrSBr show that strain can change the magnetic properties [29]. Therefore, we also check the effect of strain on the exchange constants. The main exchange constants of bulk and monolayer CrSBr as a function of strain are calculated using FPLR method, as displayed in figure 5.
For the intralayer interactions of bulk CrSBr, J 1 is strongly reduced by a compressive strain along the a axis, while J 2 and J 3 are significantly enhanced by a compressive strain along the b axis. The effects of strain on the magnetic exchange constants of monolayer CrSBr (see figures 5(g) and (h)) are similar to those of bulk CrSBr, and are in agreement with theoretical calculations for the monolayer reported in the literature [29]. As shown in figure 5(f), the interlayer interactions J z1 and J z2 are enhanced with the decrease of interlayer spacing. Remarkably, the AFM interaction J z1 is greatly enhanced by a compressive strain along the a axis, and changes from AFM to FM by a tensile strain along the a (or b) axis or compressive strain along the c axis. While the level of compressive strain along the b axis reaches 5%, T C of monolayer CrSBr increases from 211 K to 251 K. If we further increase the compressive strain level, the T C will decrease.

Conclusions
In conclusion, we present a comprehensive investigation of both bulk and monolayer CrSBr by using DFT calculations. The magnetic exchange constants are calculated using the FPLR method. The strongest terms, J 1 , J 2 , and J 3 are ferromagnetic interactions, without frustration between them, which leads to a high magnetic transition temperature for both bulk and monolayer CrSBr. In addition, J 4 and J 5 are very weak, but J 6 and J 7 cannot be neglected and J 6 is antiferromagnetic, which leads to frustration in this compound. On the other hand, although the interlayer interaction J z1 of bulk CrSBr is very weak antiferromagnetic coupling, it determines the magnetic ground state of A-type antiferromagnetism. Moreover, we have demonstrated the effect of electron doping and strain on the magnetic properties of both bulk and monolayer CrSBr. Both strategies are found to increase T C of monolayer CrSBr, and induce the antiferromagnetic to ferromagnetic phase transition of bulk CrSBr. This work demonstrate accurate calculations of magnetic exchange constants for CrSBr which will be helpful for deeply understanding their electronic and magnetic properties as well as promoting their applications in spintronics.

Data availability statement
The data that support the findings of this study are available upon request from the authors.