Electronic Bandgap and Edge Reconstruction in Phosphorene

Oct 24, 2014 - ‡Center for Nanophase Materials Sciences and §Computer Science and Mathematics Division, Oak Ridge National Laboratory, Oak Ridge, T...
7 downloads 11 Views 5MB Size
Letter pubs.acs.org/NanoLett

Electronic Bandgap and Edge Reconstruction in Phosphorene Materials Liangbo Liang,† Jun Wang,‡ Wenzhi Lin,‡ Bobby G. Sumpter,‡,§ Vincent Meunier,*,† and Minghu Pan*,‡,∥ †

Department of Physics, Applied Physics, and Astronomy, Rensselaer Polytechnic Institute, Troy, New York 12180, United States Center for Nanophase Materials Sciences and §Computer Science and Mathematics Division, Oak Ridge National Laboratory, Oak Ridge, Tennessee 37831, United States ∥ School of Physics, Huazhong University of Science and Technology, Wuhan, Hubei 430074, China ‡

S Supporting Information *

ABSTRACT: Single-layer black phosphorus (BP), or phosphorene, is a highly anisotropic two-dimensional elemental material possessing promising semiconductor properties for flexible electronics. However, the direct bandgap of single-layer black phosphorus predicted theoretically has not been directly measured, and the properties of its edges have not been considered in detail. Here we report atomic scale electronic variation related to strain-induced anisotropic deformation of the puckered honeycomb structure of freshly cleaved black phosphorus using a high-resolution scanning tunneling spectroscopy (STS) survey along the light (x) and heavy (y) effective mass directions. Through a combination of STS measurements and first-principles calculations, a model for edge reconstruction is also determined. The reconstruction is shown to self-passivate most dangling bonds by switching the coordination number of phosphorus from 3 to 5 or 3 to 4. KEYWORDS: Phosphorene, scanning tunneling microscopy/spectroscopy, direct bandgap, monatomic step edges, density functional theory, self-passivation

A

addition, we study a monostep of layered black phosphorus. The step introduces nontrivial edge states, exhibiting multiple DOS peaks in STS, which can be explained as a signature of one-dimensional (1D) Van Hove singularities corresponding to a reconstructed edge structure. Scanning tunneling microscopy (STM) measurements were performed on freshly cleaved surfaces of black phosphorus crystals inside an ultrahigh vacuum (UHV) chamber at room temperature in order to avoid sample degradation and placed in the microscope stage overnight to cool down at 80 K, as described in Materials and Methods. Figure 1a shows a representative STM topographic image of cleaved BP surface with a scanning range measured on an 18 nm × 18 nm surface area. 1D zigzag atomic rows are observed on the surface layer. In addition to an atomically resolved cleaved plane, some defects appear at the fresh-cleaved surface, as shown in Figure 1a and can be understood as either an impurity (as a bright spot) or a structural vacancy (as a dark pit).

n increasing number of alternative two-dimensional (2D) materials are being explored in the “post-graphene age”,1 such as transition metal dichalcogenides (TMDs), silicene, and germanene. These materials are studied with the hope of overcoming graphene’s deficiencies such as its zero bandgap. Among these 2D layered materials, layered black phosphorus (phosphorene) is expected to exhibit superior mechanical, electrical, and optical properties due to an intrinsic and tunable bandgap. It is, besides graphene, the only stable elemental 2D material that can be mechanically exfoliated.2,3 Its high hole mobility and direct semiconducting bandgap4 fuel hope for the development of new electronic devices in the postsilicon era. However, before phosphorene can be integrated in everyday electronics, many aspects of the material still remain to be elucidated. For example, a direct bandgap of single layer black phosphorus predicted theoretically has not been directly measured in experiments, due in part to the difficulty of isolating a single layer of phosphorene. Furthermore, the properties of its edges, which are likely different from those deep in the basal plane, have not been considered in detail. Here, we present a detailed scanning tunneling microscopy/ spectroscopy (STM/S) investigation, corroborated by ab initio calculations, to reveal the presence of a semiconducting 2 eV bandgap observed by STS, which represents the direct bandgap predicted for an isolated surface layer of the cleaved solid. In © XXXX American Chemical Society

Received: July 28, 2014 Revised: September 30, 2014

A

dx.doi.org/10.1021/nl502892t | Nano Lett. XXXX, XXX, XXX−XXX

Nano Letters

Letter

By comparing the atomic-resolved STM images in Figure 1b with the structural model of BP (Figure 1c), we find that P atoms in the topmost P atomic chains show different contrasts in the STM images. The P1 atoms are brighter in constantcurrent STM images than the P2 atoms. With no external influence, P1 and P2 atoms should be equivalent according to the symmetry of BP. Our density functional theory (DFT) calculations confirm that P1 and P2 have the same height in the defect-free and strain-free single-layer and few-layer BP. Hence, for the BP surface without height variation between P1 and P2, the simulated constant current STM shows that P1 and P2 exhibit the same brightness regardless of the applied bias. Further, the apparent distortion of the lattice, perpendicular to the P1−P2 direction observed experimentally, is not compatible with the full symmetry of the intrinsic lattice. Clearly, the symmetry has to be broken to account for the experimental results. As discussed below for the measured electronic inhomogeneity, the surface is indeed under strain. Our DFT calculations show that the height difference between P1 and P2 is proportional to the applied compressive strain. For the monolayer BP under 2% biaxial compressive strain, a geometrical 0.06 Å height variation is induced between P1 and P2. Subsequently, P1 becomes brighter in the simulated STM image than P2 (right panel in Figure 1c), consistent with the experimental STM image in Figure 1b. Bulk black phosphorus has a small semiconducting bandgap around 0.3 eV. Unlike graphene,6,7 a monolayer (1L) of black phosphorus (i.e., “phosphorene”) is a semiconductor with a predicted direct bandgap of ∼2 eV at the Γ point of the Brillouin zone,8,9 similar to monolayer transition metal dichalcogenides MoS2 and WS2.10−15 For few-layer BP, interlayer interactions reduce the bandgap for each added layer, and the bandgap eventually decreases to ∼0.3 eV corresponding to the bulk.16−19 The direct gap also moves to the Z point as a consequence of the presence of additional layers.20 Our calculations at both DFT and many-body GW levels reveal a similar trend of the bandgap of BP with the number of the layers, as shown in the right panel of Figure 2b. From 1L to bulk BP, the GW (DFT) bandgap monotonically decreases from 1.94 (0.82) eV to 0.43 (0.03) eV. As usual, DFT

Figure 1. Atomic resolved STM image of a freshly cleaved surface of black phosphorus. (a) Topographic STM image of the sample. The setup conditions for imaging were a sample-bias voltage of 1.0 V and a tunneling current of 0.12 nA. The scan size is 18 nm × 18 nm. Atomic resolved image shows a one-dimensional puckered lattice. (b) Zoomed-in image of the region corresponding to the red boxed area in a. (c) BP surface (left panel) and simulated constant current STM image (right panel) of the surface under 2% biaxial compressive strain at the bias of 1.0 V. A 0.06 Å height variation between P1 and P2 atoms is induced by the strain.

Zooming on the flat area marked by the red rectangle in Figure 1a, the close-up image (Figure 1b) shows a higher resolution STM image of BP at a few nanometer scale. According to the atomic structural model of BP (Figure 1c), as an upper atom (in darker color) sits almost directly above a lower one (in lighter color) in a puckered layer, it is hard to detect the lower atoms using STM. Therefore, the STM image of the BP surface is made up of an array of zigzag rows composed only of the upper atoms. All atoms almost stay in their original sites, as shown in this atomic model. The in-plane lattice constants along the “x” direction of a lighter effective mass and along the “y” direction of a heavier effective mass measured from the STM image are 4.4 and 3.4 Å, respectively, in good agreement with the reported values of bulk BP and previous STM results.5

Figure 2. Spectroscopic observation of the bandgap for black phosphorus. (a) Two representative tunneling spectra in the log scale measured on the black phosphorus surface show a wide bandgap with estimated size 2.05 eV. Red and blue dotted curves were measured at different locations, marked by red and blue stars on the figure inset, exhibiting the variation of measured bandgaps. (b) Left panel: calculated density of states of single-layer (1L) and bulk BP using the many-body GW approach. Right panel: the thickness dependence of calculated bandgaps of BP at both DFT and GW levels. (c) High-resolution STM image (Vbias = +1.2 V, Iset= 150 pA) with scan size of 2.4 nm × 3.6 nm. The right panel is spatial profile of LDOS at the energy of +1.15 eV. All spectral surveys were taken with a sample-bias voltage of 1.2 V, a tunneling current of 0.15 nA, and bias modulation amplitude of 3−5 mVrms at 80 K. B

dx.doi.org/10.1021/nl502892t | Nano Lett. XXXX, XXX, XXX−XXX

Nano Letters

Letter

Figure 3. Edge and edge states of a freshly cleaved surface of black phosphorus. (a) Topographic STM image of a step edge with a sample-bias voltage of 1.0 V and a tunneling current of 0.12 nA. The scan size is 24 nm × 100 nm. (b) Line profile showing the height of an atomic step edge, measured at about 5.3 Å. This indicates that the STM is visualizing a single layer edge of BP. (c) A typical dI/dV spectrum measured at the step edge, shown in part a. A series of peaks (labeled as P1 to P4) near the negative gap edge, indicating the appearance of edge states. (d) A series of dI/dV spectra measured along the “x” direction across the step edge. These unprocessed spectra were generated from a single spectral survey that was taken along the red dashed line, marked in part a. The right panel is a color plot of these dI/dV curves with the relative distance. The evolution of edge states inside of the bandgap is clearly shown in this plot. All spectral surveys were taken with a sample-bias voltage of 1.2 V, a tunneling current of 0.15 nA, and a bias modulation amplitude of 3−5 mVrms at 80 K.

systematically underestimates bandgaps of BP, while GW yields much better results compared to experiment. Furthermore, such unique thickness dependence of the BP’s bandgap indicates that it can be used as an indicator of the material thickness. At the same time, it is well-known that STM/S is a surface-sensitive technique which is only sensitive to the topmost surface atomic layer of materials. More precisely, STM probes the local density of states a few angstroms above the topmost surface. Performing tunneling spectroscopy enables the measurements of local electronic structure and bandgap of BP material, which may consist of only the contribution from the topmost surface layer to the density of states above the surface. Typical tunneling spectra (dI/dV versus V) collected on clean BP surface at locations marked by blue and red stars (in the inset of Figure 2a) represent the unperturbed electronic structure of pristine BP (see the blue and red curves in Figure 2a). The U-shaped spectra in the log scale show an intrinsic bandgap of 2.05 eV, instead of a predicted ∼0.3 eV for bulk BP. We repeated these measurements many times within different scan areas and tip states, and the dI/dV spectra which show ∼2 eV gap are very reproducible. According to our GW calculations (Figure 2b) and previous theoretical works,8,9 such a large bandgap can only occur for monolayer BP (GW gap 1.94 eV). From 2L to bulk BP, the GW gap monotonically reduces from 1.65 to 0.43 eV. Clearly, the intrinsic ∼2 eV gap measured by STS shows that it only detects the topmost layer of the BP surface. Our DFT calculations indicate that, from the surface (set at z = 0 Å) to vacuum, DOS peaks contributed by bottom layers decay much more rapidly compared to those from the top layer. At z = 3 Å, the former peaks display negligible amplitudes compared to the latter peaks. Since the STS measurements are usually performed around z = 3 Å and higher, the DOS peaks originating from bottom layers in the dI/dV curve are too weak to be detectable or are well hidden in the background noise. Therefore, the STS can only detect the

DOS peaks from the top layer and its gap (more details are provided in Figure S1 and Figure S2 in the Supporting Information). The tunneling spectra vary slightly with the spatial locations. A significant change in spectra happens at the gap edge in positive bias. Such electronic gap inhomogeneity is investigated by a spectral survey in a certain spatial area. To gain insight into the variability of the electronic structure of the top BP layer, a spectroscopic survey was systematically performed within a 2.4 nm × 3.6 nm area, with its topographic image shown in the left panel of Figure 2c. This allowed comparing a series of positiondependent differential tunneling conductance spectra to effectively visualize the variation of the bandgap (Figure 2c). The states have energy values mostly above +0.8 V or below −1.2 V, with the Fermi energy EF at 0 V. However, at some locations, the gap edge at positive energy shifts away from EF, thus reducing the DOS at the positive gap edge and slightly enlarging the bandgap. The spectra measured along the “y” direction (the direction of a heavier effective mass) exhibits homogeneous bandgap distribution. In comparison, the spectra measured along “x” direction (the direction of a lighter effective mass) exhibits more electronic inhomogeneity. A LDOS mapping obtained by plotting the spatial variation of dI/dV signal at +1.15 V (the right panel of Figure 2c), shows that such electronic inhomogeneity is more prominent along the “x” direction (i.e., it is anisotropic). Such electronic inhomogeneity can be understood based on strain-induced gap modification. Theoretical work21 shows that, by applying of uniaxial stress perpendicular to the monolayer, the lattice along the “x” direction will be more sensitive to the stress. Deformation of the lattice will modify the band structure. The energy ordering of the conduction band valleys changes with strain in such a way that it is possible to switch from a nearly direct bandgap semiconductor to an indirect semiconductor, semimetal, and metal with the compression along only one direction. The bandgap inhomogeneity observed in C

dx.doi.org/10.1021/nl502892t | Nano Lett. XXXX, XXX, XXX−XXX

Nano Letters

Letter

calculations to study the local properties of the edges. As mentioned above, the experimental edge is very rough and probably a mixture of multiple edge directions (Figure 3a). Hence, we have considered PNRs with various edge directions, including the armchair edge (1,0) (i.e., along the x direction), the zigzag edge (0,1) (i.e., along the y direction), and the edge (1,3) (which corresponds to the general edge direction determined experimentally), as shown in Figure 4a. The DFT-calculated band structure of 2D phosphorene is also presented in Figure 4a, and its valence band maximum and conduction band minimum are denoted as 2D VBM and 2D CBM, respectively. For 1D PNRs, bands corresponding to 2D VBM and 2D CBM of phosphorene are also labeled as references. For the phosphorene nanoribbon with an armchair edge (Figure 4b), DFT-band structure indicates that the doubly degenerate edge bands are conduction bands, while the PNR with zigzag edge (Figure 4c) has edge bands crossing the Fermi level (hence metallic), consistent with previous theoretical reports.3,33 Consequently, edge DOS peaks for the armchairedged PNR are above the Fermi level, while the zigzag-edged PNR has multiple edge DOS peaks both above and below the Fermi level, as shown in Figure 4f. We have also studied PNRs with (1,3) and (3,1) edges, both of which exhibit metallic characters due to the presence of edge states. This result is in stark contrast with experimental dI/dV spectra shown in Figure 3 where the edge DOS peaks are all above but close to the 2D VBM peak, and the system is always semiconducting. Clearly, all edges considered above fail to explain the experimental observations. A close examination of the edges reveals that edge P atoms have only one or two covalent bonds instead of three covalent bonds as in the bulk (more details are provided in Figure S3 of the Supporting Information). The presence of dangling bonds on the edges could render the system unstable and could be directly responsible for the energy positions of edge bands. To see how the edge self-passivation affects the edge states, we further considered edge reconstructions. For a (2 × 1) supercell of zigzag-edged PNR, edge reconstruction can happen so that some of the P atoms close to the edge have five covalent bonds and all edge P atoms have three covalent bonds (Figure 4g and Figure S3a). While bulk BP is characterized by each P having exactly three neighbors, P can also assume a coordination number of 5. This hybridization effect leads to a significant stabilization of the edge by 1.62 eV (i.e., 0.25 eV/Å normalized by the edge length). Upon edge reconstruction, the edge bands are shifted below the Fermi level and above the 2D VBM, as shown from Figure 4c to Figure 4d. Similarly, for one unit cell of the PNR with a (1,3) edge in Figure 4e (also Figure S3b), the edge reconstruction lowers the total energy by 2.35 eV (i.e., 0.22 eV/Å normalized by the edge length), and the edge bands are also shifted below the Fermi level and above the 2D VBM. Consequently, as shown in Figure 4f, for the PNR with zigzag or (1,3) reconstructed edge, the edge DOS peak is below the Fermi level and above the 2D VBM peak, in a qualitative agreement with the experimental dI/dV spectra in Figure 3. More importantly, from the zigzag to (1,3) reconstructed edge in Figure 4f, the energy separation between the edge peak and 2D VBM peak is increased. This suggests that different reconstructed edges give rise to edge peaks at different energy positions. Since the experimental edge might be a mixture of multiple edges, multiple edge peaks in experimental dI/dV spectra (Figure 3) may be due to multiple reconstructed edges. For a (2 × 1) supercell of armchair-edged PNR, edge reconstruction is also found, but the reconstructed

our STS measurement is most likely induced by subtle deformation of the lattice along the x direction. The nonequivalent STM brightness between P1 and P2 atoms in Figure 1 can thus be explained by the strain-induced surface deformation, as discussed above. Note that the surface strain has been recently theoretically predicted to induce substantial and distinct shifts of frequencies of Raman modes of phosphorene,22 and experimental Raman studies on BP have begun to be reported in the literature.23,24 High-resolution Raman measurements could probe the surface strain of phosphorene in the near future. We now turn to the investigation of the properties of phosphorene edges. In graphene, edge-states localized at the zigzag graphene edge originate from a 2-fold degenerate flat band at the Fermi energy (EF) located at one-third of the Brillouin zone away from the zone center. Such states have been predicted first in theory25−28 and later observed in monatomic step edges of graphite by using scanning tunneling microscopy,29,30 and further observed on graphene nanoribbons.31,32 These edge states (which are extended along the edge direction) decay exponentially with the distance from the edge, with decay rates depending on their momentum. The edge states hold particular interests in graphene and its nanostructures, partially for its 1D nature. For 2D phosphorene, as shown previously (Figure 2a), its electronic structure has a clear bandgap of about 2.1 eV. There is no electronic state inside the gap. On the other hand, the edge states of graphene originate from the topology of the π electron networks with a zigzag edge. This is clearly not the case for the phosphorus’ sp3 electron system. In black phosphorus, each phosphorus atom assumes a coordination number of 3. To understand what happens when the phosphorene network is broken at an edge, we carried out experimental measurements by using STS corroborated with extensive DFT calculations. A step edge on the cleaved BP surface is shown in Figure 3a. The line profile (Figure 3b) measured across this step shows that the height of the step is about 5.3 Å, indicating it is a single layer edge. As shown in the topographic image, the shape of the edge is somewhat irregular, possibly indicating a mixture of multiple edge directions. From atomic resolved STM image, we can accurately pin down the overall direction of this edge, as discussed later. A representative dI/dV curve measured at the edge exhibits the appearance of multiple nonzero peaks inside of the bandgap. To highlight these DOS peaks, only the negative part of STS curve is shown in Figure 3c. Note that all of the edges we found in the samples (about five edges) show edge peaks in the STS spectra: some display two peaks, while others have more as the one shown in Figure 3c. In order to show the evolution of these states across the edge, a line spectrum survey is measured along the red dashed line, marked in Figure 3a. A series of dI/dV spectra are shown in Figure 3d, measured with at 4 Å intervals. These peaks appear gradually inside of the gap when approaching the edge and reach their maximum peak intensities at the edge to eventually vanish after passing the edge. A color plot of these dI/dV spectra with using lateral distance as the horizontal axis is shown in the right panel of Figure 3d. Note that the penetration length of these edge states is about 10 nm, much longer than the delocalization length of graphene edge states. To understand the nature of these edge states in anisotropic BP, we carefully considered a number of possible edge structures by constructing a series of 1D phosphorene nanoribbons (PNRs). We carried out extensive DFT D

dx.doi.org/10.1021/nl502892t | Nano Lett. XXXX, XXX, XXX−XXX

Nano Letters

Letter

Figure 4. DFT and GW calculations of edge states. (a) Structure of 2D phosphorene with armchair (1,0), zigzag (0,1), and (1,3) edge directions shown. The DFT band structure of the phosphorene is also shown. (b−e) DFT band structures and band-specific charge density distributions of 1D phosphorene nanoribbons (PNRs) with various edges. The Fermi level is set at 0 eV. For 1D PNRs, the bands corresponding to 2D phosphorene’s valence band maximum and conduction band minimum are denoted as 2D VBM and 2D CBM, respectively. Edge bands around the Fermi level are also labeled. (f) DFT densities of states (DOS) for systems from parts a to e. DOS peaks for the 2D VBM and CBM are marked by blue dashed lines, while edge peaks are by red dashed lines. (g) For the PNR with zigzag reconstructed edge, the evolution of GW local densities of states along the red dashed line from the interior to the edge. (h) For the zigzag-edge reconstructed PNR supported on a 2D phosphorene layer, plane-averaged electrostatic potential with a schematic diagram of band bending across the step edge. More details on the reconstructed edges are reported in Figure S3 of the Supporting Information.

CBM gap is around 2.01 eV, consistent with the experimental value 2.05 eV. Moving toward the zigzag reconstructed edge, the 2D VBM peak gradually decreases, while the edge peak appears gradually inside of the gap and reaches the maximum intensity at the edge, in agreement with the trend in experimental dI/dV spectra. As discussed above in reference to Figure 4f, different reconstructed edges contribute to edge peaks at different energy positions, and the (1,3) reconstructed edge is expected to result in another edge peak next to the edge peak by the zigzag reconstructed edge, thus explaining the multiple edge peaks in the experimental dI/dV spectra. Note that in the experimental dI/dV spectra, the energy positions of edge peaks are shifted with the distance from the edge (Figure 3). The shift of a DOS peak approaching the step edge has been reported previously and attributed to the electrostatic potential variation induced by the charge redistribution at the step edge.34−36 Although the evolution of GW local DOSs for the zigzag-edge reconstructed PNR in Figure 4g has confirmed the appearance of the edge peak inside of the gap, it fails to show any edge peak shift. This may be due to the fact that we only considered the free-standing PNR, while experimentally it is supported on multiple-layer BP.

system is energetically less stable by 1.28 eV (i.e., 0.14 eV/Å normalized by the edge length), and the edge bands are still above the Fermi level (Figure S3c). So the armchair edge reconstruction is ruled out to account for the experimental data. In addition to edge reconstruction, we also considered edge passivation by hydrogen atoms. For both armchair- and zigzagedged PNRs, our calculations find that hydrogen passivation causes edge bands to disappear in the 2D VBM-CBM gap (more details in Figure S4 in the Supporting Information).3 Although hydrogenated armchair- and zigzag-edged PNRs are still semiconductors, edge DOS peaks are no longer present, and hence edge passivation can be safely excluded as a potential explanation of the experimental observations. To quantitatively explain the experimental dI/dV spectra in Figure 3, we have also carried out single shot G0W0 calculations on the PNR with a zigzag reconstructed edge, as shown in Figure 4g. Along the red dashed line from the PNR’s interior to edge, the evolution of GW local densities of states is presented in the right panel of Figure 4g. The 2D VBM and CBM peaks are marked by blue dashed lines, while the edge peak is highlighted by a red dashed line. Starting from the interior (the bottom DOS line), no edge peak is detected, and the 2D VBME

dx.doi.org/10.1021/nl502892t | Nano Lett. XXXX, XXX, XXX−XXX

Nano Letters

Letter

Hence, to find out the variation of electrostatic potential across the step edge, we constructed a model of zigzag-edge reconstructed PNR on top of monolayer BP, more closely representing the experimental condition. The computed planeaveraged electrostatic potential across the step edge is presented in Figure 4h, clearly indicating a downward shift approaching the edge and then an upward shift away from the edge. Accordingly, as shown in the schematic diagram of band bending in Figure 4h, the edge band bends upward approaching the edge and then downward away from the edge, which might explain the observed upshift and then downshift of edge peaks in the experimental dI/dV spectra. In conclusion, an anisotropic two-dimensional material composed of black phosphorus was studied by low temperature STM/S and theoretical calculations. A ∼2 eV bandgap was identified by STS and rationalized by calculations within the self-consistent GW approximation. A strain effect leads to the observed nonequivalence of P1 and P2 atoms in the STM images and a bandgap inhomogeneity along “x” and “y” directions. The edge structure and electronic properties were determined from a quantitative agreement between experimental dI/dV spectra and a reconstructed model where local coordination numbers are modified as a self-passivation mechanism. Materials and Methods. STM Measurements. For the STM/S measurements, BP single crystals were cleaved in situ under UHV environment at room temperature, producing a mirror-like surface. A laboratory-built low-temperature scanning tunneling microscope (STM) was used for the imaging and spectroscopic measurements. The sample was cleaved at room temperature in ultrahigh vacuum (UHV) to expose a shining surface and then loaded into the STM head for investigation at about 80 K. We obtained topographic images in constantcurrent mode and the tunneling spectra dI/dV using a lock-in technique to measure differential conductance. A commercial Pt−Ir tip was prepared by gentle field emission at a clean metal sample. The bias voltage was applied on the sample during the STM observations. The WSxM software has been used to process and analyze STM data.37 Theoretical Calculations. Plane-wave DFT calculations were performed using the VASP package38 within the generalized gradient approximation (GGA) using the Perdew−Burke− Ernzerhof (PBE) exchange-correlation functional.39 The projector augmented wave (PAW) pseudopotentials were used with a cutoff energy set at 500 eV. For bulk BP, both atoms and cell volume were allowed to relax until the residual forces were below 0.01 eV/Å, with a 9 × 12 × 4 k-point sampling in the Monkhorst−Pack scheme.40 The optimized lattice parameters of bulk BP are a = 4.52 Å, b = 3.31 Å, and c = 11.06 Å. Single- and few-layer BP systems were then modeled by a periodic slab geometry using the optimized in-plane lattice constant of the bulk. A vacuum region of 18 Å in the z direction was used to avoid spurious interactions with replicas. For the 2D slab calculations, all atoms were relaxed until the residual forces were below 0.01 eV/Å, and we used a 9 × 12 × 1 k-point sampling. Starting from the DFT ground state, quasiparticle (QP) energies were then calculated within the self-consistent GW approximation. An iteration loop was run only for the calculation of G, while W was fixed to the initial DFT obtained W0 (This procedure is called GW0 in VASP).41,42 For both bulk and layered BP systems (including 1L to 6L), QP energies were iterated three times, and the energy cutoff for response function was chosen at 150 eV.

In addition, to study edge states of BP, various 1D phosphorene nanoribbons were constructed and relaxed until the residual forces were below 0.02 eV/Å with fine k-point samplings. For zigzag-edged phosphorene nanoribbons with edge reconstruction, single shot G0W0 calculations were also carried out to obtain the density of states at the GW level. The simulated STM images were computed using converged electronic densities within the Tersoff−Hamann approximation.43 Although the PBE functional can yield a satisfactory description of BP’s electronic properties,9 the optB86b−vdW functional can more accurately reproduce the experimental geometry of BP since it accounts for vdW interactions.44 According to our calculations, optB86b−vdW yields a height difference between P1 and P2 atoms proportional to the applied compressive strain, while PBE fails to capture the straininduced surface corrugation. Hence, for the STM simulation of the BP corrugated surface, the optB86b−vdW functional is used instead of PBE.



ASSOCIATED CONTENT

* Supporting Information S

DFT calculations of charge decay perpendicular to the BP, MoS2, and graphite surface, various edge reconstructions, and hydrogenation effects to phosphorene nanoribbons. This material is available free of charge via the Internet at http:// pubs.acs.org.



AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected] (M.P.). *E-mail: [email protected] (V.M.). Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS Part of this research was conducted at the Center for Nanophase Materials Sciences (J.W., B.G.S., M.P.), which is sponsored at Oak Ridge National Laboratory by the Scientific User Facilities Division, Office of Basic Energy Sciences, U.S. Department of Energy. The work at Rensselaer Polytechnic Institute (RPI) was supported by New York State under NYSTAR program C080117 and the Office of Naval Research. The computations were performed using the resources of the Center for Computational Innovation at RPI.



REFERENCES

(1) Reich, E. S. Nature 2014, 506 (7486), 19. (2) Li, L.; Yu, Y.; Ye, G. J.; Ge, Q.; Ou, X.; Wu, H.; Feng, D.; Chen, X. H.; Zhang, Y. Nat. Nanotechnol. 2014, 9, 372−377. (3) Liu, Y.; Xu, F.; Zhang, Z.; Penev, E. S.; Yakobson, B. I. Nano Lett. 2014, 10.1021/nl5021393. (4) Liu, H.; Neal, A. T.; Zhu, Z.; Luo, Z.; Xu, X.; Tománek, D.; Ye, P. D. ACS Nano 2014, 8, 4033−4041. (5) Zhang, C. D.; Lian, J. C.; Yi, W.; Jiang, Y. H.; Liu, L. W.; Hu, H.; Xiao, W. D.; Du, S. X.; Sun, L. L.; Gao, H. J. J. Phys. Chem. C 2009, 113, 18823−18826. (6) Novoselov, K.; Geim, A.; Morozov, S.; Jiang, D.; Zhang, Y.; Dubonos, S.; Grigorieva, I.; Firsov, A. Science 2004, 306, 666−669. (7) Berger, C.; Song, Z.; Li, T.; Li, X.; Ogbazghi, A. Y.; Feng, R.; Dai, Z.; Marchenkov, A. N.; Conrad, E. H.; First, P. N.; A. de Heer, W. J. Phys. Chem. B 2004, 108, 19912−19916. (8) Takao, Y.; Morita, A. Physica B, C 1981, 105, 93−98. (9) Tran, V.; Soklaski, R.; Liang, Y. F.; Yang, L. Phys. Rev. B 2014, 89 (23), 235319. F

dx.doi.org/10.1021/nl502892t | Nano Lett. XXXX, XXX, XXX−XXX

Nano Letters

Letter

(10) Radisavljevic, B.; Radenovic, A.; Brivio, J.; Giacometti, V.; Kis, A. Nat. Nanotechnol. 2011, 6 (3), 147−150. (11) Mak, K. F.; Lee, C.; Hone, J.; Shan, J.; Heinz, T. F. Phys. Rev. Lett. 2010, 105, 136805. (12) Splendiani, A.; Sun, L.; Zhang, Y.; Li, T.; Kim, J.; Chim, C.-Y.; Galli, G.; Wang, F. Nano Lett. 2010, 10, 1271−1275. (13) Kuc, A.; Zibouche, N.; Heine, T. Phys. Rev. B 2011, 83, 245213. (14) Liang, L.; Meunier, V. Nanoscale 2014, 6, 5394−5401. (15) Huang, S.; Ling, X.; Liang, L.; Kong, J.; Terrones, H.; Meunier, V.; Dresselhaus, M. S. Nano Lett. 2014, 14, 5500−5508. (16) Keyes, R. W. Phys. Rev. 1953, 92, 580−584. (17) Warschauer, D. J. Appl. Phys. 1963, 34, 1853−1860. (18) Maruyama, Y.; Suzuki, S.; Kobayashi, K.; Tanuma, S. Physica B, C 1981, 105, 99−102. (19) Akahama, Y.; Endo, S.; Narita, S. J. Phys. Soc. Jpn. 1983, 52, 2148−2155. (20) Asahina, H.; Shindo, K.; Morita, A. J. Phys. Soc. Jpn. 1982, 51, 1193−1199. (21) Rodin, A. S.; Carvalho, A.; Neto, A. H. C. Phys. Rev. Lett. 2014, 112 (17), 176801. (22) Fei, R.; Yang, L. Appl. Phys. Lett. 2014, 105, 083120. (23) Zhang, S.; Yang, J.; Xu, R.; Wang, F.; Li, W.; Ghufran, M.; Zhang, Y.-W.; Yu, Z.; Zhang, G.; Qin, Q.; Lu, Y. ACS Nano 2014, 8, 9590−9596. (24) Xia, F.; Wang, H.; Jia, Y. arXiv preprint arXiv:1402.0270, 2014. (25) Fujita, M.; Wakabayashi, K.; Nakada, K.; Kusakabe, K. J. Phys. Soc. Jpn. 1996, 65, 1920−1923. (26) Nakada, K.; Fujita, M.; Dresselhaus, G.; Dresselhaus, M. S. Phys. Rev. B 1996, 54, 17954−17961. (27) Wakabayashi, K.; Fujita, M.; Ajiki, H.; Sigrist, M. Phys. Rev. B 1999, 59, 8271−8282. (28) Miyamoto, Y.; Nakada, K.; Fujita, M. Phys. Rev. B 1999, 59, 9858−9861. (29) Kobayashi, Y.; Fukui, K.-I.; Enoki, T.; Kusakabe, K.; Kaburagi, Y. Phys. Rev. B 2005, 71, 193406. (30) Niimi, Y.; Matsui, T.; Kambara, H.; Tagami, K.; Tsukada, M.; Fukuyama, H. Phys. Rev. B 2006, 73, 085421. (31) Tao, C.; Jiao, L.; Yazyev, O. V.; Chen, Y.-C.; Feng, J.; Zhang, X.; Capaz, R. B.; Tour, J. M.; Zettl, A.; Louie, S. G.; Dai, H.; Crommie, M. F. Nat. Phys. 2011, 7, 616−620. (32) Pan, M.; Girao, E. C.; Jia, X.; Bhaviripudi, S.; Li, Q.; Kong, J.; Meunier, V.; Dresselhaus, M. S. Nano Lett. 2012, 12, 1928−1933. (33) Guo, H.; Lu, N.; Dai, J.; Wu, X.; Zeng, X. C. J. Phys. Chem. C 2014, 118 (25), 4051−14059. (34) Ono, M.; Nishigata, Y.; Nishio, T.; Eguchi, T.; Hasegawa, Y. Phys. Rev. Lett. 2006, 96, 016801. (35) Hasegawa, Y.; Ono, M.; Nishigata, Y.; Nishio, T.; Eguchi, T. Real-Space Observation of Screened Potential and Friedel Oscillation by Scanning Tunneling Spectroscopy. J. Phys.: Conf. Series 2007, 61, 399. (36) Nakajima, Y.; Takeda, S.; Nagao, T.; Hasegawa, S.; Tong, X. Phys. Rev. B 1997, 56, 6782. (37) Horcas, I.; Fernandez, R.; Gomez-Rodriguez, J.; Colchero, J.; Gómez-Herrero, J.; Baro, A. Rev. Sci. Instrum. 2007, 78, 013705. (38) Kresse, G.; Furthmüller, J. Comput. Mater. Sci. 1996, 6, 15−50. (39) Perdew, J.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1996, 77, 3865−3868. (40) Monkhorst, H. J.; Pack, J. D. Phys. Rev. B 1976, 13, 5188−5192. (41) Liang, L.; Meunier, V. Phys. Rev. B 2012, 86, 195404. (42) Liang, L.; Girão, E. C.; Meunier, V. Phys. Rev. B 2013, 88, 035420. (43) Tersoff, J.; Hamann, D. Phys. Rev. B 1985, 31, 805. (44) Qiao, J.; Kong, X.; Hu, Z.-X.; Yang, F.; Ji, W. Nat. Commun. 2014, 5, 4475.

G

dx.doi.org/10.1021/nl502892t | Nano Lett. XXXX, XXX, XXX−XXX