Structures, stabilities, and electronic properties of defects in monolayer black phosphorus

The structures, stabilities, and electronic properties of monolayer black phosphorus (M-BP) with different kinds of defects are investigated within the frame of density-functional theory. All the possible configurations of defects in M-BP are explored, and the calculated results suggest that the stabilities of the configurations with different kinds of defects are greatly related to broken bonds, structural deformation and the character of the bonding. The configurations with two or three vacancies are energetically more favorable than the ones with a single vacancy. Meanwhile, the doping of two foreign atoms, such as sulfur, silicon or aluminum, is more stable than that of the corresponding single dopant. The electronic properties of M-BP are greatly affected by the types of defects. The single S-doped M-BP not only retains the character of a direct semiconductor, but it also can enlarge the band gap by 0.24 eV relative to the perfect one. Such results reveal that the defects not only greatly affect the electronic properties, but they also can be used as an effective way to modulate the band gap for the different applications of M-BP in electronic devices.

Two-dimensional materials, due to their novel electronic and optical properties, have attracted numerous studies for decades. As one of these materials, although graphene exhibits plenty of advantages, and it has dominated the two-dimension research field for many years [1][2][3][4] , the gapless nature of graphene greatly limits its applications in electronic devices. Thus many groups have made efforts to explore other single layer materials, such as boron nitride 5 , or transition metal dichalcogenides (TMDs) [6][7][8][9] .
One of the widely studied semiconductor materials is black phosphorus (BP) [10][11][12][13][14][15][16][17] , which possesses puckered layered structure of one elemental phosphorus allotrope. In bulk BP 18,19 , different layers interact through the weak van der Waals (vdW) force, a configuration which enables one to exfoliate a single layer from the bulk material. Recently, monolayer black phosphorus (M-BP), namely phosphorene, has been fabricated by several research groups 20,21 . Further theoretical calculations show that the band gap of BP is layer dependent with 0.30 eV in bulk BP and 0.90 eV in monolayer BP 21,22 , which is relatively smaller than the corresponding values of 0.35 eV and 2.00 eV 23 in experiments, respectively. Its tunable property of band gap and the super electronic performance largely enhance its promising applications in many electronic devices, such as lithium-ion batteries [24][25][26] , thin-film solar cells 27,28 , channel material in field-effect transistors and thermoelectric devices 29,30 .
Because of its superior electronic performance and structure-dependent novel phenomenon, BP has immediately attracted great attention and becomes one of the most research hotspots. Guo et al. 31 investigated nanotubes of BP and predicted that phosphorene was a promising candidate for future nanoelectronic devices. Several groups 32,33 studied the strain effect on the band gap of BP, and they revealed that BP possessed a rich variety of electronic and structural transformations under different pressures. Gong et al. 34 systematically investigated the possible metal contacts to monolayer BP using DFT method and proposed several future device applications of BP. Zhu et al. 35 did a detailed research on semiconducting layered blue phosphorus, and they indicated that it would be as useful as BP in functional materials. Liu et al. 36  that there were no deep gap states in GB and the deep gap states in the nanoribbons could be eliminated by passivation. Dai et al. 27 reported the effect of stacking order on BP's band gap and predicted its potential applications in thin-film solar cell. Despite all the advances in research on BP, one area that requires further studies is defects formed on the material during the fabrication process. As the defects can't be avoided in the two dimensional BP fabrication process as in other two dimensional materials 37,38 , the electronic property and carrier mobility of M-BP may be affected by the vacancies created by these defects. Meanwhile, foreign atom doping is one of the accessible ways to tailor the electronic properties of M-BP. Thus, it is greatly necessary to examine how the defects affect the structural and electronic properties of M-BP.
In this paper, we carry out the first-principles calculations to explore the intrinsic relationship between the structures, stabilities, and electronic properties of defects (vacancy, sulfur doping, silicon doping, or aluminum doping) in M-BP. Through the examination of formation energies, the most stable configuration with each kind of defect at different defect degrees is determined. The study reveals a number of interesting observations. Firstly, the broken bonds, structural deformation or the character of bonding make divacancy or two dopants M-BP even more stable than the corresponding single one. Secondly, the electronic properties can be modulated by different kinds of defects and their densities. These results provide detailed information for future applications of M-BP.

Computational details
The calculations were carried out within the frame of density functional theory (DFT), as implemented in the Vienna ab-initio simulation package (VASP) 39 . The exchange-correlation potential was described by the generalized gradient approximation (GGA) functional of Perdew, Burke, and Ernzerhof (PBE) 40 . The van der Waals (vdW) correction proposed by Grimme (DFT-D2) 41 was included in our calculations. The plane wave cut-off energy of 400 eV was chosen. A 3 × 3 × 1 k-point sampling based on Gama-centered Monkhorst-Pack scheme 42 was used for the all structural relaxation. The structural optimization process was finished until the force exerted on each atom was less than 10 −2 eV/Å and the total energy change was less than 10 −5 eV. Meanwhile, a vacuum space of 15 Å was added in the unit cell to avoid the interaction between the repeated unit cells. A 4 × 4 × 1 and a 3 × 3 × 1 supercell M-BP was adopted to accommodate vacancies and foreign atoms doping, respectively.
To explore the relative stability of the configurations, we define the formation energies, E fv , of the vacancies M-BP as the following equation, per P E dv is the total energy of M-BP with vacancies, and E per represents the total energy of perfect M-BP. n stands for the number of the removed phosphorus atoms and μ P is the chemical potential of bulk phosphorus. In order to know the relative stabilities of foreign atom substitution, the formation energies of S-doped, Si-doped or Al-doped M-BP are calculated by the following equation: Where E ds is the total energy of M-BP doped with the number m S, Si or Al atoms. E p is the energy of a single P atom in the perfect M-BP. K and m stand for the number of P, doped S, Si or Al atoms in the doped M-BP configuration, respectively. μ b is the chemical potential of bulk sulfur, silicon or aluminum.

Results and discussions
Properties of Perfect M-BP. As shown in Fig. 1, the M-BP is composed of buckled hexagons, along with a rectangular symmetry. The P atoms are subdivided into two sublayers, forming an armchair-like pattern. There are three bonds for each P atom in M-BP: two bonds are formed by the P atoms in the same layer with a distance of d 1 = 2.25 Å, and the third bond is between the interlayer P atoms with a distance of d 2 = 2.22 Å. There are two types of angles between the nearest neighbor P atoms: two angles are 103.75° between P atoms in different layers and the other one is 96.28° between P atoms in same layer. As shown in Fig. 1, the M-BP possesses a direct band gap of 0.88 eV, which is about 0.58 eV larger than the bulk, as other PBE results reported in the previous work 19,43 . The calculated band gap of M-BP is smaller than the experimental value of 2.0 eV 23 , which should be the reason that pure DFT usually underestimates the band gap of semiconductor because of self-interaction error. Both the conduction band minimum (CBM) and the valence band maximum (VBM) are located at the G point for M-BP, a situation which is very beneficial for the electron transport. Charge distributions of VBM are mainly located at the p z orbital of the P atoms in different planes, and those of CBM are mainly located at the p and s orbitals of the P atoms in the same plane 36 , another observation indicating the excellent transport properties.

Atomic structure of defects in M-BP.
To explore the characteristics of the vacancy defects or foreign atom (sulfur, silicon or aluminum) doping systems of M-BP, we consider the possible configurations. Due to the equivalency of all the phosphorus atoms, there is only one kind of single vacancy or foreign atom doping M-BP configuration. As shown in Fig. 2, the configuration with only the P 1 atom removed is similar to the one with only P 2 removed. When the number of vacancies or doping atoms becomes two, there are several possible configurations, such as P 1 -P 2 , P 1 -P 3 , P 1 -P 4 , P 1 -P 5 , P 1 -P 6 , P 1 -P 7 and P 1 -P 8 . The typical defect of P 1 -P 2 vacancy (foreign atom doped) M-BP is shown in Fig. 2. Because of the equivalent of the P 1 -P 4 and P 1 -P 7 , only P 1 -P 4 defect M-BP is considered in the following. And for three vacancies, only four kinds of configurations are considered, e.g., P 1 -P 3 -P 11 , P 1 -P 2 -P 7 , P 1 -P 6 -P 9 , P 1 -P 5 -P 10 vacancies M-BP. For simplicity, the M-BP with P n -P m vacancy is labeled as V n-m , and the M-BP with sulfur (silicon, aluminum) atoms doping is labeled as S n-m (Si n-m , Al n-m ). Here, n and m represent the index of P atoms, which are removed from the perfect one or substituted with the other atoms. For example, the M-BP with the type of P 1 -P 2 vacancy is labeled as V 1−2 M-BP.

The stabilities and electronic structures of M-BP with different vacancies. For single vacancy
in M-BP (V 1 ), an obvious structural deformation occurs. The P atom (P 2 atom) below the vacancy site (removed P 1 atom) moves toward the direction of the vacancy site by a distance of 0.37 Å, forming two bonds with nearby P atoms in the layer of vacancy with a distance of 2.47 Å. The other two with nearby P  The various configurations of the divacancy can be obtained by removing P 1 -P 2 , P 1 -P 3 , P 1 -P 4 , P 1 -P 5 , P 1 -P 6 , P 1 -P 7 , and P 1 -P 8 atoms, respectively. The corresponding doped configurations are obtained by replacing the P atoms by the doping atoms.
atoms are in another layer of P 2 atom with a distance of 2.33 Å. The formation energy of V 1 is 2.00 eV. The ground state of V 1 in M-BP is spin-polarized with a magnetic moment of 0.21 μ B . As shown in Fig. 3(a), the calculated band structures indicates that the valance band of V 1 M-BP is split into two bands near the X, M and G points, and also there is one band across the Fermi level. As shown in Fig 3(a,b), partial DOS shows that the spin polarization is located at the p z orbital, which mainly distributes around the P atoms near the vacancy. When there are two or three vacancies in M-BP, several typical configurations exist. The several typical configurations of divacancy (also two dopants in next section) are shown in Fig. 2. The relative stabilities of these configurations are examined, and the corresponding formation energies are shown in Table 1. Obviously, all divacancy configurations (except V 1−8 ) are more stable than single-vacancy configurations. Furthermore, the V 1-3 vacancies M-BP is the most stable one among all of them, with the formation energy of 0.72 eV/vacancy. The high stability of V 1−3 is greatly correlated with the large structural deformation. When two P atoms (P 1 and P 3 atoms) are removed from the perfect M-BP, the P 11 atom moves towards the P 2 atom, and the distances between the P 11 and P 2 , P 11 and P 10 , P 11 and the P 13 atom below P 11 atom become 2.39 Å, 2.27 Å, 2.25 Å, respectively (see Fig. 4(b)).
Those deformations bond every P atom to three nearby P atoms, forming new sp 3 hybridization by structural deformation, and this should also be one reason that V 1−3 vacancies M-BP is more energetic favorable than the single vacancy case. The newly formed sp 3 hybridization for P atom leads to the delocalization of the valance band state. The impurity state mainly originates from the p z orbital of the P atoms, which is far from the vacancies, as shown in Fig. 4(a,b). In addition, there is also another unoccupied impurity state below the CBM mainly distributing around the defect sites. As a result of the V 1−3 vacancies, the band gap of the defect M-BP becomes an indirect one from Y point to G point with a value of 1.02 eV.  Among all three vacancies in M-BP, the calculated results show that the configuration of V 1−3−11 M-BP, with the three vacancies in the same armchair line, is the most stable configuration. It should be noted that P 3 and P 11 are the equivalent P atoms nearest to the P 1 atom, and configuration V 1−3 is the most stable configuration of divacancy M-BP. The formation energy of V 1−3−11 vacancies is 0.68 eV, which is smaller than that of the two vacancies configuration V 1− 3 . Here, the stability of M-BP with different number of vacancies is discussed. We divide the formation of the vacancies into two parts: the P atoms are removed from perfect M-BP, and the M-BP with vacancies is transformed to a stable state. The corresponding energy changes are broken bond energy and deformation energy, respectively. Once the V 1 ,    Different from one and two vacancies case, V 1−3−11 M-BP exhibits more occupied spin up states than those that spin down below the Fermi level (see Fig. 4(c)), leading to a net spin magnetic moment of 1 μ B . And the magnetic moment mainly comes from the p z orbital of P atoms, which includes the P 2 atom and P atoms around it (see Fig. 2 and Fig. 4(d)). This distribution of spin density of V 1−3−11 is also greatly related with the structural deformation of V 1−3−11 M-BP. The P 4 atom (P 13 atom below P 11 vacancy site) moves toward the P 5 (P 10 ) atom to keep the sp 3 hybridization. However, the P 2 atom below the P 1 vacancy only bonds to the two nearest P atoms, which creates a dangling bond in the P 2 atom. In addition, there is also one unoccupied localized state below the CBM, and the V 1−3−11 M-BP owns an even larger indirect band gap of 1.20 eV (from Y point to G point).
In all, the stabilities and electronic properties of M-BP with different number of vacancies are investigated. The formation energies with different vacancies indicate that the two and three vacancies are more stable than the single vacancy. For example, the formation energy of V 1−3 M-BP is 1.28 eV lower than that of V 1 M-BP, which indicates that the second P atom (P 3 or P 11 ) can easily come out of M-BP once the first P atom (P 1 ) was removed.
The stabilities and electronic structures of sulfur, silicon or aluminum doped M-BP. One of the most effective techniques to tune the electronic property of two-dimensional systems is to dope foreign atoms. Sulfur, aluminum and silicon atoms have one more, two less and one less electron than P atom, respectively. They are expected to induce extra states and tailor the electronic properties of M-BP. In the following, the stabilities and electronic properties of M-BP with S, Si or Al doping are discussed. The electronic structures of M-BP with one dopant, sulfur or silicon, are examined firstly. In the case of the single S-doped M-BP, the calculated band structures are shown in Fig. 5(a). The S-impurity state is in the middle of the band gap, which is split into two states: one occupied spin up, the other unoccupied spin down. Thus the ground state of single S-doped M-BP is spin-polarized with a magnetic moment of 1.0 μ B . Those results indicate that S is not a n-type dopant, as a good n-type dopant should introduce filled states below the CBM closely. The S-doped M-BP remains a direct semiconductor with a band gap of 1.12 eV at G point. The impurity states of the single S-doped M-BP mainly occupy the p z orbital in partial DOS (see Fig. 5(a)). The magnetic moment is mainly localized around the P atom below the doped S atom (see Fig. 5(b)). The localized magnetism around the P atom and the occupied induced state suggest that the unpaired electrons are mainly located at around the P atom below the S atom. The magnetic distribution induced by the S atom should come from the different numbers of electrons between the S and P atoms. Since the S atom has six valance electrons, when replacing one P with S, the S atom bonds to the two P atoms in the same layer, which makes the p orbital of the S atom fully occupied. Thus the P atom below the sulfur atom is unable to form a bond with the S atom, leaving one unpaired electron. Then the magnetic moment comes from the extra valance electron of the S atom, leaving one electron of the P atom unpaired.
The electronic structure of M-BP doped with one silicon atom is also examined. Similar to the S-doped M-BP, the Si-impurity state of Si-doped M-BP is observed in the middle of the band gap and split into two states: partially occupied spin up and down, leading to a spin polarization with a magnetic moment of 0.99 μ B . Obviously, single S-doped and Si-doped M-BP exhibits different electronic structures. The single Si-doped M-BP show two states across the Fermi level, similar to that of V 1 M-BP. As shown in Fig. 5(c), the partial DOS shows that the impurity states mainly occupy the p z orbital, similar to that of S-doped M-BP. The magnetic moment is mainly located around the doped Si atom, as shown in Fig. 5(d). The main reason should come from the fact that silicon has one less electron than phosphorus, and the unpaired electron is mainly located around the silicon atom. In additional, Si is not an effective p-type dopant, as a good p-type dopant should introduce unfilled states above the VBM closely, but Si-doped M-BP does not own this characteristic as discussed above (see Fig. 5(c)).
In order to further understand how the electronic structure is affected by the doping density, M-BP with two dopant atoms is further investigated. Here, we consider the typical configurations with bi-S-doped and bi-Si-doped M-BP. Among all the typical double atom doped configurations, the most stable doping configurations for bi-S-doped and bi-Si-doped M-BP are S 1−2 and Si 1−8 , respectively. The atomic structure of S 1−2 M-BP is explored firstly. As shown in Fig. 6(a), the two doping S atoms are close to each other. Each single S atom needs to bond with two P atoms to form an sp 3 hybridization, while each P atom has to bond with three P atoms to form an sp 3 hybridization in perfect M-BP. Thus the effective utilization of the two bonds of S atom is critical for the most stable bi-S-doped M-BP. When the two S atoms are doped into the neighbor sites of the M-BP, every S atom can form two bonds with the P atom whereas there is no bonding between the S atoms. Meanwhile the distance between the S atoms in S 1−2 M-BP is enlarged to 3.31 Å to avoid the repulsive interaction between them.
Nevertheless, the most stable configuration of bi-Si-doped M-BP is Si 1−8 M-BP, which is greatly different from the one of bi-S doping, as shown in Fig. 6(b). The two doping atoms do not like to stay closely. The difference between these two dopants mainly comes from the distinct electronic structures. In addition, it should be noted that the deformation of the bi-Si-doped M-BP configurations Si 1−8 is rather large. The Si atoms move toward to the P atoms below them. Meanwhile, the P atoms move in the same direction, even breaking the two-layer structure in the process. The distance between the two Si atoms is 2.63 Å, which is obviously shorter than the one between the P atoms at the corresponding sites in perfect M-BP (3.50 Å).
We further extend our study to discuss the stabilities of the S-doped and Si-doped M-BP. The calculated formation energies of different doped M-BP are shown in Table 1. The doping of the first S atom into M-BP needs extra energy of 0.09 eV. After doping of the second S atoms, the formation energy becomes negative, which suggests that the second S atom is easier to dope. Different from those of the In the following, the electronic properties of the most stable configurations of bi-doped M-BP are also explored, as shown in Fig. 7. Unlike single S-doped or Si-doped configurations, S 1−2 and Si 1−8 M-BP do not exhibit any magnetism. As shown in Fig. 7(a), there is an impurity state at the top of the valence band. The impurity electronic state of S 1−2 M-BP is mainly delocalized around the doped-S atom. The band structures of S 1−2 M-BP show that the doped M-BP becomes an indirect semiconductor with a band gap of 1.18 eV, which is about 0.3 eV larger than the perfect one. We also take this as an example to explore how the supercell size affects the band gap. The band gap of S 1−2 M-BP with 4 × 4 × 1 supercell is 1.00 eV, about 0.18 eV smaller than that with a 3 × 3 × 1 supercell. This could be understood by the impact of the density of the S 1−2 type dopant, i.e., the smaller density of S 1−2 type dopant and the smaller band gap. This indicates the density of dopants could impact the electronic structures of M-BP in another way.
The band gap of Si 1−8 M-BP is an indirect one with a value of 0.99 eV. Si 1−8 M-BP exhibits two impurity states with one occupied near the Fermi level and another unoccupied below CBM (see Fig. 7(d)). The electronic states of Si 1−8 M-BP impurity states are shown in Fig. 7(e,f), respectively. The occupied impurity state mainly distributes at the P atoms in two lines, which contains the silicon atoms, whereas the unoccupied impurity state mainly localizes at the silicon atoms.
Finally, the stabilities and electronic properties of M-BP with Al doping are also discussed. As listed in Table 2, the formation energies indicate that the Al 1−4 M-BP is the most stable configuration of single-Al and bi-Al doped case. Furthermore, the formation energies indicate that the Al atoms are easier to dope into M-BP than the sulfur atoms. Moreover, the electronic properties of A 1 and A 1−4 M-BP are also explored, as shown in Fig. 8. As the Al atom has two less electrons than the P atom, there will be one fewer occupied state once one Al atom is doped into M-BP, as shown in Fig. 8(a,c). The unoccupied impurity state of Al 1 M-BP is below the VBM (Fig. 8(a)) and it is mainly localized around the Al atom as shown in Fig. 8(b). The band structures of Al 1 M-BP show that the doped M-BP still nearly remains a direct semiconductor with a band gap of 1.15 eV, locating at Y point. Al 1−4 M-BP has a smaller band gap, about 1.12 eV. There are also two unfilled impurity states emerging below the CBM (Fig. 8(c)), mainly localizing around the two Al dopants (Fig. 8(d)).

Conclusion
In summary, the structures, stabilities and the electronic properties of vacancy, S-doped, Si-doped, Al-doped M-BP are investigated by first-principle calculations. The most stable configurations of vacancy and doping are explored. The results reveal that the stabilities of the configurations are mainly determined by broken bonds, the structural deformation and the character of bonding. The divacancy or bi-doped configurations are more stable than single-vacancy or single-doping, and the degree of difficulty for doping atom into M-P is from Al to Si. All kinds of defects discussed could tune the electronic properties of M-BP. For single S or Si atom doping, the spin density is mainly localized around the doping atom, and there are impurity states among the band gap, occupying the p z orbital and in the middle of the band gap. Our results indicate that the band gap could be modulated by introducing vacancies or dopants. Those results should be useful for the design of novel electronic devices with M-BP.