Mechanical and Chemical Stability of Monolayer Black Phosphorous

3 hours ago - Download Citation · Email a Colleague · Order Reprints · Rights & Permissions · Citation Alerts · Add to ACS ChemWorx. SciFinder Subscri...
0 downloads 0 Views 2MB Size
Subscriber access provided by Kaohsiung Medical University

C: Surfaces, Interfaces, Porous Materials, and Catalysis

Mechanical and Chemical Stability of Monolayer Black Phosphorous Studied by Density Functional Theory Simulations Mohammad J. Eslamibidgoli, and Michael H. Eikerling J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.8b04344 • Publication Date (Web): 07 Sep 2018 Downloaded from http://pubs.acs.org on September 7, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Mechanical and Chemical Stability of Monolayer Black Phosphorous Studied by Density Functional Theory Simulations Mohammad Javad Eslamibidgolia and Michael H. Eikerlinga,* a

Department of Chemistry, Simon Fraser University, 8888 University Drive, Burnaby,

V5A 1S6, BC, Canada

Abstract Simulations based on electronic density functional theory have been employed to study the environmental stability of phosphorene under mechanical stress as well as oxidative conditions. To understand the mechanical response, bi-axial strain was applied along zigzag and arm-chair directions and the potential energy surface was generated. Poisson’s ratio and Young’s modulus were calculated along each direction indicating anisotropic response of the material. Under large strain conditions, several stable or metastable phases were identified including transformation from black phosphorus to white phosphorus and polymeric phases. To evaluate the chemical stability, surface mixing energies of phosphorene oxide were calculated as a function of oxygen coverage. Results indicate the formation of PO3 and PO4 chains at oxygen coverage above 0.5 monolayers, suggesting a multistep oxidation process that would ultimately lead to the formation of P2O5. Ab initio molecular dynamics simulations with an additional water molecule revealed the hydrophobic nature of pristine black phosphorus in comparison to the hydrophilic nature of oxidized black phopshorus.

1 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 31

Introduction Monolayer (or few-layer) black phosphorus (BP) has attracted numerous studies in recent years due to its unusual electronic and optical properties1. As the most thermodynamically stable allotrope of phosphorus, BP possesses a puckered orthorhombic structure with a thickness-dependent direct band gap, which spans a range from 0.3 eV for bulk BP to 0.9 eV for single-layer BP2-4. The material has a high carrier mobility (up to 1000 cm2V−1s−1 at room temperature)5,6 and it can be turned into either a p-type or an n-type semiconductor via using contact metal engineering7,8. BP demonstrates mechanical, electrical, optical and thermal anisotropic properties which make this material very promising for a wide range of applications, including photonic, thermoelectric, telecommunication, and photovoltaic devices9-13. Phosphorus-based composites have also recently moved into the focus as anode materials for sodium ion batteries, thanks to their high theoretical capacity of 2600 mA h g−1

14,15

.

Each P atom can react with three Na ions during electrochemical sodiation to form Na3P, leading to the highest capacity for Na ion batteris to be attained with any of the known anode materials16-20. A major drawback of BP-based materials is the lack of thermodynamic stability in oxygen-rich atmosphere, which leads to rapid compositional and structural changes and chemical degradation21,22. Huang et al. studied the role of oxygen and water in the degradation of BP21; their results confirmed that deaerated water has negligible impact on BP degradation, while O2-saturated water decomposes the material21. It was found in experiments that the oxidation of BP leads to a significant increase in surface roughness followed by the formation of large droplets of phosphorus/phosphoric acid1,22-24.

2 ACS Paragon Plus Environment

Page 3 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

However, these structural transformations of oxidized BP have not yet been captured in theoretical studies25,26. Several research groups have also explored the response of phosphorene to mechanical strain, which determines its electronic and optical properties for use in nanoelectronic and optoelectronic devices27. In a density functional theory (DFT) study, Wei et al.28 found that a monolayer of phosphorene withstands a tensile strain of up to 30% and 27% in the armchair and zigzag directions, respectively. Sha et al.29 studied the mechanical properties of phosphorene at finite temperatures using molecular dynamics (MD) simulations. They found that fracture strength and fracture strain decreased significantly as the temperature increased from 0 K to 450 K. First-principles calculations of Jiang and Park30 showed that Young’s modulus and strain responses are anisotropic due to the puckered structure of phosphorene. In the majority of the theoretical studies, only uniaxial strain was applied along each direction. Theoretical studies of both uniaxial and biaxial strain effects, in a wider range, on the possible structural transformation of phosphorene is therefore lacking. Further theoretical works are required to gain deeper understanding of the environmental stability of monolayer BP and its responses to varying mechanical and chemical conditions. In the first part of this paper, we use DFT calculations to generate the potential energy surface (PES) of phosphorene under a wide range of bi-axial strain to uncover the mechanical properties and structural evolution of this 2D material. In the second part, we use DFT to study the oxidation and consequent chemical degradation of phosphorene. Based on a wide exploration of surface energy for various possible configurations of BP oxide, we propose new surface structures for oxidized phosphorene.

3 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 31

Moreover, using ab inito molecular dynamics (AIMD) we investigate the impact of of surface oxidation on transforming the wettability of phosphorene from hydrophobic to hydrophilic. Computational Methods Electronic structure calculations were carried out with the periodic DFT package VASP31. Kohn-Sham one-electron wave functions were expanded in a plane wave basis set up with energy cut-off of 400 eV. Exchange-correlation effects were incorporated within the generalized gradient approximation (GGA), using the exchange-correlation functional by Perdew, Burke, and Ernzerhof (PBE)32 along with the projector-augmented wave (PAW) potentials33. The reciprocal space was meshed with a 12 × 8 × 1 k-point grid using the Monkhorst-Pack method34. To simulate a monolayer of phosphorene, a unit cell with periodic boundary condition was used. A vacuum space of at least 25 Å was applied to minimize the interaction between layers. In the vibrational frequency calculations, only the surface P and Oad atoms were included and the construction of the Hessian matrices was evaluated by the atomic forces after displacements of ± 0.015 Å form equilibrium. AIMD simulations were performed in the microcanonical ensemble using the Verlet algorithm with a time step of 1 fs at room temperature. Results and Discussion Part 1. Mechanical response and structural transformation of phosphorene under applied strain Black phosphorus has a unique layered structure as a result of the bonding between phosphorus atoms to three neighboring atoms by sp3 hybridized orbitals. This leads to the arrangement of phosphorus atoms in a puckered lattice structure, as shown in Figure 1.

4 ACS Paragon Plus Environment

Page 5 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

The initial structure of phosphorene in our study was obtained from bulk black phosphorus. The relaxed lattice constants for a monolayer of BP were obtained as a = 3.30 Å along the zig-zag direction, and b = 4.63 Å, along the armchair direction, illustrated in Figure 1 (b) and (c), respectively. Given the optimized lattice parameters of phosphorene, we applied a bi-axial strain along zigzag direction from 3.20 Å to 3.40 Å and along armchair directions from 4.45 Å to 4.80 Å, in increments of 0.01 Å. At each axial strain value, the lattice constant in the transverse direction was fully relaxed to ensure the force in the transverse direction is at a minimum. Figure 2 (a) illustrates the PES as a function of applied strain along zig-zag and armchair directions. The elliptical region surrounding the global minimum reveals the Poisson effect in phosphorene; the observed negative correlation in the minimum energy induced by this elliptical shape indicates that compression (or stretch) applied in axial direction causes expansion (or compression) in the transverse direction, respectively. The Poisson’s ratio, which is a measure of the Poisson effect, can be obtained from the minimum energy structures of the PES; Figure 2 (b) shows the mechanical response of the zigzag and armchair directions to the applied strain in their perpendicular axes. The absolute value of the slope is the Poisson’s ratio, which we found to be 0.71 and 0.17 along zigzag and armachair direction, respectively. These values are in agreement with the Poisson’s ratio presented in ref.28. The steeper energy surface and larger Poisson’s ratio in the zig-zag direction suggest that strain is more difficult to apply in this direction than in the armchair direction; correspondingly, the armchair direction is softer than the zig-zag direction.

5 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 31

From the PES, the strain energy density can be obtained as shown in Figure 3 (a). Under small strain, from -4% (compression) to 4% (expansion), the strain energy density behavior is symmetric and the data can be fitted to a polynomial quadratic function, i.e. y = Y x2/2 + b x + c, where Y represents the Young’s modulus. As indicated in Figure 3 (a), Young’s modulus was found to be 194.5 GPa and 46.2 GPa along zigzag and armchair direction, respectivey, showing the large anisotropy consistent with the values of the Poisson’s ratio. In prior computational studies, Y was reported as 166 GPa and 44 GPa in ref.28, and 106.4 GPa and 41.3 GPa in ref.30, for zigzag and armchair directions, respectively. Values found in the three works along armchair direction are consistent, while there is a significant variation of the Y along zigzag direction, which could be due to the different ranges of strain applied and/or the different approaches taken to calculate the Young’s modulus. To determine the region where phosphorene behaves elastically, we applied a larger uniaxial strain along the axial directions, i.e., from -15% (compression) to +40% (expansion), in increments of 5%, and generated the extended strain energy density plot, shown in Figure 3 (b). The nonlinear and assymetric behavior could be fitted well with a quartic function which is consistent with the results in ref.30. Moreover, from the slope of the curves the elastic region can be estimated to be around 25% along the zigzag direction and 30% along the armchair direction. Consistently, first-principles calculations in ref.28 suggested that phosphorene can withstand a tensile strain of up to 27% and 30% along the zigzag and armchair directions, respectively. Next, starting with the optimized phosphorene structure, the biaxial strain was varied over a wider range, i.e., from -40% (compression) to 100% (expansion) with an

6 ACS Paragon Plus Environment

Page 7 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

increment of 5% in both zigzag and armchair directions. Figure 4 shows the PES constructed from single point structural optimizations. The global minimum corresponds to the phosphorene structure, with its elliptical shape spanning a range of around -20% to +30% along each direction. As can be seen, for a small positive strain along the armchair direction (>5%), large uniaxial expansion along the zig-zag direction (>50%) induces a structural transition from black to white phosphorus forming isolated P4 clusters in the unit cell. On the other hand, uniaxial expansion applied along the armchair direction causes separation of material in this direction and transition to a polymeric form, shown as ∆ plateau in Figure 4. However, formation of polymeric phase along the zig-zag direction requires negative tensile strain in the armchair direction, shown as the Λ plateau in Figure 4. Optimized structures for the three plateaus in Figure 4 are shown in Figure S1 in the supporting information (SI). Moreover, our results predict two other metastable phases, named α and γ, which correspond to local minima in the PES. Structures for these phases are shown in Figure S2. Phase transitions between various allotropes of phosphorus, including white, red, black phosphorus, as well the polymeric form of the material under high-pressure and temperature have been reported in the literature35. However, it is also essential to note that under such conditions the size-dependent effects such as crack formation can occur which would not be captured in DFT calculations due to limits of the unit-cell size that could be considered. In order to further analyse the energetics for structural transitions of BP, we generated various potential energy diagrams from Figure 4. Figure 5 shows the potential energy diagram of phosphorene under fixed mechanical strains along the zig-zag direction for

7 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 31

varying strains along the armchair direction. The value of energy is given relative to that of the phosphorene structure, which is considered as zero (global minimum). The plateau region in Figure 5 (a)-(f) at around 0.5 eV represents the white phosphorus region, which corresponds to positive biaxial strains in both the armchair direction (>5%) and zig-zag direction (>50%). For positive strains (>80%) along the zig-zag direction and the negative strain (50%). On the other hand, the two local minima at around 0.7 eV in Figure 5 (n) correspond to the α and γ phases for positive (15%) and negative (-30%) strains along the armchair direction, respectively; in this case, the activation barrier for the transition between the α and γ phases was found to be 0.4 eV. The energy diagram of phosphorene under fixed mechanical strains along armchair direction for varying strains along the zig-zag direction is included in the supporting information (Figure S3).

8 ACS Paragon Plus Environment

Page 9 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 6 shows the band structure of phosphorene and that of the metastable γ phase. Consistent with previous theoretical studies4,25 for monolayer black phosphorus, we obtained a direct band gap of 0.92 eV at the center of Brillouin zone, while for the γ phase the band gap was found to be zero. The corresponding density of states (DOS) is plotted in the supporting information (Figure S4 and S5). Part 2. Phosphorene oxidation In this part, we discuss the interaction of oxygen with phosphorene. We used DFT to calculate the relative mixing energy of oxidized surfaces of phosphorene with oxygen coverage between zero and one monolayer (ML). Based on this energy, we assessed the relative stabilities of various surface configurations at different coverages, identified the most stable structures, and compared our results with structures proposed in the literature. The mixing energy, normalized to the surface area, were calculated with respect to the pristine phosphorene surface and the surface covered with a monolayer of oxygen from the following relation, 

ΔE  = / − 1 −    −  / ,

(1)

where A is the surface area, ΘO the oxygen coverage as the number of oxygen atoms per surface P atoms, EO/BP the total energy of the oxygen-covered surface, Eclean the energy of a clean phosphorene surface, and E1MLO/BP the energy of a phosphorene surface with a monolayer of adsorbed oxygen at on-top sites. For the oxidized BP surfaces, vibrational frequency calculations were performed to evaluate the zero point energy correction and entropy contributions to the free energy, i.e., ∆ZPE-T∆S, respectively (see eq. S1-S5 in SI). As can be seen in Table S1, the calculated correction terms are negligible compared

9 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 31

to the relative mixing energy. Moreover, no imaginary mode was found, confirming the stability of the structure. Figure 7 shows the computed ΔE  of selected configurations as a function of ΘO along with the side and top views of the optimized structures. Starting from the clean surface, as shown in Figure 7 (a), ΘO was increased in increments of 0.125 ML. Consistent with previous theoretical studies25,26, the preferred adsorption site for monoatomic oxygen is the on-top dangling position (Figure 7(b)). As ΘO increases up to 0.50 ML the ΔE  decreases almost linearly. On the other hand, from ΘO = 0.50 to ΘO = 0.75 the formation of PO3 groups is preferred, which is caused by the displacement of P atoms from their original positions in clean phosphorene. We examined various possible configurations of PO3 by slightly shifting P atoms away from the surface, as seen in Figure 7(f), or parallel to the surface, as seen in Fig. 7(g). Moreover, we changed the configuration of O atoms to form PO3 chains, cf. Figure 7 (g), vs. uniformly distributed PO3 groups, or zig-zag structures, as shown in Figure 7(h). As can be seen from the plot of ΔE  , at ΘO = 0.75 the most stable structure corresponds to uniformly distributed PO3 groups. From ΘO = 0.75 up to 1.0, the formation of PO4 is preferred, as shown in Figure 7(l). For this range of ΘO, we also examined various alternative configurations. For example at ΘO = 1.0, we considered structures with O adsorbed at the dangling site, displayed in Figure 7(i), with a PO3 chain formed, see Figure 7 (j), with a PO4 chain, see Figure 7 (k), and with a PO4 zig-zag structure shown in Figure 7 (l). It was found that the zig-zag PO4 structure has the lowest phase formation energy in comparison to the other configurations. This finding differs from the literature, where the proposed stable

10 ACS Paragon Plus Environment

Page 11 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

conformation of phosphorene oxide was the one shown in Figure 7(i)25. Thus, the proposed surface configuration in Figure 7(l) represents the most stable surface structure for oxidized BP reported in the literature to date. The formation of PO3 and then PO4 groups at high oxygen coverage suggests a multi-step oxidation process of black phosphorus that would eventually lead to the formation P2O5, in which the phosphorus atoms are bound by a tetrahedron of oxygen atoms.

A more complete list of all the

considered structures with oxygen adsorbed on one side or both sides of phosphorene is shown in Figure S6 in the supporting information. We determined the electron density difference upon oxygen adsorption on the phosphorene surfaces corresponding to the adsorption-induced charge rearrangement. The electron density difference along z direction was calculated as,   = ∬  where, 

!""

!"" #$#%,

(2)

=  −  +  , is the difference between the electron density of the

oxidized surface configurations and the sum of the electron densities of non-interacting phosphorene and oxygen atoms at the same positions. The generated 

!""

of the most

stable phosphorene oxide systems at ΘO = 0.25, 0.5, 0.75, and 1 is shown in Figure S7. Figure 8 (a) and (b) shows plots of   for varying oxygen coverages. The adsorption of oxygen involves electron transfer from P atoms to Oad atoms which creates electron accumulation around oxygen atoms and electron depletion at P sites, inducing a charge polarization at the interface. Given  , we calculated the amount of charge accumulation and depletion from the following relation, )*

' = (+*  #.

(3)

11 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 31

As can be seen from Figure 8 (c), by increasing the coverage up to 0.5 ML the amount of accumulated (positive q) and depleted (negative q) charges at O sites and P sites increases monotonically, respectively. However, from 0.5 to 1.0 ML, the amount of charge transfer for the most stable configurations exhibits non-monotonic behavior: q increases from 0.5 to 0.75 and it decreases for the most stable structure at 1 ML, i.e., the zig-zag PO4 configuration (Figure 7 (l)). In fact, the charge rearrangement caused by the displacement of P atoms and the formation of PO4 groups significantly modifies  , as shown in Figure 8 (b), compared to the reference structure at 1 ML (i.e., Figure 7 (i)). Part 3. Interaction of phosphorene and phosphorene oxide with water We calculated the binding energy of a water monomer on the clean and various oxidized surfaces of phosphorene using the following relation: ΔE -./ = -./// − / − -./  As shown in Table 1, the water molecule binds weakly to the pristine phophorene surface, with a binding energy of 80 meV if calculated at the PBE level, or 150 meV if using dispersion the corrected PBE-D3/zero functional. For comparison, it should be mentioned that the binding energy of a water molecule to pristine graphene is reported as 280 meV if calculated using the PBE functional36. On the other hand, water molecule interacts more strongly with the oxidized surfaces due to the formation of hydrogen bonds between the water molecule and Oad. The stable structures of the water monomer on (oxidized) phosphorene are shown in Figure 9. To gain insights about the wetting behavior of phosphorene oxide surfaces, AIMD simulations for water monomer interaction with clean and oxidized surfaces were performed at room temperature. Figure 10 shows the trajectory of the oxygen atom of the

12 ACS Paragon Plus Environment

Page 13 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

water molecule as sampled by AIMD. Initially, the structures were optimized with the oxygen atom in water molecule fixed away from the surface. The simulation was then initiated from the optimized structure and the trajectory was generated. As can be seen in Figure 10 (a), the water monomer is highly mobile on clean phosphorene, which indicates the weak interaction of water with this surface and thus reveals its hydrophobic nature. For the oxidized phosphorene, on the other hand, the water molecule forms a hydrogen bond with Oad, inhibiting its mobility. As shown in Figure 10 (b), (c) and (d) at coverages of 0.5, 0.75 and 1.0 ML, respectively, the water molecule spent most of the time near the surface forming a hydrogen bond with Oads. In this case, the final distances (at 8 ps) between the hydrogen atom in the water molecule and the Oads with which it forms the hydrogen bond were obtained as 2.09 Å, 2.04 Å, and 1.86 Å for 0.5, 0.75 and 1.0, respectively. This finding is consistent with the hydrophilic nature of oxidized BP, as also seen in recent experimental work21.

Conclusion The mechanical and chemical stability of phosphorene was studied by DFT-based simulations. Bi-axial strain loadings were applied along both zigzag and armchair directions and the corresponding potential energy surface was generated. Under small strain conditions, an elliptical global energy minimum was obtained, which reveals the anisotropic strain behavior of the material. Under large strain, several structural transitions were found and transformation and activation energies were determined. The global minimum corresponds to phosphorene, which can resist a tensile strain of around

13 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 31

25% and 30% along the zigzag and armchair directions, respectively. Larger bi-axial strain (> 50%) was found to induce a phase transition to white phosphorus (P4 clusters). In a separate study, the oxidation of phosphorene was investigated using DFT. The stability of surface configurations at various oxygen coverages was evaluated using the phase formation energies. At oxygen coverage > 0.5 ML, the formation of PO3 and PO4 was found to occur suggesting the transformation of black phosphorus to phosphorus pentoxide. Wettability trends were analyzed using AMID simulations of water monomer interaction with clean phosphorene and phosphorene oxide. The water molecule monomer was found to be highly mobile on the clean phosphorene surface, implying weak interaction of water with this surface and revealing its hydrophobic nature. On oxidized surfaces, the water monomer forms a hydrogen bond with Oad, restricting the water mobility and confirming the hydrophilic nature of oxidized BP.

Supporting information Supplementary figures, table and equations as refered to in this paper.

Acknowledgment The authors would like to acknowledge the financial assistance from the NSERC Strategic Project Grant (STPGP 494157-16) led by Prof. Xueliang (Andy) Sun from University of Western Ontation, ON, CANADA. Computations were performed using facilities provided by WestGrid and Compute Canada.

References

14 ACS Paragon Plus Environment

Page 15 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(1) Carvalho, A.; Wang, M.; Zhu, X.; Rodin, A. S.; Su, H.; Neto, A. H. C. Phosphorene: from theory to applications, Nat. Rev. Mater. 2016, 1(11), 16061. (2) Xia, F.; Wang, H.; Jia, Y. Rediscovering black phosphorus as an anisotropic layered material for optoelectronics and electronics. 2014, Nat. Commun., 5, 4458. (3) Tran, V.; Soklaski, R.; Liang, Y.; Yang, L; Layer-controlled band gap and anisotropic excitons in few-layer black phosphorus. 2014, Phys. Rev. B, 89(23), 235319. (4) A. Castellanos-Gomez, L. Vicarelli, E. Prada, J. O. Island, K. L. Narasimha-Acharya, S. Blanter, I.; Groenendijk, D. J.; Buscema, M.; Steele, G. A.; Alvarez, J. V.; Zandbergen, H. W. Isolation and characterization of few-layer black phosphorus, 2014, 2D Mater. 1(2), 025001. (5) Li, L.; Yu, Y.; Ye, G.J.; Ge, Q.; Ou, X.; Wu, H.; Feng, D.; Chen, X.H.; Zhang, Y. Black phosphorus field-effect transistors. 2014, Nat. Nanotech, 9(5), p.372. (6) Liu, H.; Neal, A. T.; Zhu, Z.; Luo, Z.; Xu, X.; Tománek, D.; Ye, P. D. Phosphorene: an unexplored 2D semiconductor with a high hole mobility, 2014, ACS Nano, 8.4, 40334041. (7) Perello, D. J.; Chae, S. H.; Song, S.; Lee, Y. H. High-performance n-type black phosphorus transistors with type control via thickness and contact-metal engineering, 2015, Nat. Commun., 6, 7809 (8) Koenig, S. P.; Doganov, R. A.; Schmidt, H.; Castro Neto, A. H.; Özyilmaz, B. Electric field effect in ultrathin black phosphorus, 2014, Appl. Phys. Lett., 104, 103106 (9) Wang, X.; Jones, A.M.; Seyler, K.L.; Tran, V.; Jia, Y.; Zhao, H.; Wang, H.; Yang, L.; Xu, X.; Xia, F. Highly anisotropic and robust excitons in monolayer black phosphorus, 2015, Nature Nanotech. 10.6, 517.

15 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 31

(10) Qiao, J.; Kong, X.; Hu, Z. X.; Yang, F.; Ji, W. High-mobility transport anisotropy and linear dichroism in few-layer black phosphorus. 2014, Nature Commun., 5, 4475. (11) Tao, J.; Shen, W.; Wu, S.; Liu, L.; Feng, Z.; Wang, C.; Hu, C.; Yao, P.; Zhang, H.; Pang, W.; Duan, X. Mechanical and electrical anisotropy of few-layer black phosphorus. 2015, ACS Nano, 9(11), 11362-11370. (12) Luo, Z.; Maassen, J.; Deng, Y.; Du, Y.; Garrelts, R.P.; Lundstrom, M.S.; Peide, D.Y; Xu, X. Anisotropic in-plane thermal conductivity observed in few-layer black phosphorus, 2015, Nat. Commun., 6, 8572. (13) Kou, L.; Frauenheim, T.; Chen, C. Phosphorene as a superior gas sensor: Selective adsorption and distinct I–V response. 2014, J. Phys. Chem. Lett., 5(15), 2675-2681. (14) Xu, G. L.; Chen, Z.; Zhong, G. M. ; Liu, Y.; Yang, Y.; Ma, T.; Amine, K. Nanostructured black phosphorus/Ketjenblack–multiwalled carbon nanotubes composite as high performance anode material for sodium-ion batteries, 2016, Nano Lett., 16(6), 3955-3965 (15) Sun, J.; Lee, H. W.; Pasta, M.; Yuan, H.; Zheng, G.; Sun, Y.; Li, Y.; Cui, Y.; A phosphorene–graphene hybrid material as a high-capacity anode for sodium-ion batteries, 2015, Nat. Nanotec., 10(11), 980 (16) Kim, Y.; Ha, K. H.; Oh, S. M.; Lee, K. T. High‐Capacity Anode Materials for Sodium‐Ion Batteries, 2014, Chem. Eur. J., 20, 11980-11992 (17) Sangster, J. Na-P (Sodium-Phosphorus) System, 2010, J. Phase Equilib. Diffus., 31, 62–67.

16 ACS Paragon Plus Environment

Page 17 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(18) Datta, K.; Epur, R.; Saha, P.; Kadakia, K.; Park, S. K.; Kumta, P. N. Tin and graphite based nanocomposites: Potential anode for sodium ion batteries, 2013, J. Power Sources, 225, 316–322. (19) Qian, J.; Chen, Y.; Wu, L.; Cao, Y.; Ai, X.; Yang, H. High capacity Na-storage and superior cyclability of nanocomposite Sb/C anode for Na-ion batteries, 2012, Chem. Commun., 48, 7070–7072. (20) Kim, Y.; Park, Y.; Choi, A.; Choi, N.S.; Kim, J.; Lee, J.; Ryu, J. H.; Oh, S.M.; Lee, K.T. An amorphous red phosphorus/carbon composite as a promising anode material for sodium ion batteries, 2013, Adv. Mater., 25(22), 3045-3049 (21) Huang, Y.; Qiao, J.; He, K.; Bliznakov, S.; Sutter, E.; Chen, X.; Luo, D.; Meng, F.; Su, D.; Decker, J.; Ji, W. Interaction of black phosphorus with oxygen and water, 2016, Chem. Mater., 28(22), 8330-8339 (22) Wood, J.D.; Wells, S.A.; Jariwala, D.; Chen, K.S.; Cho, E.; Sangwan, V.K.; Liu, X.; Lauhon, L.J.; Marks, T.J.; Hersam, M.C. Effective passivation of exfoliated black phosphorus transistors against ambient degradation, 2014, Nano Lett., 14(12), 6964-6970 (23) Favron, A.; Gaufrès, E.; Fossard, F.; Phaneuf-L’Heureux, A.L.; Tang, N.Y.; Lévesque, P.L.; Loiseau, A.; Leonelli, R.; Francoeur, S.; Martel, R. Photooxidation and quantum confinement effects in exfoliated black phosphorus, 2015, Nat. Mater. 14(8), 826. (24) Yau, S.L.; Moffat, T.P.; Bard, A.J.; Zhang, Z.; Lerner, M. M. STM of the (010) surface of orthorhombic phosphorus, 1992, Chem. Phys. Lett., 198, 383-388. (25) Wang, G.; Pandey, R.; Karna, S. P. Phosphorene oxide: stability and electronic properties of a novel two-dimensional material, 2015, Nanoscale, 7, 524-531.

17 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 31

(26) Ziletti, A.; Carvalho, A.; Campbell, D. K.; Coker, D. F.; Neto, A. C. Oxygen defects in phosphorene, 2015, Phys. Rev. Lett., 114(4), 046801. (27) Akhtar, M.; Anderson, G.; Zhao, R.; Alruqi, A.; Mroczkowska, J.E.; Sumanasekera, G. Jasinski, J.B. Recent advances in synthesis, properties, and applications of phosphorene, 2017, npj 2D Mater. Appl., 1, 5. (28) Wei, Q.; Peng, X. Superior mechanical flexibility of phosphorene and few-layer black phosphorus, 2014, Appl. Phys. Lett., 104, 251915 (29) Sha, Z.D.; Pei, Q.X.; Ding, Z.; Jiang, J.W.; Zhang, Y.W. Mechanical properties and fracture behavior of single-layer phosphorene at finite temperatures, 2015, J. Phys. D: Appl. Phys., 48(39), 395303. (30) Jiang, J. W.; Park, H. S. Mechanical properties of single-layer black phosphorus, 2014, J. Phys. D: Appl. Phys., 47(38), 385304. (31) Kresse, G.; Furthmüller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set, 1996, Phys. Rev. B, 54, 11169. (32) Perdew, J.; Burke, K.; Ernzerhof, M. Generalized gradient approximation made simple, 1996, Phys. Rev. Lett., 77,865. (33) Blöchl, P.E. Projector augmented-wave method, 1994, Phys. Rev. B, 50, 17953. (34) Pack, J. D.; Monkhorst, H. J. Special points for Brillouin-zone integrations, 1977, Phys. Rev. B, 16, 1748. (35) Katayama, Y.; Mizutani, T.; Utsumi, W.; Shimomura, O.; Yamakata, M.; Funakoshi, K. I. A first-order liquid–liquid phase transition in phosphorus. 2000, Nat., 403(6766), 170.

18 ACS Paragon Plus Environment

Page 19 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(36) Freitas, R. R. Q.; Rivelino, R.; Mota, F. D. B.; De Castilho, C. M. C. DFT studies of the interactions of a graphene layer with small water aggregates. 2011, J. Phys. Chemi. A, 115(44), 12348.

Table 1. Binding energy of water monomer on clean and various oxidized surfaces of phosphorene calculated at PBE and PBE-D3/zero levels. Values correspond to the structures shown in Figure 9.

∆Ewater (eV)

∆Ewater (eV)

Coverage (ml)

PBE

PBE-D3/zero

0.0

-0.08

-0.15

0.25

-0.21

-0.31

0.5

-0.19

-0.29

19 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 31

0.75

-0.15

-0.25

1.0

-0.44

-0.61

20 ACS Paragon Plus Environment

Page 21 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 1. Chemical structure of 2D phosphorene (a) Top-view, (b) side-view along armchair and (c) side-view along zigzag directions.

21 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 31

Figure 2. (a) Potential energy surface of phosphorene under strain along zigzag and armchair directions. (b) Strain response of transverse direction to axial strain. Negative of the slope represents the Poisson's ratio along zig-zag and armchair directions.

22 ACS Paragon Plus Environment

Page 23 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 3. (a) Strain energy density vs. applied axial strain along the armchair and zigzag directions fitted to y = Y x2/2 + C1 x + C2; Y represents the Young’s modulus. (b) Strain energy density vs. medium range axial strain along the armchair and zigzag directions fitted to polynomial quartic function.

23 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 31

Figure 4. Potential energy surface of phosphorene under strain along zigzag and armchair directions.

24 ACS Paragon Plus Environment

Page 25 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 5. Potential energy variation of phosphorene under applied mechanical strains along armchair derection for the varying strains along zigzag direction: (a)-(o) from 100% ZZ to -30% ZZ.

25 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 31

Figure 6. Calculated band structure of phosphorene and the meta-stable phase

26 ACS Paragon Plus Environment

Page 27 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 7. Calculated surface mixing energies for various oxygen configurations on phosphorene as a function of oxygen coverage. Energies are relative to clean phosphorene and 1 ML of oxygen adsorbed in the dangling position.

27 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 31

Figure 8. Comparison between the planar average of the induced charge density along the surface normal for phosphorene with oxygen coverages of (a) 0.125 to 0.5 ML and (b) 0.75, and 1 ML; (c) calculated chaege that is transferred between phosphorus and oxygen upon adsorption

28 ACS Paragon Plus Environment

Page 29 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 9. Most stable structures of water monomer on clean and oxidized surfaces of phosphorene. top panel shows it along zigzag direction, lower panel shows it along armchair direction

29 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 31

Figure 10. Trajectory paths of oxygen in water monomer during AIMD simulation at (a) clean phosphorene surface, (b) 0.5 ML O, (c) 0.75 ML O, and (d) 1.0 ML O coverages.

30 ACS Paragon Plus Environment

Page 31 of 31 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

TOC Graphic

31 ACS Paragon Plus Environment