High-mobility transport anisotropy and linear dichroism in few-layer black phosphorus

Abstract

Two-dimensional crystals are emerging materials for nanoelectronics. Development of the field requires candidate systems with both a high carrier mobility and, in contrast to graphene, a sufficiently large electronic bandgap. Here we present a detailed theoretical investigation of the atomic and electronic structure of few-layer black phosphorus (BP) to predict its electrical and optical properties. This system has a direct bandgap, tunable from 1.51 eV for a monolayer to 0.59 eV for a five-layer sample. We predict that the mobilities are hole-dominated, rather high and highly anisotropic. The monolayer is exceptional in having an extremely high hole mobility (of order 10,000 cm² V^-1 s^-1 ) and anomalous elastic properties which reverse the anisotropy. Light absorption spectra indicate linear dichroism between perpendicular in-plane directions, which allows optical determination of the crystalline orientation and optical activation of the anisotropic transport properties. These results make few-layer BP a promising candidate for future electronics.

ARTICLE

Received 27 Feb 2014 | Accepted 20 Jun 2014 | Published 21 Jul 2014

Two-dimensional crystals are emerging materials for nanoelectronics. Development of the eld requires candidate systems with both a high carrier mobility and, in contrast to graphene, a sufciently large electronic bandgap. Here we present a detailed theoretical investigation of the atomic and electronic structure of few-layer black phosphorus (BP) to predict its electrical and optical properties. This system has a direct bandgap, tunable from1.51 eV for a monolayer to 0.59 eV for a ve-layer sample. We predict that the mobilities are hole-dominated, rather high and highly anisotropic. The monolayer is exceptional in having an extremely high hole mobility (of order 10,000 cm2 V 1 s 1) and anomalous elastic properties which reverse the anisotropy. Light absorption spectra indicate linear dichroism between perpendicular in-plane directions, which allows optical determination of the crystalline orientation and optical activation of the anisotropic transport properties. These results make few-layer BP a promising candidate for future electronics.

DOI: 10.1038/ncomms5475 OPEN

High-mobility transport anisotropy and linear dichroism in few-layer black phosphorus

Jingsi Qiao1,2,*, Xianghua Kong1,2,*, Zhi-Xin Hu1,2, Feng Yang1,2,3 & Wei Ji1,2

1 Department of Physics, Renmin University of China, Beijing 100872, China. 2 Beijing Key Laboratory of Optoelectronic Functional Materials & Micro-nano Devices, Renmin University of China, Beijing 100872, China. 3 College of Physics and Electronic Engineering, Institute of Solid State Physics, Sichuan Normal University, Chengdu 610068, China. * These authors contributed equally to this work. Correspondence and requests for materials should be addressed to W.J. (email: mailto:[email protected]

Web End [email protected] ).

NATURE COMMUNICATIONS | 5:4475 | DOI: 10.1038/ncomms5475 | http://www.nature.com/naturecommunications

Web End =www.nature.com/naturecommunications 1

& 2014 Macmillan Publishers Limited. All rights reserved.

ARTICLE NATURE COMMUNICATIONS | DOI: 10.1038/ncomms5475

The discovery of graphene laid the foundations for many new areas of research. One of the most important foundation is the intensive investigation of

two-dimensional (2D) atomic-layer systems, including graphene itself15, transition metal dichalcogenides (TMDCs)58, silicene9,10 and germanane10,11, as candidate materials for future electronics applications1219. A high-performance device such as a eld-effect transistor (FET) requires a moderate electronic band gap, a reasonably high carrier mobility of the channel material and excellent electrode-channel contacts4,5,8,1319. Graphene offers extremely high mobilities, due to its very low carrier effective mass, and thus is considered to be a promising candidate for high-speed FET devices, but its intrinsic dispersion is gapless24,1214. Despite extensive efforts following a wide variety of approaches to the problem of opening a gap in different graphene nanostructures, all devices to date have a relatively large off current and thus a low onoff ratio13. For this reason the emergence of monolayer TMDs has attracted substantial research interest, and in particular MoS2 has recently been used to fabricate a FET15. Unlike graphene, mono-layer MoS2 is a direct-bandgap semiconductor with a carrier mobility of B200 cm2 V 1 s 1 (ref. 15), improvable up to 500 cm2 V 1 s 1 (ref. 8), a value which is not unreasonable for applications but remains orders of magnitude lower than that of graphene13,14. Thus the entire community working on atomic-layer transport continues to search for a 2D material which is semiconducting, preferably with a direct gap, has high carrier mobility and has the potential to form excellent contacts with known electrode materials.

Here we present the theoretical discovery of a new category of layered direct-bandgap semiconductor, few-layer black phosphorus (BP), with high mobility, high in-plane anisotropy and linear dichroism. BP is an allotrope of phosphorus, whose crystal structure is a strongly folded honeycomb sheet with troughs running along the y-axis (b direction). We demonstrate using density functional theory (DFT) calculations that few-layer BP systems, from the monolayer up to ve-layer structures, are thermally stable with an interlayer interaction energy of

0.46 eV and a bandgap that falls exponentially with the thickness. Carrier mobilities at room temperature are high, and furthermore exhibit strongly anisotropic behaviour. The BP monolayer has an exceptionally high hole mobility, which we predict to lie in the range from 10,000 to 26,000 cm2 V 1 s 1, and which we explain by visualizing the wavefunctions. We nd an explicit linear dichroism in the computed absorption spectra, whereby the positions of the lowest-energy absorption peaks for the two in-plane directions differ strongly. These results

demonstrate that few-layer BP is a new category of 2D semiconductor with high potential for novel applications in nanoelectronics and optoelectronics.

ResultsGeometric and electronic properties of bulk BP. Because no experimental data are available for the atomic structure and electronic properties of few-layer BP, we begin by considering bulk BP to rene the accuracy of our theoretical predictions. We have performed DFT calculations using a number of different functionals to gauge which provide the best t to experiment; a full list may be found in Supplementary Table 1 and we summarize our ndings here. If one considers only the lattice geometry then the PBE-G06 (refs 20,21) and optB86b-vdW22 methods produce the best results. However, neither local density approximation-modied BeckeJohnson (LDA-mBJ)23,24 nor HSE06 (refs 25,26) calculations based on these best geometries can provide a satisfactory value for the bulk bandgap, as shown in Supplementary Fig. 1. To optimize both sets of properties simultaneously, we found that the optB88-vdW functional22,27, combined with either the LDA-mBJ or HSE06 method to predict the electronic bandstructure, gave the best ts. Details of choosing the best-t functional are available in Supplementary Methods. Figure 1a shows the fully relaxed atomic structure of bulk BP, where the equilibrium lattice constants obtained from optB88-vdW are only 12% larger than the experiment2830. The accompanying Brillouin zone (BZ) and electronic bandstructures, labelled mBJ (optB88-vdW) and HSE06 (optB88-vdW) are shown in Fig. 1b,c. Both combinations of methods predict that bulk BP is a semiconductor with a direct bandgap at the Z point of 0.31 eV (mBJ) or 0.36 eV (HSE06), both values being fully consistent with the experimental value of 0.310.35 eV (refs 3134).

At the bandgap wavevector, which is the Z point, one valence band (VB) and one conduction band (CB) disperse rather strongly along the ZQ and ZG directions (Fig. 1c), indicating very small effective masses. A t of these bands using the nearly-free electron model gives effective carrier masses along the ZQ direction which are rather small and similar, namely 0.12 m0 for electrons and 0.11 m0 for holes, whereas these along ZG are slightly larger, taking the respective values 0.15 m0 and 0.30 m0 (all masses taken from HSE06). Considerably larger values are found along ZT0A0, where the carrier effective masses are 1.15 m0 and 0.71 m0, respectively. All of these results are fully consistent with the experimental values for bulk BP35. The corresponding mBJ results are within 0.04 m0 of the HSE06 values along ZQ, ZG and ZT0A0 except for the electron

a b c

2.0

Q

T'

1.0

0.0

1.0

Z

A'

Y

N

c

E E VBM (eV)

1.15

0.71

0.36 eV

0.31 eV

R1

R2

X

T A

0.12

0.11 0.30

0.15

[afii9835]2

G

[afii9835]1

z

y

a b

x

2.0

Z T' A(') T G X Q Z G Y N X

Figure 1 | Lattice and electronic structures of bulk black phosphorus. (a) Crystal structure of bulk BP marked with coordinate axes (x, y, z), lattice vectors (a, b, c) and structural parameters (R1, R2, y1 and y2). (b) Brillouin zone path of BP primitive cell. (c) Electronic bandstructures for bulk BP calculated with the HSE06 functional (red solid line) and the mBJ potential (blue dashed line), together with tted effective masses along the ZT0A0, ZQand ZG directions. At the right of the image, a zoomed-in plot shows the direct bandgap at Z. EVBM is the energy of valence-band maximum.

2 NATURE COMMUNICATIONS | 5:4475 | DOI: 10.1038/ncomms5475 | http://www.nature.com/naturecommunications

Web End =www.nature.com/naturecommunications

& 2014 Macmillan Publishers Limited. All rights reserved.

NATURE COMMUNICATIONS | DOI: 10.1038/ncomms5475 ARTICLE

effective mass along ZT0A0, which is larger than the HSE06 value by 0.13 m0. Values for the direction ZQ are similar to other high-mobility semiconductors, such as AlGaAs (0.1 m0)36,

and four times smaller than those of MoS2 (ref. 37).

Geometric and electronic properties of few-layer BP. The similarity between the bandstructures of bulk and monolayer BP38, and the low effective masses of bulk BP35, suggest that few-layer BP is likely to be a high-mobility, direct-bandgap 2D semiconductor. We have thus performed direct calculations of the geometric and electronic properties for few-layer BP systems. Table 1 summarizes the changes of geometrical properties as a function of the layer number from 1 to 5. The lattice parameter a increases by 0.11 on passing from bulk to monolayer BP, whereas b grows only by 0.02 . There is an abrupt reduction of a between the monolayer and the bilayer, which we believe to be the primary consequence of interlayer interactions (0.46 eV) in the bilayer. The signicant stretching of a in few-layer BP systems is caused almost exclusively by an increase of the bond angle y1 but not of the bond lengths.

Results from our bandstructure calculations for the ve few-layer BP systems are shown in Fig. 2. In the monolayer, the original Z point of the bulk BZ folds back to the Gamma (G) point, so that the original ZQ and ZT0A0 directions of the bulk

BZ project onto the GX and GY directions of the monolayer, which correspond respectively to the a and b directions in real space (Fig. 2a,b). Results obtained from the mBJ method are quantitatively the same as those obtained from HSE06. Supplementary Figure 2 and Supplementary Table 2 show that the predicted bandgap is rather less sensitive than the unit-cell dimensions to the functional used for optimizing the atomic structures, with similar values emerging for all few-layer BP systems from the best-performing functionals. However, other physical quantities, notably the location of band maximum (minimum) and the carrier effective mass (and hence the mobility) are in fact very sensitive to the choice of the functional, and hence we have focused primarily on our optB88-vdW results for the discussion of these. Monolayer BP is indeed a direct-bandgap semiconductor, as shown in Fig. 2c, with the gap value of 1.51 eV obtained at the G point. Phosphorus is a heavier element than carbon and therefore should have stronger spin-orbit coupling (SOC) in its 2D forms than graphene. We have considered full SOC effects in calculating the electronic bandstructure of bulk BP. Supplementary Figure 3 shows that inclusion of SOC terms has no appreciable effect on the primary features of the bandstructure, indicating that phosphorus is not sufciently heavy to cause any qualitative changes. Only a very small separation of the formerly fourfold degenerate bands, into two band pairs, can be found from X along N (or M for few-layer BP) to Y, which reaches a maximum of 23 meV around the N (M) point. Along GX and GY, these bands are already

separated by bandstructure effects and SOC causes no further splitting. That the spin-orbit splitting energy is so small is a consequence of the small effective nuclear charge (Zeff) of the P

atom and the weak variation of the charge gradient in an elemental system such as BP.

When two monolayers are combined to form a bilayer (Fig. 2d,e), the gap is reduced to 1.02 eV and two additional bands emerge around the gap at the G point. Together with the original VB and CB, we denote these bands as VB1, VB2, CB1 and CB2, as shown in Fig. 2f. In real space, the states from these four bands with k close to G are extended throughout the bilayer, as shown in Fig. 2g. States VB1 and VB2 differ in the interlayer region although they share the same origin in terms of the atomic orbitals. A clear bonding-like feature is visible in the interlayer region (marked by red rectangles) for VB1, which lies lower in energy, whereas VB2 shows an anti-bonding feature. Similar bonding and anti-bonding features are also found in CB1 and CB2, where they are observed not between the layers but across the troughs. These features indicate that wavefunction overlap, rather than van der Waals effects, plays the primary role in mediating the interlayer interaction, which explains the abrupt reduction of the lattice constant a from monolayer to bilayer. This interlayer interaction introduces band dispersions of VBs and CBs in the direction normal to the layers, which leads to the reduction of the bandgap (by 0.5 eV) from monolayer to bilayer.

The thicker the few-layer BP system, the stronger is the interlayer interaction and thus the larger is the dispersion of VBs (CBs), resulting in a smaller bandgap. Thus the gap falls continuously on adding more layers to the bilayer, reaching0.59 eV in ve-layer BP, as illustrated in Fig. 2h. We tted these values by an exponential decay relation and found that the corresponding bulk gap obtained by extrapolation is 0.53 eV, a value 0.17 eV larger than the one we calculated for bulk BP. We ascribe this difference to the elongated lattice parameter a in few-layer BP and results calculated with a series of constrained a values conrm our expectation; these results imply that the bandgap is very sensitive to the lateral strain along the a (x) direction and could be modulated by varying interatomic separations and angles. During the preparation of this manuscript, we became aware that a strain-induced bandgap modulation has recently been demonstrated theoretically39.

The small effective masses remain in all of the few-layer BP systems (Table 2). For the GX direction in the monolayer, equivalent to ZQ in the bulk, the carrier effective masses are 0.15 m0 (hole) and 0.17 m0 (electron), only 0.040.05 m0 larger than those in bulk BP. It is remarkable that in the GY direction, equivalent to ZT0A0 in the bulk, the VB appears to be nearly at close to the G point, with an effective mass of 6.35 m0 (see

Fig. 2c), nine times its bulk value of 0.71 m0. By contrast, the effective mass for the CB is 1.12 m0, very close to its bulk value of1.15 m0. Unlike the bandgap, only the hole effective mass along

Table 1 | Structural information for few-layer BP.

NL a () b () Dc () R1 () R2 () h1/h10 () h2/h20 ()1 4.58 3.32 3.20 2.28 2.24 103.51 96.002 4.52 3.33 3.20 2.28 2.24 102.96 96.21/95.92 3 4.51 3.33 3.20 2.28 2.24 102.81/102.74 96.30/95.99 4 4.50 3.34 3.20/3.21 2.28 2.24 102.76/102.67 96.34/96.01 5 4.49 3.34 3.20/3.21 2.28 2.24 102.71/102.63 96.37/96.05 Bulk 4.47 3.34 3.20 2.28 2.25 102.42 96.16

BP, black phosphorus.

Lattice constants a, b and Dc (the interlayer spacing between two adjacent BP layers) and structural parameters R1, R2, y1 and y2 of few-layer and bulk BP calculated using optB88-vdW. There are slight differences in y1 and y2 between the outermost and inner layers.

NATURE COMMUNICATIONS | 5:4475 | DOI: 10.1038/ncomms5475 | http://www.nature.com/naturecommunications

Web End =www.nature.com/naturecommunications 3

& 2014 Macmillan Publishers Limited. All rights reserved.

ARTICLE NATURE COMMUNICATIONS | DOI: 10.1038/ncomms5475

c

Monolayer

a

d e

f

h

b

b

Band gap (eV)E E VBM(eV)E E VBM(eV)

y

a

x

a

3.0

2.0

1.0 CB2

CB1

VB2

VB1

0.0

1.0

2.0

3.0 G

1.5

1.2

0.9 0.91

z

z

z

z

x

y

g

X M Y G

CB1

b

Y M

G X

VB1 VB2 CB2

3.0

2.0

1.0

0.0

1.0

2.0

3.0 G X M Y G

Bilayer

1.53

HSE06 (vdWDF) mBJ (vdWDF)

PBE

HSE06 (PBEG06)

1.51

1.41

x

0.49 0.36

0.53

0.45

0.31

0.6

0.27

0.3

0.0

0.24 0.15

1 2 3 4 5 Inf Bulk

y

Number of layers

Figure 2 | Electronic structures of few-layer BP. (a,b) Top view of the atomic structure of the monolayer and the associated Brillouin zone. (d,e) Side views of the atomic structure of the bilayer. (c,f) Bandstructures of monolayer and bilayer BP calculated with the HSE06 functional (red solid lines) and the mBJ potential (blue dashed lines), respectively. Two valence (VB1 and VB2) and two conduction states (CB1 and CB2) are marked in panel f. (g) Spatial structure of wavefunctions for the four marked states illustrated in the xz and yz planes using an isosurface of 0.0025 e 3. (h) Evolution of the direct bandgaps as a function of the sample thickness. Functionals used for structural optimization are shown in parentheses. Bandgap values are markedfor the monolayer system, for the extrapolation of our results and for real bulk BP.

Table 2 | Predicted carrier mobility.

Carrier type NL mx*/m0 my*/m0 E1x E1y Cx_2D Cy_2D lx_2D ly_2D

GX GY (eV) (J m 2) (103 cm2 V 1 s 1)e 1 0.17 1.12 2.720.02 7.110.02 28.94 101.60 1.101.14 B0.08

2 0.18 1.13 5.020.02 7.350.16 57.48 194.62 B0.60 0.140.16

3 0.16 1.15 5.850.09 7.630.18 85.86 287.20 0.760.80 0.200.22

4 0.16 1.16 5.920.18 7.580.13 114.66 379.58 0.961.08 0.260.30

5 0.15 1.18 5.790.22 7.350.26 146.58 479.82 1.361.58 0.360.40

h 1 0.15 6.35 2.500.06 0.150.03 28.94 101.60 0.640.70 1026

2 0.15 1.81 2.450.05 1.630.16 57.48 194.62 2.62.8 1.32.2

3 0.15 1.12 2.490.12 2.240.18 85.86 287.20 4.45.2 2.23.2

4 0.14 0.97 3.160.12 2.790.13 114.66 379.58 4.45.2 2.63.2

5 0.14 0.89 3.400.25 2.970.18 146.58 479.82 4.86.4 3.04.6

2D, two-dimensional.

Carrier types e and h denote electron and hole, respectively. N represents the number of layers, m * and m * are carrier effective masses for directions x and y, respectively, E (E ) and C (C ) are the deformation potential and 2D elastic modulus for the x (y) direction. Mobilities m and m were calculated using equation (1) with the temperature T set to 300 K.

4 NATURE COMMUNICATIONS | 5:4475 | DOI: 10.1038/ncomms5475 | http://www.nature.com/naturecommunications

Web End =www.nature.com/naturecommunications

& 2014 Macmillan Publishers Limited. All rights reserved.

NATURE COMMUNICATIONS | DOI: 10.1038/ncomms5475 ARTICLE

GY shows a strongly layer-dependent evolution, in that it decreases from 6.35 m0 for the monolayer to 1.81 m0 for the bilayer and eventually to 0.89 m0 in the ve-layer system, close to the bulk value of 0.71 m0.

Carrier mobility. The electronic properties of few-layer BP are governed to a large extent by the carrier mobilities, which are in turn strongly inuenced (but not solely determined) by their effective masses. Thus we also provide theoretical predictions for the carrier mobilities, in both x and y directions, in few-layer BP systems by applying a standard 2D model as discussed under Methods below. We apply a phonon-limited scattering model, in which the primary mechanism limiting carrier mobility is scattering due to phonons36,4042. Two properties of the few-layer BP lattice, namely the deformation potential E1 and the elastic modulus C2D in the propagation direction of the longitudinal acoustic wave, are then the most relevant factors determining the mobility43,44. These quantities were calculated using the optB88-vdW functional and the data are presented in Table 2. If few-layer BP is cut into ribbons, a one-dimensional model should be adopted for predicting the carrier mobilities. We present the results and discussion for a 10-nm-width BP ribbon in the Supplementary Discussion and Supplementary Table 3.

Our predicted mobilities for the ve few-layer BP systems (Table 2) are in general moderately large (hundreds to thousands of cm2 V 1 s 1), moderately anisotropic and asymmetric between electrons and holes, with the holes being more mobile in both directions. Concerning directional anisotropy, electron and hole mobilities for the x direction are nearly twice the values along y for holes and four times as large for electrons. However, the monolayer is exceptional in every way. The electron mobility along x is nearly 14, rather than four, times the value along y, that is, 1,1001,140 cm2 V 1 s 1 compared with B80 cm2 V 1 s 1, which is consistent with the trend of the other few-layer systems.

To our great surprise, however, the hole mobility along x is 1638 times smaller than that along y, that is, 640700 cm2 V 1 s 1 versus 10,00026,000 cm2 V 1 s 1, making y the direction of higher hole conductivity. This extraordinarily large value for monolayer BP is a consequence of the extremely small deformation potential, E1y 0.150.03 eV, and occurs despite

the fact that the carriers are very heavy (6.35 m0). This monolayer exception should be valuable in applications such as the separation of electrons and holes.

The value of E1y for holes in a monolayer is very striking in that it is an order of magnitude smaller than typical values of E1,

which are 5.0 eV for graphene44, 3.9 eV for MoS2 (ref. 37), and3.7 eV for h-BN40 and 7 eV for AlGaAs quantum wells36. In few-layer BP systems, E1y for holes has a conventional value2.970.18 eV in a ve-layer structure, decreasing smoothly to1.630.16 eV in a bilayer. These results can be explained rather well by the form of the VB wavefunctions shown in Fig. 2g, where interlayer and stacking-induced intra-layer overlap increases with layer number. These effects are completely absent in monolayer BP, making the VB wavefunction quite isolated along y, as it is also for the VB2 (anti-bonding) state of a bilayer highlighted by the green rectangle in Fig. 2g. Thus small structural deformations due to longitudinal phonon oscillations along y have almost no effect on this wavefunction and cause little change in its energy, resulting in the very small deformation potential. The situation is quite different for the CB (electron) wavefunction, which has E1x 2.720.02 eV in the monolayer and 5.020.02 eV

in the bilayer. The in-plane bonding and anti-bonding features of the CB wavefunctions in few-layer BP systems cause a substantial enhancement of the energy-level perturbations due to structural deformations along y. The monolayer exception aside, all E1 values for CB states are rather large (5.07.7 eV), for both x and y directions, and are almost independent of the number of layers, a result which can be explained by the strong CB wavefunction overlap, especially along y (Fig. 2g). For the VB states, the potentials are somewhat smaller (23 eV), precisely because their wavefunctions have less overlap between different P atoms.

Optical absorption spectra and linear dichroism. We have also predicted the optical absorption spectra of few-layer BP systems by computing the dielectric function. Two absorption spectra are shown in Fig. 3a,b for light incident along the z direction and linearly polarized in the x and y directions, respectively; while the results for incident light polarized in the z direction may be found in Supplementary Fig. 4. These results demonstrate a strong linear dichroism: for a dielectric polarization in the x direction, the band edge of the rst absorption peak is found at the bandgap and thus falls rapidly with the thickness of the sample (Fig. 3a). By contrast, with y-polarized light this peak is found at 3.14 eV in the monolayer and its position falls only slightly with thickness, remaining at 2.76 eV in bulk BP (Fig. 3b). From the symmetries of

a b c

0.0

1-layer 2.00.0

x y

2.0

1.55

1-layer 3.14

Polarized light

Detector

yy

2-layer

1.05 2-layer 3.00

1.50.0

5.00.0

Absorbance (%)

2.93

2.00.0

3-layer 5.00.0

3-layer

4-layer 4-layer 2.92

5-layer 0.60

2.00.0

0.78

0.67 1.5

0.0

5.00.0

BP

D

S

xx

0.1

10.00.0

5-layer 2.90

Bulk Bulk

0.46

1.2

2.00.0

Oxides

0.0

0.4 0.8

1.6 2.0 2.0 2.4

2.76

2.8 3.2 3.6

Back gate

Energy of incident light (eV)

Energy of incident light (eV)

Figure 3 | Optical absorption spectra. (a,b) Optical absorption spectra of few-layer BP for light incident in the c (z) direction and polarized along the a (x) and b (y) directions, respectively. Black dashed lines show an approximate linear t used to estimate the band edges for the rst absorption peak, which are highly anisotropic between x and y. The band edge drops rapidly with sample thickness for x (from 1.55 eV for the monolayer (red) to 0.60 eV for ve layers (orange)) but only slightly for y (from 3.14 eV to 2.90 eV). (c) Schematic illustration of a proposed experimental geometry to determine the orientation of few-layer BP structures using optical absorption spectroscopy, and thus to utilize the anisotropic electronic properties of BP. The light is linearly polarized in a chosen orientation and near-normally incident on the sample. The sample should be rotated to identify the a or b direction by monitoring the absorption signal, after which the source and drain electrodes, denoted as S and D respectively, may be deposited to fabricate an FET device.

NATURE COMMUNICATIONS | 5:4475 | DOI: 10.1038/ncomms5475 | http://www.nature.com/naturecommunications

Web End =www.nature.com/naturecommunications 5

& 2014 Macmillan Publishers Limited. All rights reserved.

ARTICLE NATURE COMMUNICATIONS | DOI: 10.1038/ncomms5475

the wavefunctions in Fig. 2g, clearly the (odd) dipole operator connects the VB and CB states for x-polarization, allowing the direct-bandgap process, but this is symmetry-forbidden for ypolarization and the transition occurs between VB and CB states elsewhere in the BZ.

DiscussionHigh carrier mobilities would be expected in BP systems due to the small effective masses of 0.10.2 m0. However, our calculations reveal that the key issue is not the effective mass but the very small elastic modulus C2D: all the BP systems we have studied are very soft in the x direction, due to the fact that bond angles are changed easily (with little effect on bond lengths). When compared to other monolayer systems, the modulus of BP for the x direction is 29.0 J m 2, a value respectively 2, 10 and 11 times smaller than MoS2 (ref. 37), h-BN40 and graphene44. These results indicate that much higher mobilities could be most likely achieved if few-layer BP systems could be modied to become stiffer, that is, to have a larger elastic modulus for the x direction, for example by external stress, doping or substrate connement. However, possible changes of the deformation potential as a result of modifying the elastic moduli are not yet clear and require further investigation. The modulus for the y direction is signicantly larger than that for x, but the higher carrier effective masses, result in carrier mobilities being smaller than those along x. Both the bandgap and the mobility depend strongly on the structural properties and this offers a number of possibilities for further control over the electrical and mechanical properties of few-layer BP.

Further, the linear dichroism we nd allows an electronic determination of the sample orientation. We propose an experimental setup to deduce the orientation of few-layer BP systems using optical spectroscopy, as illustrated in Fig. 3c. If the incident light is linearly polarized in a chosen orientation and is near-normally incident, the absorption spectrum should vary when the sample is rotated and such that the a and b directions can be identied by monitoring the absorption signal. Once the sample orientation is determined, it is far easier to fabricate electrodes or gates utilizing the highest-mobility direction of a few-layer BP sample in an FET-type device.

In summary, we have shown theoretically that few-layer BP is a novel category of 2D semiconductor offering a direct bandgap, high carrier mobility and high transport anisotropy, all of which are tunable by controlling the layer thickness. The bandgap decreases from 1.51 eV for a monolayer to 0.59 eV for ve-layer BP, a set of values tting well into the gap between graphene nanoribbons and TMDCs. The anisotropy of electric and optical properties is another unique feature distinguishing BP from other 2D materials, such as silicene9,10, h-BN5,45 and germanane10,11. Higher conductivity is generally found in the direction oriented perpendicular to the troughs (the x direction) and holes are usually more mobile than electrons. For holes, the mobility in the x direction increases from 600 cm2 V 1 s 1 for a monolayer to over 4,000 cm2 V 1 s 1 for ve-layer BP. We found that monolayer BP is a very special material, in which we predict an extremely high hole mobility of 10,00026,000 cm2 V 1 s 1. All of our mobility results can be understood rather well from the form of the real-space wavefunctions, which show a much stronger interlayer coupling in few-layer BP than in graphene or TMDCs. All these results make few-layer BP a very promising candidate material for future applications in electronics and optoelectronics.

Note added in proof: during the review process for this manuscript, we became aware that few-layer BP samples have

very recently been exfoliated from bulk BP and fabricated into electronic devices46,47. Measurements on these samples verify the high hole mobilities46,47 and transport anisotropy46 we predict. Although preliminary optical characterizations have been performed46,47, linear dichroism has not yet been explored.

Methods

DFT calculation. DFT calculations were performed using the generalized gradient approximation for the exchange-correlation potential, the projector augmented wave method48,49 and a plane-wave basis set as implemented in the Vienna abinitio simulation package50. The energy cutoff for the plane-wave basis was set to 500 eV for all calculations. Two k-meshes of 10 8 4 and 10 8 1 were

adopted to sample the rst BZ of the conventional unit cell of bulk and few-layer BP, and the mesh density of k points was kept xed when calculating bandstructures using primitive cells. In optimizing the system geometry, van der Waals interactions were considered by the vdW-DF level with the optB88 exchange functional (optB88-vdW)22,27. The shape and volume of each supercell were optimized fully and all atoms in the supercell were allowed to relax until the residual force per atom was o0.001 eV 1. Electronic bandstructures were calculated by the mBJ23,24 and hybrid functional (HSE06)25,26 methods based on the atomic structures obtained from the full optimization by optB88-vdW. The charge gradient used in the mBJ calculations, c 1.1574, was extracted from the

bulk value.

Carrier mobility calculation. In 2D the carrier mobility is given by the expression4042

m2D

e 3C2DkBTm emdEi12

1

p . The term E1 represents the

deformation potential constant of the valence-band minimum for hole or conduction-band maximum for electron along the transport direction, dened by Ei1 DVi=Dl=l0. Here DVi is the energy change of the ith band under proper cell

compression and dilatation (calculated using a step of 0.5%), l0 is the lattice constant in the transport direction and Dl is the deformation of l0. The elastic modulus C2D of the longitudinal strain in the propagation directions (both x and y) of the longitudinal acoustic wave is derived from (E E0)/S0 C(Dl/l0)2/2, where E

is the total energy and S0 is the lattice volume at equilibrium for a 2D system. All structural properties in the calculation of carrier mobilities were obtained from optB88-vdW and properties related to the electronic structure were computed with the HSE06 functional. The temperature used for the mobility calculations was 300 K.

Uncertainty of deformation potential. In calculations of the electron and hole mobilities, the deformation potential is derived from linear ts to the respective energies of the conduction-band maximum and the valence-band minimum as functions of the lattice dilation or compression. Typical tting results for E1y

(electron state) and E1x (hole state) are shown in Supplementary Fig. 5ac.

The slope represents the deformation potential, which usually has an unavoidable uncertainty accounted for in the standard tting error. Supplementary Figure 5d summarizes the relative errors for the deformation potentials of electrons and holes in few-layer BP for the x and y directions. For a given direction and type of carrier, the error is generally smallest for monolayer BP, while that for ve-layer BP is the largest, except for the hole deformation potential along y. The error inE1y for hole states appears rather larger than the others, especially for monolayer

BP where it exceeds 20%, presumably due to the small absolute value of E1y

for this case. Supplementary Fig. 5c shows details of the t. The absolute value of E1y for hole states is only 0.15 eV, an order of magnitude smaller than theother energies; an error of 1 meV, believed to be the accuracy limit of DFT,may then lead to a tting error of more than 10%. In this respect, a ttingerror of 20% is denitely reasonable, and suggests that our calculation is highly accurate.

Absorption spectra calculation. Absorption spectra were calculated from the dielectric function using the expression A(o) 1 e a(o) Dz, where ao oecn is the absorption coefcient, n

where m e is the effective mass in the transport direction and mdis the average effective mass determined by md

m xm y

q is the index of refraction, e1 and e2 are the real and imaginary parts of the dielectric function, o is the light frequency, c is the speed of light in vacuo and Dz represents the unit-cell size in the c direction. The electronic structures were obtained from HSE06 results and the k-mesh was doubled in calculating dielectric functions. Excitonic contributions were not considered in our calculations. The total number of bands considered was set to be twice that used in the total-energy and bandstructure calculations. Because the dielectric function is a tensor, the absorption spectra along the three directions a (x), b (y) and c (z) were obtained separately.

e e

p e

2

6 NATURE COMMUNICATIONS | 5:4475 | DOI: 10.1038/ncomms5475 | http://www.nature.com/naturecommunications

Web End =www.nature.com/naturecommunications

& 2014 Macmillan Publishers Limited. All rights reserved.

NATURE COMMUNICATIONS | DOI: 10.1038/ncomms5475 ARTICLE

References

1. Katsnelson, M. I. Carbon in Two Dimensions (Cambridge University Press,

2012).

2. Novoselov, K. S. et al. Two-dimensional gas of massless Dirac fermions in graphene. Nature 438, 197200 (2005).

3. Geim, A. K. & Novoselov, K. S. The rise of graphene. Nat. Mater. 6, 183191 (2007).

4. Hwang, E. & Das Sarma, S. Acoustic phonon scattering limited carrier mobility in two-dimensional extrinsic graphene. Phys. Rev. B 77, 115449 (2008).

5. Geim, A. K. & Grigorieva, I. V. Van der Waals heterostructures. Nature 499, 419425 (2013).

6. Wang, Z. M. MoS2: Materials, Physics, and Devices (Springer International Publishing, 2014).

7. Wang, Q. H. et al. Electronics and optoelectronics of two-dimensional transition metal dichalcogenides. Nat. Nanotech. 7, 699712 (2012).

8. Fivaz, R. & Mooser, E. Mobility of charge carriers in semiconducting layer structures. Phys. Rev. 163, 743755 (1967).

9. Vogt, P. et al. Silicene: compelling experimental evidence for graphenelike two-dimensional silicon. Phy. Rev. Lett. 108, 155501 (2012).

10. Houssa, M. et al. Electronic properties of hydrogenated silicene and germanene. Appl. Phys. Lett. 98, 223107 (2011).

11. Bianco, E. et al. Stability and exfoliation of germanane: a germanium graphane analogue. ACS Nano 7, 44144421 (2013).

12. Berger, C. et al. Ultrathin epitaxial graphite: 2D electron gas properties and a route toward graphene-based nanoelectronics. J. Phys. Chem. B 108, 1991219916 (2004).

13. Liao, L. et al. High-speed graphene transistors with a self-aligned nanowire gate. Nature 467, 305308 (2010).

14. Schwierz, F. Graphene transistors. Nat. Nanotechnol. 5, 487496 (2010).15. Radisavljevic, B. et al. Single-layer MoS2 transistors. Nat. Nanotech. 6, 147150 (2011).

16. Wang, H. et al. Integrated circuits based on bilayer MoS2 transistors. Nano Lett. 12, 46744680 (2012).

17. Yoon, Y., Ganapathi, K. & Salahuddin, S. How good can monolayer MoS2 transistors be? Nano Lett. 11, 37683773 (2011).

18. Popov, I., Seifert, G. & Tomnek, D. Designing electrical contacts to MoS2 monolayers: a computational study. Phys. Rev. Lett. 108, 156802 (2012).

19. Xia, F., Farmer, D. B., Lin, Y. M. & Avouris, P. Graphene eld-effect transistors with high on/off current ratio and large transport band gap at room temperature. Nano Lett. 10, 715718 (2010).

20. Perdew, J. P., Burke, K. & Ernzerhof, M. Generalized gradient approximation made simple. Phy. Rev. Lett. 77, 38653868 (1996).

21. Grimme, S. Semiempirical GGA-type density functional constructed with a long-range dispersion correction. J. Comput. Chem. 27, 17871799 (2006).

22. Klime, J., Bowler, D. R. & Michaelides, A. Van der Waals density functionals applied to solids. Phys. Rev. B 83, 195131 (2011).

23. Becke, A. D. & Johnson, E. R. A simple effective potential for exchange.J. Chem. Phys. 124, 221101 (2006).24. Tran, F. & Blaha, P. Accurate band gaps of semiconductors and insulators with a semilocal exchange-correlation potential. Phys. Rev. Lett. 102, 226401 (2009).

25. Heyd, J., Scuseria, G. E. & Ernzerhof, M. Hybrid functionals based on a screened Coulomb potential. J. Chem. Phys. 118, 82078215 (2003).

26. Heyd, J., Scuseria, G. E. & Ernzerhof, M. Erratum: Hybrid functionals based on a screened Coulomb potential. J. Chem. Phys. 124, 219906 (2006).

27. Klimes, J., Bowler, D. R. & Michaelides, A. Chemical accuracy for the van der Waals density functional. Condens. Matter Phys. 22, 022201 (2010).

28. Hultgren, R., Gingrich, N. S. & Warren, B. E. The atomic distribution in red and black phosphorus and the crystal structure of black phosphorus. J. Chem. Phys. 3, 351355 (1935).

29. Brown, A. & Rundqvist, S. Renement of the crystal structure of black phosphorus. Acta Cryst. 19, 684685 (1965).

30. Cartz, L. et al. Effect of pressure on bonding in black phosphorus. J. Chem. Phys. 71, 17181721 (1979).

31. Keyes, R. The electrical properties of black phosphorus. Phys. Rev. 92, 580584 (1953).

32. Warschauer, D. Electrical and optical properties of crystalline black phosphorus. J. Appl. Phys. 34, 18531860 (1963).

33. Maruyama, Y., Suzuki, S., Kobayashi, K. & Tanuma, S. Synthesis and some properties of black phosphorus single crystals. Physica B C 105, 99102

(1981).34. Akahama, Y., Endo, S. & Narita, S.-I. Electrical properties of black phosphorus single crystals. J. Phys. Soc. Jpn 52, 21482155 (1983).

35. Morita, A., A., H., Kaneta, C. & Sasaki, T. Anisotropic mobility and electron-phonon interaction in black phosphorus. Proc. 17th Int. Conf. Phys. Semicond. 13201324 (1985).

36. Walukiewicz, W., Ruda, H., Lagowski, J. & Gatos, H. Electron mobility in modulation-doped heterostructures. Phys. Rev. B 30, 45714582 (1984).

37. Kaasbjerg, K., Thygesen, K. S. & Jacobsen, K. W. Phonon-limited mobility in n-type single-layer MoS2 from rst principles. Phys. Rev. B 85, 115317 (2012).

38. Du, Y., Ouyang, C., Shi, S. & Lei, M. Ab initio studies on atomic and electronic structures of black phosphorus. J. Appl. Phys. 107, 093718 (2010).

39. Rodin, A. S., Carvalho, A. & Castro Neto, A. H. Strain-induced gap modication in black phosphorus. Phys. Rev. Lett. 112, 176801 (2014).

40. Bruzzone, S. & Fiori, G. Ab-initio simulations of deformation potentials and electron mobility in chemically modied graphene and two-dimensional hexagonal boron-nitride. Appl. Phys. Lett. 99, 2108 (2011).

41. Takagi, S.-i., Toriumi, A., Iwase, M. & Tango, H. On the universality of inversion layer mobility in Si MOSFETs: part I-effects of substrate impurity concentration. IEEE Trans. Electr. Dev. 41, 23572362 (1994).

42. Fiori, G. & Iannaccone, G. Multiscale modeling for graphene-based nanoscale transistors. Proc. IEEE 101, 16531669 (2013).

43. Bardeen, J. & Shockley, W. Deformation potentials and mobilities in non-polar crystals. Phys. Rev. 80, 7280 (1950).

44. Xi, J. et al. First-principles prediction of charge mobility in carbon and organic nanomaterials. Nanoscale 4, 43484369 (2012).

45. Liu, L., Feng, Y. & Shen, Z. Structural and electronic properties of h-BN. Phy. Rev. B 68, 104102 (2003).

46. Liu, H. et al. Phosphorene: an unexplored 2D semiconductor with a high hole mobility. ACS Nano 8, 40334041 (2014).

47. Li, L. et al. Black phosphorus eld-effect transistors. Nat. Nanotech. 9, 372377 (2014).

48. Blchl, P. E. Projector augmented-wave method. Phys. Rev. B 50, 1795317979 (1994).

49. Kresse, G. & Joubert, D. From ultrasoft pseudopotentials to the projector augmented-wave method. Phys. Rev. B 59, 17581775 (1999).

50. Kresse, G. & Furthmller, J. Efcient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys. Rev. B 54, 1116911186 (1996).

Acknowledgements

We thank B. Normand for a critical reading of the manuscript. This work was supported by the National Natural Science Foundation of China (NSFC) under Grant Nos 11004244 and 11274380, the Ministry of Science and Technology (MOST) of China under Grant No. 2012CB932704, the Beijing Natural Science Foundation (BNSF) under Grant No. 2112019 and the Basic Research Funds of Renmin University of China from the Central Government under Grant Nos. 12XNLJ03, 14XNH060 and 14XNH062. W.J. was supported by the Program for New Century Excellent Talents in Universities. Calculations were performed at the Physics Laboratory for High-Performance Computing of Renmin University of China and at the Shanghai Supercomputer Center.

Author contributions

W.J. conceived this research. J.Q., Z-X.H. and W.J. performed atomic and electronic structure calculations. X.K. calculated deformation potentials, elastic moduli and carrier mobilities. F.Y. and W.J. predicted the optical absorption spectra. J.Q., X.K. and W.J. wrote the manuscript and all the authors commented on it.

Additional information

Supplementary Information accompanies this paper at http://www.nature.com/naturecommunications

Web End =http://www.nature.com/ http://www.nature.com/naturecommunications

Web End =naturecommunications

Competing nancial interests: The authors declare no competing nancial interests.

Reprints and permission information is available online at http://npg.nature.com/reprintsandpermissions/

Web End =http://npg.nature.com/ http://npg.nature.com/reprintsandpermissions/

Web End =reprintsandpermissions /

How to cite this article: Qiao, J. et al. High-mobility transport anisotropy and linear dichroism in few-layer black phosphorus. Nat. Commun. 5:4475doi: 10.1038/ncomms5475 (2014).

This work is licensed under a Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License. The images or other third party material in this article are included in the articles Creative Commons license, unless indicated otherwise in the credit line; if the material is not included under the Creative Commons license, users will need to obtain permission from the license holder to reproduce the material. To view a copy of this license, visit http://creativecommons.org/licenses/by-nc-nd/4.0/

Web End =http:// http://creativecommons.org/licenses/by-nc-nd/4.0/

Web End =creativecommons.org/licenses/by-nc-nd/4.0/

NATURE COMMUNICATIONS | 5:4475 | DOI: 10.1038/ncomms5475 | http://www.nature.com/naturecommunications

Web End =www.nature.com/naturecommunications 7

& 2014 Macmillan Publishers Limited. All rights reserved.