⇦Chinese Physics Letters, 2023, Vol. 40, No. 7, Article code 077301 Superexchange Interactions and Magnetic Anisotropy in MnPSe$_3$ Monolayer Guangyu Wang (王光宇)1,2†, Ke Yang (杨柯)3,1†, Yaozhenghang Ma (马曜峥行)1,2, Lu Liu (刘禄)1,2, Di Lu (芦地)1,2, Yuxuan Zhou (周宇轩)1,2, and Hua Wu (吴骅)1,2,4* Affiliations 1Laboratory for Computational Physical Sciences (MOE), State Key Laboratory of Surface Physics, and Department of Physics, Fudan University, Shanghai 200433, China 2Shanghai Qi Zhi Institute, Shanghai 200232, China 3College of Science, University of Shanghai for Science and Technology, Shanghai 200093, China 4Collaborative Innovation Center of Advanced Microstructures, Nanjing 210093, China Received 6 April 2023; accepted manuscript online 30 May 2023; published online 26 June 2023 †These authors contributed equally to this work.
*Corresponding author. Email: wuh@fudan.edu.cn Citation Text: Wang G Y, Yang K, Ma Y Z X et al. 2023 Chin. Phys. Lett. 40 077301    Abstract Two-dimensional van der Waals magnetic materials are of great current interest for their promising applications in spintronics. Using density functional theory calculations in combination with the maximally localized Wannier functions method and the magnetic anisotropy analyses, we study the electronic and magnetic properties of MnPSe$_3$ monolayer. Our results show that it is a charge transfer antiferromagnetic (AF) insulator. For this Mn$^{2+}$ $3d^5$ system, although it seems straightforward to explain the AF ground state using the direct exchange, we find that the nearly 90$^\circ$ Mn–Se–Mn charge transfer type superexchange plays a dominant role in stabilizing the AF ground state. Moreover, our results indicate that, although the shape anisotropy favors an out-of-plane spin orientation, the spin-orbit coupling (SOC) leads to the experimentally observed in-plane spin orientation. We prove that the actual dominant contribution to the magnetic anisotropy comes from the second-order perturbation of the SOC, by analyzing its distribution over the reciprocal space. Using the AF exchange and anisotropy parameters obtained from our calculations, our Monte Carlo simulations give the Néel temperature $T_{\rm N}=47$ K for MnPSe$_3$ monolayer, which agrees with the experimental 40 K. Furthermore, our calculations show that under a uniaxial tensile (compressive) strain, Néel vector would be parallel (perpendicular) to the strain direction, which well reproduces the recent experiments. We also predict that $T_{\rm N}$ would be increased by a compressive strain.
cpl-40-7-077301-fig1.png
cpl-40-7-077301-fig2.png
cpl-40-7-077301-fig3.png
cpl-40-7-077301-fig4.png
cpl-40-7-077301-fig5.png
cpl-40-7-077301-fig6.png
cpl-40-7-077301-fig7.png
cpl-40-7-077301-fig8.png
cpl-40-7-077301-fig9.png
 DOI:10.1088/0256-307X/40/7/077301 © 2023 Chinese Physics Society Article Text Two-dimensional (2D) van der Waals (vdW) materials have attracted widespread attention for several decades due to their highly tunable physical properties and potential applications in multifunctional electronic devices.[1-3] Recently, with the 2D magnetism observed in atomically thin CrI$_3$ and Cr$_2$Ge$_2$Te$_6$, the 2D vdW magnetic materials have become one of the hot topics in condensed matter research.[4,5] These cleavable materials provide an ideal platform for exploring magnetism in the 2D limit, where new physics and novel properties are expected. Moreover, as the vdW interaction between the magnetic layers is usually weak, one may combine these 2D vdW magnetic materials to construct heterostructures, which could produce new properties that are not found in the individual constituents and thus offer immense potential for spintronic applications.[6-13] One practical strategy to search 2D magnetism is to explore the single layers exfoliated from a bulk vdW material that exhibits robust magnetic order. The transition metal phosphorous trichalcogenides $M$P$X_3$ ($M$ = Mn, Ni, Fe, Cd; $X$ = S, Se) are an important 2D vdW magnetic material platform.[14-18] Among them, MnPSe$_3$ has been extensively studied due to its fascinating and promising properties.[19-27] MnPSe$_3$ crystallizes in the space group $R\bar{3}$ (No. 148), and the honeycomb layers form the rhombohedral ABC stacking along the $c$ axis. The P–P dimer is located vertically across the center of the honeycomb formed by Mn ions, as seen in Fig. 1. MnPSe$_3$ bulk is a vdW layered antiferromagnetic (AF) insulator having the in-plane spins with $XY$ anisotropy and the Néel temperature of $T_{\rm N}=74$ K,[19] and it has the charge-transfer (CT) type optical gap of 2.27 eV.[20] A very recent experimental study showed that MnPSe$_3$ monolayer is dynamically stable and can exist as a 2D freestanding crystal,[25] and that it retains the AF ground state with the decreasing $T_{\rm N}=40$ K. Moreover, it was found that under a uniaxial tensile strain, the Néel vector is parallel to the strain direction.[25] So far, however, there has been no systematic theoretical study concerning the magnetic and anisotropic properties of MnPSe$_3$ monolayer, which stimulates our present work. In this Letter, we provide insights into the exchange interactions, magnetic anisotropy, and strain effect of MnPSe$_3$ monolayer, using density functional calculations in combination with the maximally localized Wannier functions (MLWFs) method and the magnetic anisotropy analyses. Our calculations show that MnPSe$_3$ monolayer is a CT-type insulator with the AF ground state. We demonstrate that the superexchange interaction plays a more important role than the direct exchange in the stabilization of the AF ground state. As for the magnetic anisotropy energy (MAE), which is crucial for 2D magnetism according to the Mermin–Wagner theorem,[28] our results show that the SOC contributes much more via the second-order perturbation effect to the MAE than the shape anisotropy. We also carefully analyze the distribution of the SOC perturbation energy in the reciprocal space. Using the AF exchange and anisotropy parameters obtained from our calculations, our Monte Carlo simulations show that MnPSe$_3$ monolayer has $T_{\rm N}=47$ K, in agreement with the experimental 40 K.[25] Moreover, we find that under a uniaxial tensile (compressive) strain, the Néel vector would be parallel (perpendicular) to the strain direction, and the $T_{\rm N}$ would be increased by a compressive strain.
cpl-40-7-077301-fig1.png
Fig. 1. Top view (a) and side view (b) of MnPSe$_3$ monolayer drawn by the VESTA code,[13] showing the honeycomb Mn ions with up (down) spins, the edge-sharing MnSe$_6$ octahedra, and the P–P dimers.
Computational Details. Density functional theory (DFT) calculations were carried out using the Vienna ab initio simulation package.[29] The projector augmented-wave method with a plane-wave basis set was used, and the cut-off energy was set to 400 eV.[30] The exchange and correlation energy was described by the generalized gradient approximation (GGA) with the Perdew–Burke–Ernzerhof functional.[31] The lattice parameters of MnPSe$_3$ monolayer are optimized to be $a=b=6.393$ Å, which are almost the same as the experimental bulk values of $a=b=6.387$ Å.[19] In the structural optimization, all the atoms are fully relaxed till the total energy converges to $10^{-6}$ eV and forces to 0.01 eV/Å. To avoid an interaction between the layer and its periodic images, a vacuum of 20 Å was set along the $c$-axis. The Monkhorst–Pack $k$-mesh of $12\times12\times1$ was used for integration over the Brillouin zone.[32] To describe the on-site Coulomb interactions of Mn 3$d$ electrons, the GGA+$U$ method[33] was used with a typical value of the Hubbard $U =4.0$ eV.[34-36] As seen below, our calculations with $U = 4.0$ eV well reproduce the experimental insulating gap, and moreover, the calculated exchange and anisotropy parameters account for the experimental $T_{\rm N}$ reasonably well. We also included the SOC effect in our GGA$+U+$SOC calculations to study magnetic anisotropy. The tight-binding parameters were obtained from MLWFs analysis using the Wannier90 package,[37,38] and they were used to check the relative importance of the direct exchange and the superexchange interactions. Moreover, we carried out Monte Carlo simulations to estimate the $T_{\rm N}$ of MnPSe$_3$ monolayer, using the Metropolis method[39] and a $20 \times 20$ spin lattice with a periodic boundary condition. Such a $20 \times 20$ spin lattice turns out to be sufficiently large to ensure a reliable simulation of the magnetic transition temperature. During the simulation steps, each spin was rotated randomly in all directions. At each temperature, we used $3 \times 10^{8}$ Monte Carlo steps to reach an equilibrium and $6 \times 10^{8}$ steps for a statistical average. The magnetic specific heat at a given temperature was calculated according to $C_v=(\langle E^2 \rangle - \langle E \rangle^2) / (k_{\scriptscriptstyle{\rm B}} T^2)$, where $E$ is the total energy of the spin system.
cpl-40-7-077301-fig2.png
Fig. 2. Density of states (DOS) of MnPSe$_3$ monolayer by GGA+$U$. It is an AF insulator with a CT-type gap from Se $4p$ to Mn $3d$. The blue (red) curves in (b) stand for the $t_{\rm 2g}$ ($e_{\rm g}$) orbitals. Fermi level is set at zero energy.
Results and Discussion—Electronic Structure and Magnetic Exchange Interactions. We performed GGA+$U$ calculations to investigate the electronic structure of MnPSe$_3$ monolayer in the AF insulating ground state. The Mn 3$d$ orbitals in the local MnSe$_6$ octahedral crystal field split into the $t_{\rm 2g}$ triplet and $e_{\rm g}$ doublet. Our calculations show that each of Mn ions has a local spin moment of 4.6$\mu_{\scriptscriptstyle{\rm B}}$, indicating the formal Mn$^{2+}$ state with the $t_{\rm 2g}^3e_{\rm g}^2$ ($S=5/2$) configuration. The fully occupied up-spin $t_{\rm 2g}^3e_{\rm g}^2$ states are also clearly seen in Fig. 2(b). A reduction of the spin moment is due to a strong covalency between Mn $3d$ and Se $4p$. Figure 2 shows that the topmost valence bands arise from the Se $4p$ states, and the Mn $3d$ conduction bands lie above the Fermi level by about 2 eV. The strong Mn $3d$–Se $4p$ covalency was also proposed by the recent x-ray photoemission and x-ray absorption spectroscopy,[40] and the CT-type excitation gap of 2.27 eV from Se $4p$ to Mn $3d$ in the previous optical absorption study[20] is well reproduced here by our present calculations. Therefore, our calculations suggest that MnPSe$_3$ monolayer is a CT insulator. Now we study the magnetic properties of MnPSe$_3$ monolayer. Our calculations show that the AF Néel state is the ground state and it has a lower total energy than the ferromagnetic (FM) state by 26.4 meV per formula unit (f.u.). The AF exchange parameter may be estimated to be 1.41 meV (26.4/3$S^2$) by mapping the energy difference to a Heisenberg spin Hamiltonian ($H = \sum_{i,j} \frac{J}{2} {{\boldsymbol S}_{i}} \cdot {{\boldsymbol S}_{j}}$) for the Mn$^{2+}$ ions ($3d^5$, $S = 5/2$, counting each pair by $JS^2$) [see more information in the Supplemental Materials (SMs)].[41] As the Mn$^{2+}$ ions are fully spin polarized, their possible direct exchange can only be AF. Then it may be straightforward to explain the AF ground state of MnPSe$_3$ monolayer simply by the direct exchange. However, the superexchange of neighboring magnetic ions through shared ligands usually plays an important role in determining the magnetic order of transition metal compounds.[42] Therefore, we need to check the relative importance of the direct exchange and superexchange using the MLWFs analyses as given in detail below.
cpl-40-7-077301-fig3.png
Fig. 3. The GGA+$U$ calculated DFT (black solid lines) and Wannier-interpolated (red dashed lines) band structures of MnPSe$_3$ monolayer.
cpl-40-7-077301-fig4.png
Fig. 4. Contour-surface plots of the Mn 3$d$ [(a)–(e)] and Se 4$p$ [(f)–(h)] Wannier functions. For all plots we choose an isosurface level of $\pm 3.0$ (yellow for positive values and cyan for negative values) using the VESTA visualization program.[13]
Direct exchange and superexchange originate from the electronic hopping of Mn–Mn and Mn–Se–Mn, respectively. To obtain the hopping parameters, we constructed the MLWFs using Mn 3$d$, P 3$s$ and 3$p$, and Se 4$s$ and 4$p$ orbitals as initial guesses. Figure 3 shows the comparison between DFT and Wannier-interpolated band structures of MnPSe$_3$ monolayer, indicating that the chosen MLWFs perfectly reproduce the ab initio electronic states. The spatial localization is the most important property of the MLWFs, and it is intimately related to the accuracy of the hopping parameters. To verify the localized features of the MLWFs, we plot the Wannier functions for Mn 3$d$ and Se 4$p$ orbitals in real space. As shown in Fig. 4, the 3$d$ and 4$p$ Wannier functions are centered on Mn and Se ions, respectively. Note that we choose the local octahedral $xyz$ coordinate system, with the $xyz$ axes being directed from Mn to neighboring Se ions. In addition, the shape of these Wannier functions is consistent with the $d$ and $p$ orbitals in isolated atoms. These results suggest that the Wannier functions are well localized in real space. Therefore, with the help of these Wannier functions, the hopping integrals of different ions and orbitals can be obtained, and they are crucial for checking the relative importance of the direct exchange and superexchange interactions. We first investigate the direct exchange in MnPSe$_3$ monolayer. Figure 5(a) shows the schematic diagram of the direct exchange. The direct hopping between two neighboring Mn$^{2+}$ ions is allowed only for spin-antiparallel electrons. We summarize in Fig. 5(b) the hopping parameters from the MLWFs calculations, and find that the $d_{x^2-y^2}$ orbitals of two neighboring Mn ions have the largest hopping integral of 0.11 eV and the $d_{xy}$ orbitals the second largest one of 0.10 eV. All other hopping integrals are about half of them in strength or even less, thus their contribution to the direct exchange is of less concern.
cpl-40-7-077301-fig5.png
Fig. 5. (a) Schematic diagram of the AF direct exchange channel via $d_{x^2-y^2}$ orbitals of Mn1 and Mn2. (b) Hopping integrals of spin-down channels between 3$d$ orbitals of Mn1 and Mn2 calculated by MLWFs basis set. See Fig. S3(a) for the hopping integrals of spin-up channels in the SMs.[41] (c) The CT-type AF superexchange channel via the same $p_x$ orbital. (d) Hopping integrals of spin-down channels between Mn 3$d$ and Se 4$p$ orbitals. See Fig. S3(b) in the SMs[41] for the hopping integrals of spin-up channels. For brevity, we use the notation $d_{x^2}$ for $d_{x^2-y^2}$, and $d_{z^2}$ for $d_{3z^2-r^2}$.
As for the superexchange via Se 4$p$ orbitals in this CT-type insulator, here we focus on the following virtual hopping processes with double holes of Se 4$p$ orbitals in the intermediate excited state after two CT excitations from Se $4p$ to Mn $3d$: $d_i^n p^6 d_j^n \xrightarrow[]{\boldsymbol{1}} d_i^{n+1} p^5 d_j^{n} \xrightarrow[]{\boldsymbol{2}} d_i^{n+1} p^4 d_j^{n+1} \xrightarrow[]{\boldsymbol{3}} d_i^{n+1} p^5 d_j^{n} \xrightarrow[]{\boldsymbol{4}} d_i^n p^6 d_j^n$. The superexchange channels in MnPSe$_3$ monolayer can be divided into two categories: via different $p$ orbitals or via same $p$ orbitals but with different spins. As the Se $4p$ electrons on their fat orbitals have negligibly weak Hund exchange, the intermediate double hole states of two different $p$ orbitals either in the spin = 1 state or in the spin = 0 one practically have no energy difference, and therefore have no net contribution to the magnetic couplings. In strong contrast, the superexchange channels via same $p$ orbitals allow only AF coupling in MnPSe$_3$. As seen in Fig. 5(c): first, one spin-down $p$-electron from the intermediate Se ion hops to Mn1, then the other spin-up $p$-electron on the same orbital hops to Mn2, and finally these electrons hop back. The former $p_x-d_{x^2-y^2}$ hopping is $pd\sigma$ type, and the latter $p_x-d_{xy}$ hopping is $pd\pi$ type. Thus, it is unsurprising that the former is larger than the latter, see Fig. 5(d). Moreover, due to the shorter Mn–Se distance of 2.72 Å (than the Mn–Mn one of 3.70 Å) and the more expansive Se $4p$ orbital (than the Mn $3d$ one, see Fig. 4), it is evident that the hoppings between Mn and Se ions ($t_{\rm pd}$) in the superexchange are much larger than the above direct $t_{\rm dd}$ hoppings, as seen in Figs. 5(b) and 5(d). As the strength of the Mn–Mn direct exchange is proportional to $t_{\rm dd}^2$/($U+4J_{\scriptscriptstyle{\rm H}}$), and that of Mn–Se–Mn superexchange is proportional to $t_{\rm pd}^2 t_{\rm pd}'^{2}/\varDelta ^3$,[42] and considering the Hubbard $U=4$ eV and Hund exchange $J_{\scriptscriptstyle{\rm H}}\approx 0.8$ eV of the Mn $3d$ electrons and the charge transfer energy $\varDelta\approx 4.5$ eV estimated from our MvLWFs calculations, we can infer that the superexchange is much stronger than the direct exchange in MnPSe$_3$ monolayer, for example, the $(x^2-y^2)\longleftrightarrow p_x-xy$ superexchange is about 5.0 meV and the $(x^2-y^2)\longleftrightarrow (x^2-y^2)$ direct exchange is about 1.7 meV. To summarize, our above results show that MnPSe$_3$ monolayer is a CT-type AF insulator, which is in good agreement with the previous experiments.[20,25] Although the direct exchange, being absolutely AF for those fully spin-polarized $S=5/2$ Mn$^{2+}$ ions, seems straightforward to explain the AF ground state, we also care about the importance of the common superexchange. Then using the MLWFs calculations of the Mn–Mn and Se–Mn hopping parameters, we find that the Mn–Se–Mn superexchange via the strong Mn–Se covalency is more important than the direct Mn–Mn exchange in determining the AF ground state. Thus, our present results clarify the CT-type AF insulating behavior of MnPSe$_3$ monolayer and prove that the superexchange rather than the direct exchange plays a major role in stabilizing the AF ground state. The MAE and its Distributions in Reciprocal Space. Now we explore the magnetic anisotropy in MnPSe$_3$ monolayer. In 2D materials, the magnetic anisotropy breaks the spin rotational invariance and is crucial for stabilizing a long-range magnetic order against a thermal fluctuation. The SOC-induced magnetic anisotropy and the magnetic shape anisotropy are two common ones appearing in 2D magnetic materials.[43-45] The SOC is included via the Hamiltonian term $H_{\scriptscriptstyle{\rm SOC}} = \xi \sigma \cdot L$ ($\xi$ is the coupling constant) into our GGA+$U$+SOC calculations. We assume the spin orientation along the $a$, $b$, and $c$ axes, respectively, and then run their self-consistent calculations till the energy difference converges within 2 µeV/f.u. The results show that MnPSe$_3$ monolayer has the easy $ab$ plane and the $c$ axis spin orientation has a higher energy by 241 µeV/f.u. The shape anisotropy energy resulting from the dipole-dipole interaction is \begin{align} E_{\rm dipole-dipole}=\,&\frac{1}{2}\frac{\mu_0}{4\pi}\sum_{i \neq j}^{N} \frac{1}{r^{3}_{ij}}\Big[\boldsymbol{M}_{i} \cdot \boldsymbol{M}_{j}\notag\\ &-\frac{3}{r^{2}_{ij}}(\boldsymbol{M}_{i} \cdot \boldsymbol{r}_{ij})(\boldsymbol{M}_{j} \cdot \boldsymbol{r}_{ij})\Big], \tag {1} \end{align} where $\boldsymbol{M}_{i}$ represents the Mn$^{2+}$ magnetic moments and $\boldsymbol{r}_{ij}$ is a vector connecting the Mn-sites $i$ and $j$. Our results show that the shape anisotropy prefers the easy $c$ axis, and the $ab$ plane spin orientation has a higher energy by 46 µeV/f.u. Taking into account both the contributions of SOC and shape anisotropy, we find that $c$ axis is the hard magnetization axis with the MAE of 195 µeV/f.u. relative to the easy $ab$ plane. Obviously, the SOC is dominant and it determines the easy plane spin orientation of MnPSe$_3$ monolayer as experimentally measured.[25] As the SOC induced MAE is dominant, here we have a close look at it. The energy gain from SOC according to the perturbation theory is \begin{align} E_{\scriptscriptstyle{\rm SOC}} = \langle o | H_{\scriptscriptstyle{\rm SOC}} | o \rangle + \frac{|\langle o | H_{\scriptscriptstyle{\rm SOC}} | u \rangle |^2}{E_o - E_u}, \tag {2} \end{align} where $| o \rangle$ and $| u \rangle$ stand for occupied and unoccupied states, $E_o$ and $E_u$ are the corresponding energies. In MnPSe$_3$, the Mn$^{2+} 3d^5$ state has only the pure spin $S=5/2$ and leaves no room for orbital degree of freedom. This indicates that the expectation value of the orbital operator is zero, and the first-order perturbation $\langle o | H_{\scriptscriptstyle{\rm SOC}} | o \rangle$ in Eq. (2) has no contribution to the MAE. Then we care about the second-order perturbation energy $|\langle o | H_{\scriptscriptstyle{\rm SOC}} | u \rangle |^2 / (E_o - E_u)$. The MAE caused by the second-order perturbation of SOC is[46-49] \begin{align} {\rm MAE}&= E_{\scriptscriptstyle{\rm SOC}}^{\perp} - E_{\scriptscriptstyle{\rm SOC}}^{\parallel}\notag\\ &\approx \frac{1}{2}\sum_{{i,k}} \left[ \varepsilon_{i}^{\perp}(k) \cdot n_{i}^{\perp}(k) - \varepsilon_{i}^{\parallel}(k) \cdot n_{i}^{\parallel}(k) \right] \cdot \omega(k), \tag {3} \end{align} where $\varepsilon_{i}^{\perp}(k)$ and $\varepsilon_{i}^{\parallel}(k)$ are the $i$th band energy at each $k$-point for MnPSe$_3$ monolayer with the $c$ axis and in-plane spin orientations, respectively; $n_i^{\perp}(k)$ ($n_i^{\parallel}(k)$) and $\omega(k)$ are the corresponding occupation number and the weight of the $k$ point. If one plots the distributions of MAE in the reciprocal space, as shown in Fig. 6(a), the MAE contributions at the $K$ and $K'$ points ($\sim$ $200$ meV in strength) would be far from the real value of the MAE (being here only 0.2 meV as seen above), which seems surprising at a first glance. Actually, the ‘giant’ MAE contributions at $K$ and $K'$ are due to the first-order perturbation of SOC: for systems without an inversion symmetry, the first-order perturbation energy is nonzero at specific $k$-points.[46-48] As the inversion symmetry of MnPSe$_3$ monolayer is broken by its AF magnetic ordering, the MAE contributions at $K$ and $K'$ points are dominated by the first-order perturbation energy, which reaches up to $\sim$ $200$ meV as seen in Fig. 6(a). However, these contributions of first-order perturbation at the specified $k$-points ($k_x$, $k_y$) will be canceled by those at the ($-k_x$, $-k_y$) points with the same spin orientation, or by those at the same $k$-points but with the opposite spin orientations as seen in Fig. 6(b). In other words, the first-order perturbation has no net contribution to the MAE upon an integration over the whole Brillouin zone. Moreover, one would clearly see the contributions of second-order perturbation to the MAE in the reciprocal space via the expression $[(E^{c}-E^{a})+(E^{-c}-E^{-a})]/2$ as plotted in Fig. 6(c). Here the MAE contribution at each $k$-point is indeed tenths of meV in the energy scale, which is in line with the above computed SOC MAE of 241 µeV/f.u.
cpl-40-7-077301-fig6.png
Fig. 6. The distributions of MAE in the reciprocal space: (a) $E^{c}-E^{a}$, (b) $E^{-c}-E^{-a}$, (c) $[(E^{c}-E^{a})+(E^{-c}-E^{-a})]/2$. See more in the main text.
cpl-40-7-077301-fig7.png
Fig. 7. Monte Carlo simulations of the magnetic specific heat of MnPSe$_3$ monolayer. The inset shows the Néel temperatures under different uniaxial strains.
For comparison, we also plot the distribution of the MAE of MnPSe$_3$ monolayer in the FM state over the reciprocal space, see Fig. S4 in the SMs.[41] As the inversion symmetry is preserved in the FM state, the first-order perturbation effect of the SOC is zero at each $k$ points for the MAE,[47] and the distributions of MAE in the reciprocal space show only the contributions from the second-order perturbation effects of the SOC, which have the same energy scale of tenths of meV as in the above AF ground state. To summarize, in MnPSe$_3$ monolayer, SOC has a significantly greater contribution to MAE than the shape anisotropy, and the SOC determines the experimental in-plane spin orientation.[25] In addition, as the $S=5/2$ Mn$^{2+}$ ion in MnPSe$_3$ monolayer has no orbital moment in its ground state, the MAE arises from the second-order perturbation effects of the SOC. Moreover, in order to clearly see the distribution of the MAE over the reciprocal space, one needs to eliminate the first-order perturbation effect of the SOC, which is nonzero at specific $k$ points (but counteracted at $-k$ points) due to the inversion symmetry breaking by the ground state AF order, as shown in Fig. 6. Néel Temperature and Spin Orientation under Uniaxial Strain. We carry out Monte Carlo simulations to estimate the Néel temperature of MnPSe$_3$ monolayer, using the following spin Hamiltonian: \begin{align} \!H=\sum_{i,j} \frac{J}{2} {{\boldsymbol S}_{i}}\!\cdot\!{{\boldsymbol S}_{j}}\!+\! \sum_{i}\Big\lbrace D(S_{i}^{c})^{2}\!+\!E_{n} \big[(S_{i}^{a})^{2}\! -\! (S_{i}^{b})^{2}\big] \Big\rbrace, \tag {4} \end{align} where the first term stands for the isotropic Heisenberg exchange (AF when $J > 0$), the second term $D$ describes the longitudinal ($c$ axis) anisotropy, and the last term $E_n$ represents the transverse ($ab$ plane) anisotropy which would appear under uniaxial strain (see below). Using the above GGA+$U$ computed AF exchange parameter $J=1.41$ meV and the anisotropy parameter $D=31.2$ µeV estimated by the MAE of 195 µeV per $S=5/2$ Mn$^{2+}$ ion ($D$ = MAE/$S^2$), the calculated magnetic specific heat from our Monte Carlo simulations shows $T_{\rm N} = 47$ K for MnPSe$_3$ monolayer (see Fig. 7), which agrees with the experimental 40 K.[25] Strain is an effective tool for tuning various properties of 2D materials. As seen in Fig. 8, our GGA+$U$+SOC calculations find that the AF exchange parameter $J$ is enhanced under uniaxial compressive strains but is reduced under uniaxial tensile strains, the anisotropy parameter $E_n$ shows a similar tendency, but the anisotropy parameter $D$ shows an opposite tendency. Using those parameters, our Monte Carlo simulations show that the $T_{\rm N}$ of MnPSe$_3$ monolayer increases with the compressive strains but decreases with the tensile strains, as seen in the inset of Fig. 7. Moreover, we find that the spin orientation of the Mn$^{2+}$ ions could be tuned by the strains. As shown in Figs. 9(a) and 9(b), when a 2% uniaxial tensile strain is applied to MnPSe$_3$ monolayer, the spin orientation (i.e., the Néel vector) would be along the direction of the uniaxial strain, which is in good agreement with the experiment.[25] In contrast, the Néel vector is perpendicular to the strain direction under 2% uniaxial compressive strain, as seen in Figs. 9(c) and 9(d). When the uniaxial tensile (compressive) strain is along 45$^{\circ}$ with respect to the $a$ axis, the Néel vector would also be parallel (perpendicular) to the strain direction, see Fig. S5 in the SMs.[41] The results of 5% uniaxial strain are similar to the 2% case except that the in-plane magnetic anisotropy is increased from 10 µeV to 25 µeV, see Fig. S6 in the SMs.[41] To summarize, our GGA+$U$+SOC calculations and Monte Carlo simulations well reproduce the experimental $T_{\rm N}$ of the AF MnPSe$_3$ monolayer and the in-plane spin orientation.[25] Moreover, we find that the uniaxial strains can be used to tune the spin orientation, which agrees with the experiment,[25] and we predict that the $T_{\rm N}$ would be increased by a compressive strain.
cpl-40-7-077301-fig8.png
Fig. 8. The AF exchange parameter $J$ and the magnetic anisotropy parameters $D$ and $E_n$ for MnPSe$_3$ monolayer under different uniaxial strains along the $a$ axis.
cpl-40-7-077301-fig9.png
Fig. 9. Polar diagrams of the MAE of MnPSe$_3$ monolayer in the $ab$ plane: under 2% uniaxial tensile strain along 0$^\circ$ (a) and 90$^\circ$ (b) with respect to the $a$ axis, and under 2% uniaxial compressive strain along 0$^\circ$ (c) and 90$^\circ$ (d) with respect to the $a$ axis.
In summary, through density functional theory calculations, magnetic exchange analyses based on MLWFs, and Monte Carlo simulations, we have studied the electronic structure and magnetic properties of MnPSe$_3$ monolayer. Our results show that MnPSe$_3$ monolayer is a CT-type AF insulator. It is proved that, although the direct exchange between two neighboring $S=5/2$ Mn$^{2+}$ ions is certainly an AF type, the experimental AF order is actually dominated by the nearly 90$^\circ$ Mn–Se–Mn superexchange via the strong Mn $3d$–Se $4p$ covalency. We also find that the SOC via the second-order perturbation effects overwhelms the shape anisotropy in yielding the MAE and thus determines the experimental in-plane spin orientation. Moreover, to clearly see the distribution of the MAE over the reciprocal space, one needs to eliminate the first-order perturbation effect of the SOC, which is nonzero at specific $k$ points (but counteracted at $-k$ points) due to the inversion symmetry breaking by the AF order. Using the computed AF exchange parameter and magnetic anisotropy parameters, our Monte Carlo simulations well reproduce the experimental $T_{\rm N}$. Furthermore, we find that a uniaxial strain can be used to tune the spin orientation, which is parallel (perpendicular) to the tensile (compressive) strain direction as observed experimentally, and we also predict that the $T_{\rm N}$ would be increased by a compressive strain. Thus, this work sheds light on the CT-type AF insulating behavior of MnPSe$_3$ monolayer, the major superexchange interactions, magnetic anisotropy and its tuning by strain. Acknowledgments. This work was supported by the National Natural Science Foundation of China (Grant Nos. 12174062, 12241402, and 12104307). References Electric Field Effect in Atomically Thin Carbon FilmsThe rise of grapheneThe electronic properties of grapheneLayer-dependent ferromagnetism in a van der Waals crystal down to the monolayer limitDiscovery of intrinsic ferromagnetism in two-dimensional van der Waals crystalsTunneling Spin Valves Based on Fe3 GeTe2 /hBN/Fe3 GeTe2 van der Waals HeterostructuresMagnetic 2D materials and heterostructuresIntrinsic Van Der Waals Magnetic Materials from Bulk to the 2D Limit: New Frontiers of SpintronicsVan der Waals Heterostructures for High‐Performance Device Applications: Challenges and OpportunitiesNéel-type skyrmions and their current-induced motion in van der Waals ferromagnet-based heterostructuresExotic Dielectric Behaviors Induced by Pseudo-Spin Texture in Magnetic Twisted BilayerTwo-Dimensional Electron Gas in MoSi2 N4 /VSi2 N4 Heterojunction by First Principles CalculationVESTA 3 for three-dimensional visualization of crystal, volumetric and morphology dataMagnetic structure of the quasi-two-dimensional antiferromagnet NiPS 3 Magnetic structure and magnon dynamics of the quasi-two-dimensional antiferromagnet FePS 3 Weak Van der Waals Stacking, Wide-Range Band Gap, and Raman Study on Ultrathin Layers of Metal Phosphorus TrichalcogenidesEmergent superconductivity in an iron-based honeycomb lattice initiated by pressure-driven spin-crossoverMagnetic order in XY-type antiferromagnetic monolayer CoPS 3 revealed by Raman spectroscopyNeutron diffraction study of the layered compounds MnPSe3 and FePSe3Optical absorption and reflectivity study of the layered MnPSe_3 seleniophosphateHalf-Metallicity in MnPSe3 Exfoliated Nanosheet with Carrier DopingPressure-Driven Cooperative Spin-Crossover, Large-Volume Collapse, and Semiconductor-to-Metal Transition in Manganese(II) Honeycomb LatticesGate-Controllable Magneto-optic Kerr Effect in Layered Collinear AntiferromagnetsAntiferromagnet–Semiconductor Van Der Waals Heterostructures: Interlayer Interplay of Exciton with Magnetic OrderingImaging the Néel vector switching in the monolayer antiferromagnet MnPSe3 with strain-controlled Ising orderMagnetic exchange interactions in the van der Waals layered antiferromagnet Mn P Se 3 Strong anisotropy in the mixed antiferromagnetic system Mn 1 − x Fe x PSe 3 Absence of Ferromagnetism or Antiferromagnetism in One- or Two-Dimensional Isotropic Heisenberg ModelsEfficient iterative schemes for ab initio total-energy calculations using a plane-wave basis setFrom ultrasoft pseudopotentials to the projector augmented-wave methodGeneralized Gradient Approximation Made SimpleSpecial points for Brillouin-zone integrationsElectron-energy-loss spectra and the structural stability of nickel oxide: An LSDA+U studyIntrinsic magnetic topological insulators in van der Waals layered MnBi2 Te4 -family materialsNanodevices engineering and spin transport properties of MnBi2Te4 monolayerIsing Magnetism and Ferroelectricity in Ca 3 CoMnO 6 wannier90: A tool for obtaining maximally-localised Wannier functionsMaximally localized Wannier functions: Theory and applicationsThe Monte Carlo MethodBulk and surface electronic structure of MnPSe 3 revealed by photoemission and x-ray absorption spectroscopyTransition Metal Compounds VI 3 : A two-dimensional Ising ferromagnetTriaxial magnetic anisotropy in the two-dimensional ferromagnetic semiconductor CrSBr CrSbSe 3 : A pseudo one-dimensional ferromagnetic semiconductorState-tracking first-principles determination of magnetocrystalline anisotropyMicroscopic origin of magnetocrystalline anisotropy in transition metal thin filmsUnconventional magnetic anisotropy in one-dimensional Rashba system realized by adsorbing Gd atom on zigzag graphene nanoribbonsMagnetic anisotropies of Mn-, Fe-, and Co-doped monolayer MoS2
[1] Novoselov K S, Geim A K, Morozov S V, Jiang D, Zhang Y, Dubonos S V, Grigorieva I V, and Firsov A A 2004 Science 306 666
[2] Geim A K and Novoselov K S 2007 Nat. Mater. 6 183
[3] Castro N A H, Guinea F, Peres N M R, Novoselov K S, and Geim A K 2009 Rev. Mod. Phys. 81 109
[4] Huang B V, Clark G, Navarro-Moratalla E, Klein D R, Cheng R, Seyler K L, Zhong D, Schmidgall E, McGuire M A, Cobden D H, Yao W, Xiao D, Jarillo-Herrero P, and Xu X D 2017 Nature 546 270
[5] Gong C, Li L, Li Z, Ji H, Stern A, Xia Y, Cao T, Bao W, Wang C, Wang Y, Qiu Z Q, Cava R J, Louie S G, Xia J, and Zhang X 2017 Nature 546 265
[6] Wang Z, Sapkota D, Taniguchi T, Watanabe K, Mandrus D, and Morpurgo A F 2018 Nano Lett. 18 4303
[7] Gibertini M, Koperski M, Morpurgo A F, and Novoselov K S 2019 Nat. Nanotechnol. 14 408
[8] Li H, Ruan S, and Zeng Y J 2019 Adv. Mater. 31 1900065
[9] Liang S J, Cheng B, Cui X, and Miao F 2020 Adv. Mater. 32 1903800
[10] Park T E, Peng L, Liang J, Hallal A, Yasin F S, Zhang X, Song K M, Kim S J, Kim K, Weigand M, Schütz G, Finizio S, Raabe J, Garcia K, Xia J, Zhou Y, Ezawa M, Liu X, Chang J, Koo H C, Kim Y D, Chshiev M, Fert A, Yang H, Yu X, and Woo S 2021 Phys. Rev. B 103 104410
[11] Shen Y H, Tong W Y, Hu H, Zheng J D, and Duan C G 2021 Chin. Phys. Lett. 38 037501
[12] Gao R L, Liu C, Fang L, Yao B X, Wu W, Xiao Q L, Hu S B, Liu Y, Gao H, Cao S X, Song G S, Meng X J, Chen X S, and Ren W 2022 Chin. Phys. Lett. 39 127301
[13] Momma K and Izumi F 2011 J. Appl. Crystallogr. 44 1272
[14] Wildes A R, Simonet V, Ressouche E, McIntyre G J, Avdeev M, Suard E, Kimber S A J, Lançon D, Pepe G, Moubaraki B, and Hicks T J 2015 Phys. Rev. B 92 224408
[15] Lançon D, Walker H C, Ressouche E, Ouladdiaf B, Rule K C, McIntyre G J, Hicks T J, Rønnow H M, and Wildes A R 2016 Phys. Rev. B 94 214407
[16] Du K Z, Wang X Z, Liu Y, Hu P, Utama M I B, Gan C K, Xiong Q, and Kloc C 2016 ACS Nano 10 1738
[17] Wang Y G, Ying J J, Zhou Z Y, Sun J L, Wen T, Zhou Y N, Li N N, Zhang Q, Han F, Xiao Y M, Chow P, Yang W, Struzhkin V V, Zhao Y S, and Mao H K 2018 Nat. Commun. 9 1914
[18] Liu Q Y, Wang L, Fu Y, Zhang X, Huang L L, Su H M, Lin J H, Chen X B, Yu D P, Cui X D, Mei J W, and Dai J F 2021 Phys. Rev. B 103 235411
[19] Wiedenmann A, Rossat-Mignod J, Louisy A, Brec R, and Rouxel J 1981 Solid State Commun. 40 1067
[20] Grasso V and Silipigni L 1999 J. Opt. Soc. Am. B 16 132
[21] Li X X, Wu X J, and Yang J L 2014 J. Am. Chem. Soc. 136 11065
[22] Wang Y G, Zhou Z Y, Wen T, Zhou Y H, Li N N, Han F, Xiao Y M, Chow P, Sun J L, Pravica M, Cornelius A L, Yang W G, and Zhao Y S 2016 J. Am. Chem. Soc. 138 15751
[23] Sivadas N, Okamoto S, and Xiao D 2016 Phys. Rev. Lett. 117 267203
[24] Onga M, Sugita Y, Ideue T, Nakagawa Y, Suzuki R, Motome Y, and Iwasa Y 2020 Nano Lett. 20 4625
[25] Ni Z L, Haglund A V, Wang H, Xu B, Bernhard C, Mandrus D G, Qian X, Mele E J, Kane C L, and Wu L 2021 Nat. Nanotechnol. 16 782
[26] Calder S, Haglund A V, Kolesnikov A I, and Mandrus D 2021 Phys. Rev. B 103 024414
[27] Bhutani A, Zuo J L, McAuliffe R D, dela C C R, and Shoemaker D P 2020 Phys. Rev. Mater. 4 034411
[28] Mermin N D and Wagner H 1966 Phys. Rev. Lett. 17 1133
[29] Kresse G and Furthmüller J 1996 Phys. Rev. B 54 11169
[30] Kresse G and Joubert D 1999 Phys. Rev. B 59 1758
[31] Perdew J P, Burke K, and Ernzerhof M 1996 Phys. Rev. Lett. 77 3865
[32] Monkhorst H J and Pack J D 1976 Phys. Rev. B 13 5188
[33] Dudarev S L, Botton G A, Savrasov S Y, Humphreys C J, and Sutton A P 1998 Phys. Rev. B 57 1505
[34] Li J H, Li Y, Du S Q, Wang Z, Gu B L, Zhang S C, He K, Duan W, and Xu Y 2019 Sci. Adv. 5 eaaw5685
[35] An Y P, Wang K, Gong S J, Hou Y S, Ma C L, Zhu M F, Zhao C X, Wang T X, Ma S H, Wang H Y, Wu R Q, and Liu W M 2021 npj Comput. Mater. 7 45
[36] Wu H, Burnus T, Hu Z, Martin C, Maignan A, Cezar J C, Tanaka A, Brookes N B, Khomskii D I, and Tjeng L H 2009 Phys. Rev. Lett. 102 026404
[37] Mostofi A A, Yates J R, Lee Y S, Souza I, Vanderbilt D, and Marzari N 2008 Comput. Phys. Commun. 178 685
[38] Marzari N, Mostofi A A, Yates J R, Souza I, and Vanderbilt D 2012 Rev. Mod. Phys. 84 1419
[39] Metropolis N 1949 J. Am. Stat. Assoc. 44 335
[40] Fujii M, Yamaguchi T, Ohkochi T, De C D, Cheong S W, and Mizokawa T 2022 Phys. Rev. B 106 035118
[41]See the Supplemental Materials for the calculation of the magnetic exchange parameters, the hopping integrals of spin-up channels, the distributions of the MAE in the FM state over the reciprocal space, and the polar diagrams of the MAE under strains.
[42] Khomskii D I 2014 Transition Metal Compounds (Cambridge: Cambridge University Press)
[43] Yang K, Fan F, Wang H, Khomskii D I, and Wu H 2020 Phys. Rev. B 101 100402(R)
[44] Yang K, Wang G, Liu L, Lu D, and Wu H 2021 Phys. Rev. B 104 144416
[45] Wang G Y, Liu L, Yang K, and Wu H 2021 Phys. Rev. Mater. 5 124412
[46] Wang D S, Wu R, and Freeman A J 1993 Phys. Rev. Lett. 70 869
[47] van der Laan G 1998 J. Phys.: Condens. Matter 10 3239
[48] Qin Z Z, Qin G Z, Shao B, and Zuo X 2017 Nanoscale 9 11657
[49] Yue Y L, Jiang C, Han Y L, Wang M, Ren J, and Wu Y K 2020 J. Magn. Magn. Mater. 496 165929