Spin–orbit coupling in buckled monolayer nitrogene

Buckled monolayer nitrogene has been recently predicted to be stable above the room temperature. The low atomic number of nitrogen atom suggests, that spin–orbit coupling in nitrogene is weak, similar to graphene or silicene. We employ first principles calculations and perform a systematic study of the intrinsic and extrinsic spin–orbit coupling in this material. We calculate the spin mixing parameter b2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$b^2$$\end{document}, reflecting the strength of the intrinsic spin–orbit coupling and find, that b2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$b^2$$\end{document} is relatively small, on the order of 10-6\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$10^{-6}$$\end{document}. It also displays a weak anisotropy, opposite for electrons and holes. To study extrinsic effects of spin–orbit coupling we apply a transverse electric field enabling spin–orbit fields Ω\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Omega$$\end{document}. We find, that Ω\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Omega$$\end{document} are on the order of a single μ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu$$\end{document}eV in the valence band, and tens to a hundred of μ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu$$\end{document}eV in the conduction band, depending on the applied electric field. Similar to b2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$b^2$$\end{document}, Ω\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Omega$$\end{document} is also anisotropic, in particular for the conduction electrons.

parameter b 229 , while the extrinsic SOC, arising due to breaking of the space inversion symmetry, is quantified by the amplitude of spin splittings so and spin-orbit fields . All these quantities provide essential information about the strength and anisotropy of SOC in the band structure.
We find, that b 2 is of the order of 10 −6 , both for electrons and holes, and displays a weak in-plane to out-of plane anisotropy. The extrinsic SOC shows significant diversity between the valence and conduction bands. In the former, the values of spin-orbit fields vary between a few to a dozen of µ eV for the considered amplitudes of an external electric field, and are rather isotropic. In the latter, the values of are roughy ten times bigger and are strongly anisotropic.

Results and discussion
Intrinsic spin-orbit coupling. We begin our considerations from the examination of the electronic band structure of nitrogene. The calculated first principles relativistic band structure is shown in Fig. 1. Similar to other buckled 2D materials made of group V elements, such as blue phosphorene, arsenene or antimonene 22 , the band gap of monolayer nitrogene is indirect. The size of the indirect gap for the PBE exchange correlation functional is E g = 4.06, eV, and is in agreement with the values reported by others 13 . The top-most valence band has two local maxima lying between the Ŵ K and Ŵ M paths in the First Brillouin Zone (FBZ) (Fig. 1c). The maxima are located approximately 0.5 eV above the saddle point at the Ŵ point, and differ in energy by 35 meV. The edge of the conduction band is located in the middle of the Ŵ M path in the FBZ.
As already stated by Lee et al. 10 , the effects of spin-orbit coupling on the band structure on buckled nitrogene are weak. Indeed, the inclusion of relativistic effects in the calculation makes no significant effects on the non-relativistic band structure. The most prominent ones are spin-orbital gaps opened at high symmetry points of the FBZ. In the valence band, the spin-orbital gap at the Ŵ point (between the bands marked 2 and 3) is � Ŵ SO = 17.3 meV, while at the K point (between the bands marked 1 and 2') it is K SO = 2.1 meV. For comparison, the same splittings for graphene are � Ŵ SO = 9 meV and K SO = 24 µeV 23 , and for blue phosphorene � Ŵ SO = 48 meV, K SO = 10 µeV 22 . To characterize the intrinsic SOC in the band structure away from the high symmetry points we calculate the spin mixing parameter b 2 . This parameter measures the amplitude of the spin component being admixed to the Bloch state of opposite spin by the SOC. Importantly, the parameter b 2 can be easily accessed experimentally from the Elliott relation connecting b 2 with the spin τ s and momentum τ p relaxation times: b 2 = τ p (4τ s ) −1 , provided the Elliott-Yafet mechanism dominates spin relaxation 29,30 . This is usually the case, when the spin lifetime follows the same characteristics as the momentum lifetime [31][32][33][34] . Knowing τ p and τ s from the experiment, one can extract the sample independent parameter b 2 , and compare it with the theoretical values. Recently we have successfully applied this strategy to characterize b 2 and spin relaxation in black phosphorus 12,34 . In numerical simulations the parameter b 2 can be calculated directly from the wave functions, provided the spin subbands of a Bloch state ψ k,n are degenerate. This requirement is met if the time reversal symmetry and space inversion symmetry of the sample are present simultaneously. In such a case the two Bloch subbands are or, alternatively by calculating the deviation of the expectation value of the spin operator ŝ i from its nominal value 1/2 (in units of ) 36 where σ = {⇑, ⇓} and ψ σ ,i k,n (r) is the eigenstate of ŝ i . For normalized states 0 ≤ b 2 k,n,i ≤ 0.5 , where b 2 k,n,i = 0.5 corresponds to fully spin mixed states and b 2 k,n,i = 0 to pure spin up (down) spinors. The calculated Fermi contour averaged spin mixing parameter b 2 plotted versus the position of the Fermi level is shown in Fig. 2. The value of b 2 is in the range of 10 −6 for both, the valence and the conduction band. This is roughly ten times the value of b 2 for graphene and ten times less than for black and blue phosphorus 22,34,37 . At E F ≈ 32 meV the second valence band maximum (along the ŴM path in Fig. 1) crosses the Fermi level and k-points from this wedge of the FBZ start contributing to b 2 , slightly modifying its slope. The parameter b 2 is almost doping independent. Such behavior is typical for bands being energetically separated from other lower and higher lying bands 22 .
The parameter b 2 shows weak anisotropy with respect to the spin quantization axis. In the valence band ( Fig. 2a) b 2 for SQA=z (out of plane polarization) is roughly twice as large as for SQA=x, This result is surprising, since most of buckled elemental monolayer materials display the opposite trend 22 . In this context nitrogene is similar to graphene, which exhibits similar anisotropy of b 2 for holes away from the Dirac point. For conduction electrons (Fig. 2b) we observe the opposite trend in b 2 , namely, b 2 SQA=x/y /b 2 SQA=z ≈ 1.6 . Although the anisotropy of b 2 is not very high, this result deviates from other buckled pnictogens, for which b 2 is mostly isotropic.
The results presented Fig. 2 have been obtained for the PBE exchange-correlation functional. Since for hybrid functionals the band gap increases by 2 eV 13 , we have checked how the band gap correction influences the values of the spin mixing parameter and performed calculations for the HSE 38 functional. Although the band gap increased to 6.34 eV (see Supplementary Fig. S1) the corresponding spin mixing parameter ( Supplementary  Fig. S2) stays almost unaffected and varies by at most a few percent. Therefore, the results obtained for the PBE functional can be taken as conclusive.
Extrinsic spin-orbit coupling. Monolayer nitrogene has been predicted to be structurally stable on metal surfaces. Even though the interaction with the substrate is weak and makes no significant changes to the band structure of nitrogene 13 , the crystal potential at the interface breaks the inversion symmetry of the nitrogene lattice and enables the extrinsic Bychkov-Rashba spin-orbit coupling 28 . The extrinsic SOC has two main effects on the electron spin. First, it removes the degeneracy of spin states, and second, it induces crystal momentum where is the Planck constant, and σ is the vector of Pauli matrices. Instead of placing nitrogene on a particular substrate we model its presence by applying an external transverse electric field E E E = (0, 0, E) , whose amplitude can be precisely controlled in numerical calculations. This approach, allowing us to simulate different substrates, is justified due to weak hybridization of states of nitrogene and the substrate 13 .
In Fig. 3 we show the Fermi contour averaged spin splitting so calculated for several values of the electric field. In contrast to the spin mixing parameter b 2 , so differs significantly between the valence and conduction band. In the former (Fig. 3a), so takes the values from a few to a dozen of µ eV for the considered values of electric fields and grows by approximately 3 µ eV per 1Vnm −1 . In the conduction band (Fig. 3b) the corresponding values are almost ten times bigger than in the valence band, and vary between 20 µ eV to 100 µeV, for E = 1 Vnm −1 and E = 4 Vnm −1 respectively, with a linear growth of 33 µ eV per 1 Vnm −1 . Similarly to b 2 , the slope of so in the valence band changes slightly when the Fermi level E F ≈ 32 meV and the second valence band maximum (between the Ŵ and M points in Fig. 1) enters the Fermi contour.
The spin splitting so characterizes the strength of the extrinsic SOC in a band. The anisotropy of SOC can be accessed through the components of the spin-orbit field k , k,i , where S i , i = {x, y, z} is the expectation value of spin one-half operator at a given k point, and S k = S 2 k,x + S 2 k,y + S 2 k,z 39 . The calculated Fermi surface averaged components i are presented in Fig. 4. In the valence band (Fig. 4a,b), x and y take similar values, while z (Fig. 4c) is almost twice smaller and displays different monotonicity. In the conduction band, the in-plane spin components x , y (Fig. 4d,e) are large, tens of µeV, and are doping independent, while z (Fig. 4f) takes the values of a few µeV, and is doping dependent. Such big differences between � x/y and z lead to a sizeable, up to � x/y /� z ≈ 40 , doping-dependent anisotropy.
To understand these results one needs to look at the components of the electron spin which shape i . We show them in Fig. 5a,b. The in-plane spin components S x and S y form the typical circulating Rashba spin texture (Fig. 5a). Within the considered doping range the length of S x and S y is approximately constant, as shown in Fig. 5c. The S z component displays much bigger diversity, and takes small values in the center of FBZ, while for bigger crystal momenta we observe the spin-valley locking effect with a strong spin polarization. In effect, when doping increases the average value of |S z | in the valence band first decreases from the value |S z | ≈ 0.1 to |S z | ≈ 0.03 , and for E F ≤ −32 meV its starts increasing and saturates at the value |S z | ≈ 0.075 (Fig. 5c). The qualitative change to S z takes place when the Fermi contour reaches the k-points close to the anticrossings marked by the green rectangle in Fig. 1a, what happens exactly at E F ≈ −32 meV. For k in the range from the Ŵ point to the anticrossing at E − E F ≈ −10 eV, S z ≈ 0 (see the Supplementary Fig. S3). Increasing k towards the K-point, S z starts growing and reaches the maximum S z ≈ 0.5 at k above the anticrossing lying at E − E F ≈ −9 eV. Since the valence band maximum lies in between of the two anticrossings, we first observe a decrease of |S z | followed by its increase.
In the conduction band (Fig. 5b,d), the in-plane spins also form the Rashba texture, similar to the valence band, but the z component of spin is very small in the wedge of the BZ corresponding to doping range  www.nature.com/scientificreports/ (represented by black ellipses). In effect, x and y take the values close to so , while z is of the order of a few µ eV (Fig. 4d-f), generating a large anisotropy of the extrinsic SOC. As can be seen from the results shown above, SOC in nitrogene is weak. The calculated spin mixing parameter, b 2 ≈ 10 −6 , is approximately ten times bigger than for graphene 22 . This suggests, that the intrinsic SOC should not contribute much to spin scattering in nitrogene. Taking the typical momentum lifetime τ p = 100 fs, we can make a rough estimate of spin lifetime from the Elliott-Yafet mechanism 29,31 ,τ EY s ≈ τ p /4b 2 ≈ 50 ns. On the other hand, several factors may affect SOC and spin dynamics in a material. For instance in graphene, the out-of-plane lattice distortions strongly affect both, the strength of the SOC and momentum scattering, dramatically reducing spin lifetime [40][41][42][43] . A transition from the flat to rippled graphene results in the emergence of an additional (intrinsic) SO term H curv 42,44 . For typical ripples of radii in the range R ∼ 50 nm-100 nm the characteristic energy of this interaction is curv ≈ 0.2 K, and can exceed the energy of the intrinsic SOC, int ≈ 0.01 K − 0.1 K 23,42 , making the relevant spin-flip processes important. Additionally, rippling breaks the mirror symmetry ( z → −z ) protecting the electron spin, and enables additional spin relaxation channels by flexural phonons and random spin-orbit fields 41,43 .
In contrast to graphene, nitrogene is naturally buckled. Thus, the effects of buckling are embedded in the intrinsic SOC and b 2 , and are much stronger than effects of rippling discussed above. It was shown, that int at the K-point in graphene grows quadratically with the buckling height δ , � int ∼ (δ/a) 2 , a being the lattice constant, and for δ/a ≈ 0.08 it jumps to int ≈ 1 K 23 . In the case of nitrogene, the buckling height is δ ≈ 0.9 Å, a = 2.3 Å, which gives δ/a ≈ 0.4 . Following the quadratic dependence for graphene, one gets for nitrogene int = 25 K, close to the value K so = 23 K extracted from our first principles calculations. For comparison, static ripples of radii R ∼ 50 nm-100 nm, give the correction to the intrinsic SOC of the order of � curv ∼ a/R ≈ 0.2 K 42 ,-negligible in the case of nitrogene.
Although lattice ripples should not significantly affect the intrinsic SOC in nitrogene, they can generate a small random spin-orbit fields leading to faster spin decoherence, in a similar way as it takes place graphene 41 . Fortunately, lattice rippling can be to a large extent eliminated by encapsulation of the host layer by, e.g., hexagonal boron nitride 45,46 .
spin-orbit coupling plays an important role in contemporary solid state physics, spintronics and topological quantum computing. Apart from these areas of physics, it has also been intensively studied in cold atoms systems. Therefore, in the next few lines we briefly compare SOC in crystalline solids with SOC in Bose gases.
In crystalline solids, SOC originates from the crystal potential, in which an electron is moving. In the rest frame of the electron, the electric field induced by the crystal potential acts as an effective, momentum dependent Zeeman field acting on the electron spin. The strength of SOC in a band is determined by the chemical www.nature.com/scientificreports/ composition of the crystal (materials made of heavier elements display stronger SOC effects than those made of light elements) and the topology of the band structure, and to some extent may be modified externally, e.g., by electric fields 47 , strain 48 , proximity effects 49 or twisting 50 . The form of the spin-orbit interaction is dictated by the symmetry of the crystal; for instance, the famous Bychkov-Rashba 28 and Dresselhaus 51 types of SOC result from broken structure and bulk inversion symmetry, respectively. In cold bosons systems, SOC is realized by coupling the motion of an atom to its internal (hyperfine) pesudospin states, corresponding to the electron's spin up and spin down states [52][53][54][55][56] . Implementation of such coupling, called synthetic SOC 54 , requires an extra effort of dressing atomic states with lasers, but offers a full control of this interaction at will 57 . By a proper combination of laser fields and atomic pseudo-spin states, a variety of SOCs can be created and dynamically modified 54,57,58 . For example, Lin et al. 54 , realized one dimensional SOC, corresponding to equal contributions of Bychkov-Rashba and Dresselhaus SOC in conventional systems; Wu et al. 57 implemented a scheme allowing for a controllable transition from 1D to 2D SOC. More exotic forms of SOC having no counterparts in real materials, such as, a 3D analogue of Rashba SOC 59,60 , are also possible, making cold atoms systems a powerful platform for exploring spin-orbit coupling and many body physics. Like in conventional materials, in bose gases, SOC is essential to the emergence of fascinating physical phenomena, e.g., a degenerate ground state of spin-orbit coupled Bose-Einstein condensates 61

Conclusions
We have investigated the fundamental spin-orbit coupling in buckled monolayer nitrogene. Based on first principles calculations we found that spin mixing parameter characterising the intrinsic SOC is small, on the order of b 2 ≈ 10 −6 , and displays weak anisotropy. The extrinsic SOC, characterized by the Rashba spin-orbit fields , is also weak, on the order of µ eV in the valence and tens to a hundred of µ eV in the conduction band. Similar to the intrinsic SOC, the extrinsic SOC is also anisotropic. The anisotropy is particularly strong in the conduction band. Weak spin-orbit coupling in nitrogene suggests, that doped nitrogene or nitrogene nanoribbons may be attractive materials for spintronics applications.

Methods
First-principles calculations were performed with the Quantum Espresso package 70,71 . The norm-conserving pseudopotential with the Perdew-Burke-Ernzerhof (PBE) 72,73 version of the generalized gradient approximation (GGA) exchange-correlation functional was chosen, taking the kinetic energy cutoff of the plane-wave basis 80 Ry for the wave function and 320 Ry for charge density respectively. These parameters were found to give well converged results. Calculations with a hybrid functional were done with the Heyd-Scuseria-Ernzerhof (HSE06) 38 functional, with the Fock exchange contribution 20%. Monolayer nitrogene was placed in a vacuum of 21 Å to minimize interactions between periodic copies of the system. Self-consistency was achieved with 21k x × 21k y × 1k z Monkhorst-Pack grid 74 . The optimized lattice constant a was determined by minimizing the total energy followed by fitting a parabolic function. In each step the positions of atoms were fully relaxed by the quasi-Newton scheme as implemented in Quantum Espresso, assuming force convergence threshold 10 −4 Ry/ bohr. We found a = 2.29 Å and the out of plane buckling of the lattice δ = 0.96 Å, similar to those reported by other authors 13 . Calculations with the electric field were carried out with the dipole correction 75 . Fermi contour averages of spin mixing parameter b 2 , spin-orbit field 2 and spin components s i were calculated according to the formula where A k stands for b 2 k , 2 k and s k,i , S BZ is the area of the Fermi surface, ρ(E F ) is the density of states per spin at the Fermi level, v F (k) is the Fermi velocity.