Oxidation-enhanced thermoelectric efficiency in a two-dimensional phosphorene oxide

We performed density functional theory calculations to investigate the thermoelectric properties of phosphorene oxide (PO) expected to form by spontaneous oxidation of phosphorene. Since thermoelectric features by nature arise from the consequences of the electron-phonon interaction, we computed the phonon-mediated electron relaxation time, which was fed into the semiclassical Boltzmann transport equation to be solved for various thermoelectric-related quantities. It was found that PO exhibits superior thermoelectric performance compared with its pristine counterpart, which has been proposed to be a candidate for the use of future thermoelectric applications. We revealed that spontaneous oxidation of phosphorene leads to a significant enhancement in the thermoelectric properties of n-doped phosphorene oxide, which is attributed to the considerable reduction of lattice thermal conductivity albeit a small decrease in electrical conductivity. Our results suggest that controlling oxidation may be utilized to improve thermoelectric performance in nanostructures, and PO can be a promising candidate for low-dimensional thermoelectric devices.

Over the last few decades, much attention has been paid to some particular materials that efficiently convert heat into electricity or vice versa for their potential applications in thermoelectric power generating or cooling systems 1,2 . Their thermoelectric efficiencies can be assessed at temperature T by a dimensionless figure of merit where S, σ , and κ are Seebeck coefficient, electrical and thermal conductivities. The thermal conductivity can be further decomposed into the electronic and lattice thermal conductivities κ e and κ l , or κ = κ e + κ l . A promising thermoelectric candidate with a large ZT value should concurrently possess incompatible characteristics: high electrical conduction with a high σ and efficient thermal insulation with a low κ , not to mention high Seebeck coefficient, which is usually observed in semiconductors. Intriguingly, such a conflicting condition has recently been realized in some chalcogenides such as bismuth telluride ( Bi 2 Te 3 ) 1 , antimony telluride ( Sb 2 Te 3 ) 3 , tin selenide (SnSe) 4,5 , and copper selenide ( Cu 2 Se) 6,7 . Although some of these materials are currently being used in thermoelectric devices, their toxicity and/or expensive materials cost will hinder future thermoelectric applications 8 .
As alternatives that would provide inexpensive and environmentally-friendly thermoelectricity, two-dimensional (2D) materials have drawn much attention since, in general, they tend to have high values of Seebeck coefficient originating from the effect of the low-dimensionality, such as an abrupt change in the density of states (DOS) 9 and a quantum confinement effect 10,11 . Nevertheless, it turned out that some 2D materials failed to satisfy such expectations. For instance, the Seebeck coefficient of graphene is small because of its metallic nature 12 , and MoS 2 exhibited a low ZT value around 0.05 attributed to its relatively large thermal conductivity 13 .
After the discovery of a single layer of black phosphorus or phosphorene 14 , there have been several theoretical studies proposing that phosphorene would be a candidate material for thermoelectric applications [15][16][17][18] . Based on numerical calculations performed within the constant relaxation time approximation using the deformation potential theory, they reported that phosphorene, which is environmentally friendly or non-toxic, has not www.nature.com/scientificreports/ only a large Seebeck coefficient originating from its semiconductivity, but also highly anisotropic transport behaviors of high electrical and low thermal conductivities along a certain direction caused by its intrinsic puckered structure 15,16 . An explicit inclusion of electron-phonon interactions, however, revealed that those reported room-temperature ZT values were an order of magnitude overestimated 19 , implying that phosphorene may not be a suitable choice for thermoelectric applications. To make it worse, phosphorene is known to be structurally unstable under the atmospheric conditions and thus can be easily oxidized due to the lone pair electrons of each phosphorus atom [20][21][22][23][24] .
Here we propose that spontaneous oxidation of phosphorene may be a resolution of the issues in phosphorene for a thermoelectric purpose. To correctly describe the thermoelectric properties of both phosphorene and PO or to calculate their ZT values accurately, we explicitly evaluated the electron relaxation times τ n,k for given electronic band indices n and electron wave vectors k in the presence of electron-phonon interactions. Then a set of {τ n,k } was used to solve the Boltzmann transport equation based on the first-principles density functional theory to obtain the electronic contributions to the quantities in Eq. (1). Combining them with their lattice thermal conductivity calculated in our earlier study 24 , we determined their thermoelectric-related characteristics. It was found that the p-type phosphorene exhibits better thermoelectricity along the armchair direction but worse along the zigzag direction than PO, which shows a much smaller directional dependence. Furthermore, we found that the n-type PO is thermoelectrically superior not only to the n-type phosphorene regardless of their directions but also to the p-type phosphorene. This result is attributed to a significant reduction in κ l originating from spontaneous oxidation 24 albeit a small decrease in power factor S 2 σ . ZT values of PO were estimated to range from 0.1 to 0.5 at T = 300 to 700 K in the electron doping concentration. Therefore, n-doped PO would be a promising and eco-friendly thermoelectric material, and controlling the degree of oxidation may be a useful strategy to improve the thermoelectric properties of the material.

Results and discussion
Various oxidation processes of phosphorene were reported earlier [20][21][22][23][24] . Our previous studies 23,24 showed that fully-oxidized phosphorene or PO possesses highly-anisotopic configuration similar to its pristine counterpart, as shown in Fig. 1. We revealed that PO with a direct bandgap of 0.83 eV possesses two nonsymmorphic symmetries with the inversion symmetry broken, guaranteeing symmetry-protected band structures including the band degeneracy and four-fold degenerate Dirac points 23 . It was also shown that the flexible bonds between phosphorus and oxygen atoms in PO lead to a significant decrease in the lattice thermal conductivity κ PO l , which was evaluated to be in the range between 2.4 and 7.1 W/mK depending on the transport direction, compared to that of pristine phosphorene κ P l in the range of 65 and 146 W/mK 24 . This reduction implies that PO would provide a favorable condition for thermoelectric performance.
To consider the effect of the electron-phonon interaction on transport properties such as σ and κ , we evaluated the electron-phonon scattering rate and the carrier relaxation time, which is the inverse of the former, of both phosphorene and PO for the all available initial states |n, k� . According to the Fermi's golden rule, the electron scattering rate is proportional to the number of available final states satisfying conservation conditions. As seen in Fig. 2a, the electron-phonon scattering rates are somewhat roughly proportional to their corresponding total DOSs, implying the modification of DOS near band edges would be a viable strategy for enhancing the transport properties. We further observed that the carrier relaxation time of PO is larger near the lowest conduction band but much smaller near the highest valance band than that of phosphorene, as shown in Fig. 2b. This indicates that n-type PO would be a better choice for a thermoelectric material than p-type PO, because the carrier relaxation time is closely associated with the spectral conductivity tensor as described in Eq. (3), and thus with the conductivity σ.
We used the calculated relaxation time to evaluate the carrier mobility tensor µ of phosphorene and PO as a function of carrier density concentration at various temperatures, by dividing the electrical conductivity tensor σ by carrier density. Figure 3 shows the evaluated carrier mobility, µ a M,X , of phosphorene (M = P) and PO (M = PO) as a function of the carrier, either hole ( a = p ) or electron ( a = n ), density ( ρ p or ρ n ) along the armchair (X = A) and zigzag (X = Z) direction at three different temperature values. As expected, the carrier mobility clearly exhibits anisotropic behaviors for both phosphorene and PO due to their anisotropic puckered geometries. For phosphorene, we observed that both hole (p-type) and electron (n-type) mobilities, µ p P,Z and µ n P,Z , are diminutive  www.nature.com/scientificreports/ along the zigzag direction, independent of the type of carrier, whereas those along the armchair direction µ p P,A and µ n P,A are substantially high, as shown in Fig. 3a and b. For instance, we evaluated µ p P,A ≈ 150 cm 2 /Vs and µ n P,A ≈ 60 cm 2 /Vs for their respective carrier density values of 1 × 10 12 cm −2 at room temperature (T = 300 K). These results are in good agreement with previously evaluated values 19,25 . For PO, on the other hand, the hole mobility µ p PO is quite low varying, from 0.25 to 2.6 cm 2 /Vs, as shown in Fig. 3c, due to short carrier relaxation time near the valence band edge compared to that of phosphorene presented in Fig. 2b, whereas the electron mobility exhibits much higher values, for instance, µ n PO,Z ≈ 129 cm 2 /Vs with ρ n ≈ 10 11 cm −2 at T = 300 K, as shown in Fig. 3d originating from a significant reduction in the electron-phonon scattering. Note that, in any cases, the mobility tends to decrease as temperature increases because of the enhancement of electron scattering.
The relaxation time was also fed into Eq. (3) to evaluate the spectral conductivity tensor (E, T) and thus the moments L (α) (µ, T) for various transport coefficients given in Eq. (4). Such moments were employed to calculate the electric conductivity, the Seebeck coefficient, and the electronic thermal conductivity tensors using Eqs. (5) - (7). The components corresponding to a specific transport direction, either the armchair or zigzag direction, were extracted from the transport tensors to assess the thermoelectric figure of merit ZT given in Eq. (1) along the selected transport direction. Figure 4 shows the computed power factors ( S 2 σ ) and the electronic thermal conductivities ( κ e ) of both p-and n-type PO as a function of the carrier density at three different temperatures of 300, 500, and 700 K. Both S 2 σ and κ e are not only highly anisotropic similar to the mobility presented in Fig. 3, but they also depend strongly on the type of carriers. Particularly, the power factors evaluated at T = 300 K have a peak at ρ p ≈ 10 13 cm −2 and ρ n ≈ 10 12 cm −2 for the p-and n-type PO, respectively. These carrier density values can be experimentally feasible since similar values were already achieved in other 2D thermoelectric materials, such as graphene 26   www.nature.com/scientificreports/ T = 300 K is about 3.4 mW/K 2 m , which is more than three times higher than that for the p-type one, as shown in Fig. 4a and b. It is noted that this value is still smaller than values of ∼ 29 mW/K 2 m and ∼ 6.5 mW/K 2 m for p-and n-type phosphorene, respectively, evaluated at T = 300 K , as shown in Fig. S1a and b in Supplementary Information. For the n-doped case, however, since the reduction in the lattice thermal conductivity by ∼ 1/26 times 24 is much more dominant than that in power factor, the n-doped PO may still have high ZT values in spite of such low power factors. In most semiconductors, the electronic contribution to the thermal conductivity κ e is usually much smaller than the lattice contribution κ l , and thus the former may provide a negligible effect on the thermoelectric efficiency or ZT. For example, phosphorene has κ e < 10 W/mK at T = 300 K with a moderate doping concentration range from 1 × 10 12 to 1 × 10 13 cm −2 , which is much smaller than its lattice thermal conductivity 24 , as shown in Fig. S1c and d in Supplementary Information. κ e of PO, however, plays a crucial role in determining its ZT value, due to a significant reduction of κ l in PO 24 . Figure 4c and d shows the evaluated κ e of PO as a function of the carrier density of holes and electrons. Similar to its mobility and power factor, κ e of PO has higher values along the zigzag direction than along the armchair direction, which satisfies the Wiedemann-Franz law. Especially for electron-doped cases, the electronic thermal conductivity κ e of n-type PO reached ∼ 4 W/mK , which is comparable to the κ PO l , at around a carrier density of ∼ 10 13 cm −2 . More interestingly, κ e becomes even greater than κ l at higher temperatures because κ l is inversely proportional to T 24 . It is, thus, inevitable to consider the electronic contribution to the thermal conductivity for an accurate estimation of the ZT value of PO.
With all these computed quantities, we are now ready to evaluate the thermoelectric figures of merit ZT of PO and phosphorene for comparison. Figure 5 summaries the temperature dependence of ZT of p-and 10 11 10 12 10 14 10 13 Carrier density (cm -2 ) 10 11 10 12 10 14 10 13 Carrier density (cm -2 ) 10 11 10 12 10 14 10 13 Carrier density (cm -2 )   -and (b) n-type phosphorene (red solid lines) and PO (blue dashed lines) as a function of temperature. Square and triangle symbols denote ZT computed along the zigzag and armchair transport directions. We used ρ p = 6.4 × 10 12 cm −2 and ρ n = 4.8 × 10 12 cm −2 for phosphorene and ρ p = 1.0 × 10 14 cm −2 and ρ n = 1.2 × 10 12 cm −2 for PO, with which the respective figures of merit was estimated to reach their corresponding maxima at T = 700 K , as shown in Fig. S2  www.nature.com/scientificreports/ n-type phosphorene and PO along the zigzag and armchair transport directions computed with each particular value of hole and electron carrier density, with which their ZT values reached to its maximum at T = 700 K as shown in Fig. S2 in Supplementary Information. The p-type phosphorene exhibits a better thermoelectric performance along the armchair direction than not only along the zigzag direction but also along any directions for the n-type counterpart, due mainly to its higher power factors shown in Fig. S1 in Supplementary Information. Its ZT value reaches ∼0.4 at 700 K along the armchair direction in spite of its large κ l . On the other hand, PO shows better thermoelectric performance in the electron-doped case than in the hole-doped case. Particularly, the ZT ≈ 0.5 evaluated at 700 K along the zigzag direction is significantly larger than ZT ≈ 0.03 of phosphorene. It interestingly turned out that the thermoelectric efficiency in the electron doped PO becomes almost isotropic with the ZT value of ∼ 0.45 along the armchair direction. It is suggested that the oxidation in phosphorene diminishes thermoelectric anisotropy observed in pristine phosphorene. This diminution in anisotropy is explained by the opposite directional trends in the power factor and lattice thermal conductivity 24 . Furthermore, the optimal electron doping concentration ( ρ n = 1.2 × 10 12 cm −2 ) for PO is much smaller than that ( ρ p = 6.4 × 10 12 cm −2 , ρ n = 4.8 × 10 12 cm −2 ) for phosphorene. This is a very good advantage because high concentration doping is not only hardly achievable in reality but also in heavily doped semiconductors; there would be significant impurity scattering effects, which are difficult to be treated directly in theoretical calculations.

Summary and conclusions
We presented the thermoelectric properties of a new 2D material, phosphorene oxide, investigated by performing first-principles calculations. For a precise evaluation of its thermoelectric efficiency, we accurately calculated phonon-mediated electrical conductivity and lattice thermal conductivity with anharmonic phonon effects. We found that the highly anisotropic structure of the PO results in anisotropic electrical and thermal transport properties. Nevertheless the thermoelectric efficiency of the n-doped PO is almost isotropic or independent of the transport directions. Our results suggest that the n-type PO is superior in the thermoelectric behaviors to the p-type phosphorene due mainly to a significant reduction in lattice thermal conductivity albeit a small decrease in power factor. It is worth noting that spontaneous oxidation of phosphorene improves overall thermoelectric properties as well as structural stability, and thus the n-type PO would be a potential candidate for environmentally-friendly and cost-effective thermoelectric applications.

Methods
To investigate thermoelectric properties of phosphorene and PO, we first calculated the electronic structures by carrying out first-principles calculations based on the density functional theory 29 using Quantum Espresso code 30 . Plane wave basis was used to expand the electronic wavefunctions with a kinetic energy cutoff of 80 Ry. The core and valence electrons were represented by the norm-conserving pseudopotentials 31,32 , and the exchangecorrelation functional was treated within the generalized gradient approximation of Perdew-Burke-Ernzerhof 33 .
To avoid interactions from layers in neighboring cells, we added a vacuum slab of 20 Å parallel to the layer in the unit cell. The Brillouin zone (BZ) was sampled using a Ŵ-centered 10 × 10 × 1 k-points mesh for the integration over the BZ. Then, the phonon dispersion was calculated using the density functional perturbation theory (DFPT) 34 with a 10 × 10 × 1 q-points grid. As a next step, we explicitly evaluated the energy-and temperature-dependent relaxation time τ n,k (E, T) for given n and k defined as 35,36 in atomic units ( = e = 1 ). Here m,ν and dq/� BZ represent respectively the summation over both electronic band and phonon branch indices m and ν , and the normalized q-space integration over the whole BZ with its volume BZ . The Fermi-Dirac and Bose-Einstein distributions, f (E m,k+q , T) and g(ω ν,q , T) represent the electronic and phononic occupations for their corresponding energy eigenvalues E m,k+q and ω ν,q at temperature T. δ(ǫ) is the usual Dirac delta function. G mn,ν (k, q) is the first-order electron-phonon matrix element given by �m, k + q|∂ ν,q V |n, k�/ 2ω ν,q with ∂ ν,q V the derivative of the self-consistent potential associated with a phonon mode ω ν,q 34,37 . It is related to the probability amplitude of phonon-mediated scattering process from a state |n, k� to another |m, k + q� . G mn,ν (k, q) was computed by employing DFPT on 10 × 10 × 1 coarse k-and q-points meshes and then using maximally localized Wannier functions [38][39][40] through the EPW package 35 to interpolate over 100 × 100 × 1 fine k-and q-points grids.
For various transport-related quantities, we solved the semi-classical Boltzmann transport equation 41,42 . The computed electron relaxation times τ n,k given in Eq. (2) together with the group velocity v n,k = ∇ k E n,k calculated from the electronic structure were used to computed the spectral conductivity tensor given as where is the unit cell volume, and the symbol ⊗ denotes the standard outer product. Since phosphorene and PO are 2D systems, the evaluation of the conductivity tensor requires the thicknesses of the systems, which were (2) 1 τ n,k (E, T) = 2π m,ν dq � BZ |G mn,ν (k, q)| 2 × f (E m,k+q , T) + g(ω ν,q , T) δ(E − E m,k+q − ω ν,q ) × 1 − f (E m,k+q , T) + g(ω ν,q , T) δ(E − E m,k+q + ω ν,q ) ,  www.nature.com/scientificreports/