† Corresponding author. E-mail:

Project supported by the State Key Program for Basic Researches of China (Grant No. 2014CB921103) and the Natural Science Foundation of Jiangsu Province, China (Grant Nos. BK20141441 and BK2010012).

According to first-principles calculations, it is our prediction that bilayer phosphorene (BLP) will become a quasi-two-dimensional superconductor under a certain degree of interlayer compression. A decreasing interlayer distance may realize the transition in the BLP from a semiconducting phase to a metallic phase. On the other hand, a severe vertical compression may make the BLP lattice become dynamically unstable. It is found that in the stable metallic phase of the BLP, interlayer phonon modes dominate the electron-phonon coupling *λ*. The obtained *λ* can be greater than 1 and the superconducting temperature *T*_{c} can be higher than 10 K.

Initiated by the isolation of graphene,^{[1]} research on two-dimensional (2D) materials has aroused wide interests due to the potential for significant and wide-ranging novel applications. By exfoliating a layered crystal into single layers, many 2D materials^{[2–4]} have been obtained. The structure of black phosphorus (BP) consists of puckered double layers, which are held together by weak van der Waals (vdW) force.^{[5]} Phosphorene (few-layer BP) was a newly found 2D crystal by the exfoliation from a bulk crystal and used to create field-effect transistors.^{[6–8]} Large current on-off ratios, high mobilities up to 3900 cm^{2}/V·s–4000 cm^{2}/V·s, and quantum oscillations have been demonstrated in few-layer phosphorene devices.^{[9,10]} It was reported that the phosphorenes possess superior mechanical flexibility and can sustain strain up to 30% and even 32%.^{[11]} First-principles calculations suggested that strain could induce a direct-indirect band gap semiconductor transition or a semiconductor–semimetal transition in phosphorene.^{[11–14]} Furthermore, Tan *et al.*^{[15]} showed that the doped phosphorenes with carbon, oxygen, and sulfur are ferromagnetism, which may play an important role in developing spin electric devices.

The investigation about the superconductivity of 2D materials^{[16,17]} has stirred much interest due to possible applications, such as nanoscale superconducting transistors and nanoscale superconducting quantum interference devices. It was reported early that the single crystal of bulk BP showed superconductivity with *T*_{c} higher than 10 K under high pressures.^{[18,19]} Recently, the monolayer phosphorene by a theoretically assumed electron doping with a compensating uniform positive background was suggested^{[20]} to be a potential monolayer superconductor by first-principles simulations. Very recently, we predicted that the superconductivity can be induced in bilayer phosphorene (BLP) by Li intercalation, and the transition temperature can be increased up to 16.5 K by first-principles calculations.^{[21]} The early experiment^{[22]} showed that with increasing pressure, the distance between puckered layers decreases much faster than the in-plane lattice constants, for the layered BP is a strong anisotropic structure. As a result, it is interesting to see if few-layers BP can become superconducting under vertical compression (i.e., decreasing the interlayer distance). The BLP is the thinnest few-layer BP system that can be used to study the effect of changing interlayer distance.

In this work, electronic structure, lattice dynamics and electron-phonon (EP) coupling of BLP systems under vertical stress are studied via first-principles density functional theory (DFT). The calculated results show that the vertical stress on the BLP can induce a transition from a semiconductor to a BCS superconductor with the superconducting temperature *T*_{c} up to 10 K. The interlayer phonon modes play an important role in enhancing the EP coupling of BLP with decreasing the interlayer distance.

The calculations are performed in a plane-wave pseudopotential representation through the PWSCF program of the Quantum-ESPRESSO distribution.^{[23]} The ultrasoft pseudopotentials are used to model the electron-ion interactions. The valence electronic wave functions are expanded in a plane-wave basis set with a kinetic energy cutoff of 30 Ry (1 Ry = 13.6056923 eV). The general gradient approximation (GGA-PBE) is used for the exchange and correlation energy functional. The BLP is modeled using a series of slabs, which are separated from each other by a vacuum layer of thickness of about 35 Å. The vdW correction proposed by Grimme (DFT-D2) is included to model the interaction between the two layers of BLP. For the electronic-structure calculations, the Brillouin-zone integrations are performed by using the Gaussian smearing technique with a width of 0.02 Ry. Within the framework of the linear-response theory, the dynamical matrixes and the electron-phonon interaction coefficients are calculated.

BLP consists of two puckered layers, which are interacted by weak vdW forces and are stacked following an *AB* stacking order. The top and side views of the atomic structure of BLP are presented in Figs. *x* and *y* directions, respectively. For bulk BP, the early experiment^{[22]} showed that the vdW bonds are shortened while the covalent bonds change slightly with increasing pressure. As a result, the BLP under vertical pressure is modeled by decreasing the interlayer distance. We first take the interlayer distance *D* as fixed, and then relax the unit cell and other atomic positions to obtain their optimized values. With the decrease of *D*, the calculated lattice constants, puckered distance *d*_{z}, bond length and bond angles are shown in Figs. *b*, bond length *L*, and bond angle *θ*_{1} within the layer hardly change with *D*, reflecting the rigidity of the strong covalent bonding along the zigzag direction. For lattice constant *a* and bond angle *θ*_{2}, however, they show different responses to the vertical compression due to the anisotropic geometric structure. For 0.75*D*_{0} ≤ *D* ≤ *D*_{0}, the lattice constant *a* and the bond angle *θ*_{2} hardly change with *D*, but for *D* ≤ 0.75*D*_{0}, they gradually increase with decreasing the interlayer distance. The puckered distance *d*_{z} shows little change, indicating that the puckered character of BLP remains basically unchanged in the range of interlayer distance studied here.

The calculated band structures of BLP for several interlayer distances are shown in Fig. ^{3} hybridization. The top of the valence band has *π*-like character, while the bottom of the conduction band has *π**-like character. For the free BLP (Fig. *D* = 0.95*D*_{0}. The minimum of the conduction band moves from point *Ȳ* and *D*, the indirect band gap decreases gradually and it goes to zero at *D* = 0.85*D*_{0} (see Fig. *D* will result in the transition from a semiconductor to a semimetal (see Fig. *N*(*E*_{F}) listed in Table *D*.

It is well known that DFT underestimates the band gap of semiconductors. Other methods incorporated many-body corrections, such as hybrid functionals,^{[24]} GW approximation^{[25]} seem to correct the systematic DFT underestimation. We use the hybrid functional Gau-PBE^{[26]} to calculate the electronic structure of free BLP, the band gap is estimated to be about 1.0 eV, which is larger than the band gap in the conventional DFT. Although the hybrid functional approach results in an appreciably larger band gap, the observed trend of electronic structure under vertical compression with DFT investigations above is in line with the hybrid functional approach.

In order to study lattice dynamical stability of the metallic BLP, we calculate its phonon spectra. The calculated phonon dispersions and the phonon DOS for several *D* are shown in Fig. *D* = 0.8*D*_{0}, 0.75*D*_{0}, and 0.7*D*_{0}, their phonon spectra and phonon DOS are very similar, reflecting that the change of vdW force has only a little effect on the lattice dynamics. For *D* = 0.65*D*_{0}, however, there appear imaginary frequencies in the two acoustic branches, as shown by the dotted lines of Fig. *D* is smaller than 0.65*D*_{0}, the bilayer structure will become more unstable. As a result, the range of *D* for the stable metallic BLP should be within the range of 0.65*D*_{0} < *D* < 0.85*D*_{0}.

In what follows we discuss the EP interaction. According to the Migdal–Eliashberg theory,^{[26]} the Eliashberg spectral function *α*^{2}*F*(*ω*) is given by

*N*(

*E*

_{F}) is the electronic DOS at Fermi level,

^{[23]}By plotting

*α*

^{2}

*F*(

*ω*), we can estimate the relative strength of the EP coupling. The frequency-dependent EP coupling is given by

*α*

^{2}

*F*(

*ω*) and frequency-dependent EP coupling

*λ*(

*ω*) are shown in Fig.

*α*

^{2}

*F*(

*ω*) has a significantly enhanced peak at the frequency below 55 cm

^{−1}, which will dominate the EP coupling. For example, at

*D*= 0.7

*D*

_{0},

*λ*(55 cm

^{−1}) ≃ 1.08 is beyond 70% of the total EP coupling

*λ*(∞)≃ 1.46. The characteristic phonon modes which dominate the EP coupling are mainly from the interlayer modes. According to the previous first-principles calculation by Du

*et al.,*

^{[27]}the interlayer interaction is the vdW Keesom force, and the permanent dipole is formed within each puckered layer. It is speculated that this interlayer attractive interaction plays an important role in enhancing the EP coupling of BLP for smaller interlayer distance.

For three different *D*, the obtained EP coupling constant *λ* and logarithmically averaged frequency *ω*_{ln}, as well as the estimated *T*_{c} with *μ** = 0.1 taken are summarized in Table *λ* increases dramatically with the decrease of *D*. The increase of both *N*(*E*_{F}) and *λ* is favorable to enhancing *T*_{c}. When *D* = 0.7*D*_{0}, *λ* reaches as high as 1.46 and the estimated *T*_{c} is about 10 K. The superconductivity of bulk BP was reported early under high pressure.^{[18,19]} Our results suggest that BLP can become a BCS superconductor under a certain degree of vertical compression. Due to the decrease of the van der Waals gap, the conduction band and the valence band are overlapped, resulting in the semiconductor-metal transition. In our previous work,^{[21]} we predicted that the superconductivity could be induced in BLP by Li intercalation, and the transition temperature could be increased up to 16.5 K. The intercalated Li atoms act as a charge reservoir of the BLP and the charge transfer from Li to P atoms leads to a semiconducting to metallic transition. The prediction of superconductivity in BLP induced by vertical compression or by Li intercalation is an interesting subject for future experimental study.

In conclusion, we predict that the BLP can become superconducting by adjusting the interlayer distance. A decreasing interlayer distance may realize the transition in the BLP from a semiconducting phase to a metallic phase. Our calculations of phonon structures show that the stable metallic BLP under vertical compression is within the range of 0.65*D*_{0} < *D* < 0.85*D*_{0}. The interlayer phonon modes dominate the EP coupling, and so the interlayer vdW attractive interaction may play an important role in enhancing the EP coupling of BLP at smaller interlayer distance. The obtained EP coupling *λ* can be greater than 1 and the superconducting temperature *T*_{c} can be higher than 10K, suggesting that the phosphorene is a good candidate for low-dimensional superconductor.

**Reference**

1 | Science 306 666 |

2 | Chin. Phys. B 24 097103 |

3 | Chin. Phys. B 24 087201 |

4 | Chin. Phys. Lett. 32 037302 |

5 | J. Chem. Phys. 3 351 |

6 | Nat. Nanotech. 9 372 |

7 | ACS Nano 8 4033 |

8 | Appl. Phys. Lett. 104 103106 |

9 | Nat. Nanotechnol. 10 608 |

10 | 2D Materials 2 011001 |

11 | Appl. Phys. Lett. 104 251915 |

12 | Nanotechnology 25 455703 |

13 | Phys. Rev. Lett. 112 176801 |

14 | J. Phys.: Condens. Matter 27 175006 |

15 | |

16 | Chin. Phys. B 24 117404 |

17 | Chin. Phys. Lett. 31 017401 |

18 | Solid State Commun. 49 879 |

19 | Science 160 994 |

20 | Europhys. Lett. 108 67004 |

21 | Appl. Phys. Lett. 106 113107 |

22 | Phys. Rev. B 33 6177 |

23 | J. Phys.: Condens. Matter 21 395502 |

24 | Phys. Chem. Chem. Phys. 10 6615 |

25 | Phys. Rev. Lett. 93 126406 |

26 | |

27 | J. Appl. Phys. 107 093718 |