Density functional theory + U modeling of polarons in organohalide lead perovskites

We investigate the possible formation of polarons in four organic perovskites (CH3NH3PbI3, CH3NH3PbBr3, CH3NH3PbCl3, and CH3NH3PbI2Cl1) using a density functional theory (DFT) calculations with local potentials and hybrid functionals. We show that DFT+U method with U = 8 eV predicts a correct band-gap and matches the forces on ions from hybrid calculations. We then use the DFT + U approach to study the effect of polarons, i.e. to search the configuration space and locate the lowest energy localized band gap state self-trapped hole (STH). STH configurations were found for three pure halides and one mixed halide system. Spin orbit coupling (SOC) was also taken into account and the results may be found in the supplementary material. This study focuses on the +U method; however, SOC corrections added to the DFT+U calculations also resulted in STH states in all four systems. © 2016 Author(s). All article content, except where otherwise noted, is licensed under a Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/). [http://dx.doi.org/10.1063/1.4972341]

performed using the Vienna Ab initio Simulation package (VASP), together with the projector augmented wave method (VASP-PAW) [25][26][27][28] and GGA-PBE for exchange and correlation. 29First, we present the results of calculations for the perovskite CH 3 NH 3 PbI 3 , which is a benchmark absorber for solar cell devices.A pseudocubic unit cell taken from F. Brivio et al 30 was relaxed using the conjugate gradient method until all forces fell below 2 meV/Å.Lead 5d orbitals were treated as valence throughout.During self-consistent calculations, an explicit 256-point k-mesh (determined from a single 8x8x8 automatic Monkhorst-Pack 31 k-mesh calculation) was used with a 400 eV kinetic-energy cutoff.Gaussian smearing was implemented with a width of 1 meV.During non-self-consistent calculations, a linear k point scheme was used to extensively sample the Brillioun zone in all Cartesian directions away from the band gap (R point).To study the polaronic system, one valence electron was removed during relaxation and the LDA+U method was used with U = 8eV applied to the halide atom (iodine and chlorine in the I 2 Cl molecule).Dudarev's rotationally invariant LDA+U method was employed in the p orbitals of the halide atom.Simulations were carried out using 2×2×2 supercells.Dielectric constants were calculated using a total energy/wavefunction cutoff of 10 -9 eV and the finite difference method to calculate the derivatives of the cell-periodic part of the orbitals.Gaussian smearing was used with a 0.1 eV width in the optical properties calculations.Ten thousand grid points were used for the static dielectric constant for density of states (DOS) evaluation to improve the resolution of the DOS plots.For the high frequency constants, the same total energy/wavefunction cutoff of 10 -9 was used, as was Gaussian smearing with 0.1eV width and ten thousand DOS grid points.The number of bands used in the high frequency calculation was increased from the number of bands in the static frequency calculation from 32 to 120.Effective masses were calculated using parabolic band fitting in each of the seven symmetry directions away from the R point.
To determine the appropriate Hubbard correction parameter, comparisons of multiple identical DFT+U calculations were done varying only U, both for systems with and without SOC.Fig. 1 shows optimization of the U parameter in terms of band gap energy, the calculated band diagram and the total density of states (TDOS) with projected density of states (PDOS) for CH 3 NH 3 PbI 3 using the local density approximation (LDA+U) with the Hubbard-like term of U = 8eV; this method of optimizing U was also shown to work by E. Kioupakis et al. 32 We obtained this parameter by ionic force matching between LDA+U and hybrid exchange-correlation functionals (HSE06) and checked that both methods result in good piecewise linearity.Figs.1(a) and 1(b) shows that U = 8 eV with full valence occupancy yields the band-gap energy value, which matches well with the experimental value of 1.6 ± 0.1 eV; 33,34 Fig. 1(c) shows that with one electron removed from valence, the LDA+U method realizes a polaron as a bound state in the band gap with an upward shift of the conduction band minimum (CBM) from 1.62 eV to 1.90 eV.The polaron is seen in Fig. 1(c) as an asymmetry between spin-up and spin-down electrons in the TDOS (black line) and in the PDOS (color lines) as mostly halide contribution with a small contribution from lead.Thus, the hole polaron resides between the lead and halide atoms.
As is sometimes the case, a polaron may be very large, spreading across many unit cells.Therefore, the exact location of the polaron is unknown in the iodine system due to computational limitations of this study.However, distortions are seen as ionic shifts in the halide and Pb atoms.In the iodine system, the four I atoms at the [010] plane shifted outward 1Å along 110 , the four I atoms at the [100] plane shifted outward 1Å along 011 and the four Pb atoms at the [001] plane shifted 1Å along 101 .This indicates that the large polaron may exist in the 111 direction in a supercell that is larger than 2x2x2.
The LDA+U method is also successful in realizing polarons in the other pure halide systems using U = 8 eV.In each system, the CBM is shifted upward and the polaron is seen as a partially bound state in the band gap.The full valence occupancy band gaps and polaron energies in relation to valence band maxima (VBM) for each halide are summarized in Table I, along with dielectric constants and effective masses for both full occupancy and less one electron systems; SOC results for the values are include in the supplementary material.One notices that as the size of the halide atom increases, the band gap energy decreases with the pure iodine system having the smallest band gap as well as the smallest polaron energy offset from the valence band maximum (VBM); the mixed halide has a band gap energy of 2.0 eV, close to that of the pure bromine system, which is roughly the average size between chlorine and iodine.The pure bromine system band gap compares favorably with the experimental results of R.A. Jishi et al 35 and the pure chlorine with the experimental results of N. Kitazawa et al. 36 Both the static and high frequency dielectric constants increase with increasing atomic size of the halide atom.The high frequency results are similar to those of F. Brivio et al 37 but as the static frequencies are very sensitive to the theory used, large variations may be seen between competing theoretical methods.For all three pure halide systems, when going from a fully occupied system to one with one less electron, the effective mass of free electrons with respect to the rest mass TABLE I. Calculated material properties of perovskite structures.Effective masses are in units of m 0 , the rest mass of an electron and the static dielectric constants ε 0 are the average of the trace of the static dielectric tensor.The effective masses are each the numerical average of the masses for all seven symmetry directions. of the electron decreases while the effective mass of free holes with respect to the rest mass of the electron increases (except in I 3 ); both effective masses increase for the mixed halide.The effective mass of free electrons decreasing in the pure halides relates to an increase in free electron mobility.

Full
As the use of parabolic band fitting is sensitive to the distance from each k point as well as the number of k points used, results may vary; 0.5% deviations away from the band gap and 15 total k points were used here and results were comparable to those of D. Li et al 38 for both pure iodine as well as the mixed halide system.Fig. 2 shows that a polaron is also realized in the mixed halide system CH 3 NH 3 PbI 2 Cl.The U = 8 eV correction was added to the p orbitals of I and Cl which results in a localized band gap state as well as an asymmetry in the DOS with one electron removed from valence.
We conclude that the DFT+U method is appropriate for realizing polarons in halide perovskite materials.The Hubbard correction term may be determined by matching forces between GGA+U and hybrid functionals and looking to experimental values of the band gap and comparing this to the results of identical self-consistent relaxations over a range of U values.Polarons are then realized as shifts in the CBM of the band structure, localized states in the band gap and asymmetry in the TDOS in the band gap.With the DFT+U method realizing polarons in the halides, one may argue that polaron repulsion is a strong candidate for explaining the mobility phenomena seen experimentally.

FIG. 1 .
FIG. 1.(a) Band gap energy vs U parameter for LDA+U calculations in MAPbI 3 showing correct band gap energy of E g = 1.6 eV for U = 8 eV; full valence was used to optimize the U parameter.(b) Band structure of full valence LDA+U calculation for MAPbI 3 with U = 8eV to show correct band gap energy.The Fermi level is set at the energy zero.(c) Partial density of states (PDOS) of MAPbI 3 using LDA+U method with one electron removed from valence; polaron is seen as an asymmetry in the TDOS and as a shifting upward of the conduction band minimum from 1.6 eV to 1.9 eV.Black lines are the total density of states showing spin up and down electrons, while color lines relate to the projected density of states of each atomic species.

FIG. 2 .
FIG. 2. (a)Band structure and (b) total and projected density of states of MAPbI 2 Cl mixed halide system using LDA+U method with one electron removed.No polaron is realized in the mixed halide system.The Fermi level is set at the energy zero.