Bimolecular recombination in methylammonium lead triiodide perovskite is an inverse absorption process

Photovoltaic devices based on metal halide perovskites are rapidly improving in efficiency. Once the Shockley–Queisser limit is reached, charge-carrier extraction will be limited only by radiative bimolecular recombination of electrons with holes. Yet, this fundamental process, and its link with material stoichiometry, is still poorly understood. Here we show that bimolecular charge-carrier recombination in methylammonium lead triiodide perovskite can be fully explained as the inverse process of absorption. By correctly accounting for contributions to the absorption from excitons and electron-hole continuum states, we are able to utilise the van Roosbroeck–Shockley relation to determine bimolecular recombination rate constants from absorption spectra. We show that the sharpening of photon, electron and hole distribution functions significantly enhances bimolecular charge recombination as the temperature is lowered, mirroring trends in transient spectroscopy. Our findings provide vital understanding of band-to-band recombination processes in this hybrid perovskite, which comprise direct, fully radiative transitions between thermalized electrons and holes.

Reflection and transmission spectra of CH 3 NH 3 PbI 3 at 4 K.

Supplementary Note 3: Experimentally Determined Absorption Spectra
Absorption coefficient spectra for the full temperature range of CH 3 NH 3 PbI 3 are shown in Supplementary Figure 3. A phase transition occurs in the range 140-170 K between the lowtemperature orthorhombic and the tetragonal structures. At 155 K there is clear evidence of a coexistence of the two phases. 1 The excitonic absorption peak is most pronounced at low temperatures. At higher temperatures the peak becomes smeared out because of thermal broadening. The absorption coefficient (α) was obtained from the reflectance (R s ) and transmittance (T s ) by, where d is the thickness of the sample. Supplementary Equation 1 is valid for samples on non-absorbing substrates such as the quartz substrates used in our study.  Figure 3: Absorption spectra of CH 3 NH 3 PbI 3 for the full temperature range from 4 K to 295 K as a function of energy (eV) on the bottom axis and wavelength (nm) on the top axis. The absorption coefficient was determined using the reflection and transmission spectra and Supplementary Equation 1.

Supplementary Note 4: GW and DFT First Principle Calculations
All calculations are performed for the low-temperature orthorhombic crystal structure of CH 3 NH 3 PbI 3 , as reported by Baikie et al. 2 This choice is motivated by the fixed orientation of the CH 3 NH + 3 cations in the orthorhombic phase.
We calculate the electronic structure of CH 3 NH 3 PbI 3 within the local density approximation to density functional theory (DFT/LDA) 3, 4 including spin-orbit coupling, as implemented in Quantum Espresso. 5 Furthermore, we use the GW approximation 6 as implemented in the Yambo code, 7 and a self-consistent scissor GW scheme 8,9 to calculate the quasiparticle band structures. We use Wannier interpolation as implemented in the Wannier90 code [10][11][12][13] in order to obtain the quasiparticle eigenvalues at arbitrary k-points in the Brillouin zone.
In the following we give the computational details of every calculation step. and others, 16 self-consistency is important in the correct calculation of the quasiparticle band gap and effective masses of lead-halide perovskites. Here we take self-consistency into account by using a self-consistent scissor scheme described in Refs. 8, 9. Notably, in this calculation we are using a denser k-point mesh in the calculation of the quasiparticle self-energy than that reported in Refs. 8, 9, and we obtain a band gap of 1.57 eV. The band gap is 0.14 eV smaller than that obtained in our previous calculations, while the effective masses are the same as we reported in Ref. 9.
Wannier interpolation. The Wannier interpolation follows the same steps as described in Ref. 9. We use the maximally-localized Wannier functions calculated at the DFT/LDA level to interpolate the quasiparticle eigenvalues, as described in Ref. 9. The details of this calculation are reproduced here for convenience. We consider the electronic states around the band gap, using the Pb-p and I-p orbitals as initial guesses, and extracting 72 and 24 maximally-localized Wannier functions for the valence and conduction bands, respectively.
We optimize the Wannier functions for the valence and conduction band simultaneously.
Joint density of states (JDOS). For the calculation of the JDOS we use a 100 × 100 × 100 uniform k-point mesh centered at Γ to discretize the Brilloiun zone. The energy eigenvalues for each of these k-points are obtained from Wannier interpolation. For the calculation of the JDOS we sum over all available transition energies calculated for the 100 × 100 × 100 uniform k-point mesh. Each transition is represented by a Gaussian function centered at the transition energy. The width of the transition is calculated as the imaginary part of the electron-phonon self-energy using the EPW code, 17 as described in Ref. 18.
Independent particle absorption spectrum. We calculate the imaginary part of the independent particle relative dielectric function using the expression: 19, 20 where e is the electronic charge, ε 0 is the dielectric permitivity of vacuum, Ω is the unit cell volume, m e is the electron mass,ê is the polarization vector of the electromagnetic field, p is the electron momentum operator, ψ ck and ψ vk are the Kohn-Sham wave functions for the valence and conduction band states respectively, ck and vk are the corresponding Kohn-Sham eigenvalues, andhω is the photon energy. The integral is performed over the Brillouin zone and Ω BZ is the Brillouin zone volume. The summation is performed over all conduction and valence band states. In our calculation we include only the first conduction band, which is sufficient to capture all optical transitions in the spectrum up to 0.6 eV above the onset.
To calculate the optical absorption using the GW band structure we neglect the offdiagonal terms of the self-energy, as described in Ref. 19, and calculate the dipole matrix elements (shown in Figure 2(c) of the main manuscript) as: 19 where E ck and E vk are the quasiparticle eigenvalues for the conduction and valence band respectively. For comparison, we calculate the GW optical absorption spectrum using either the scaled dipole matrix elements (Supplementary Equation 3) or the dipole matrix elements calculated directly from DFT, using the following two expressions: ε GW not scaled We calculate the dipole matrix elements | ψ ck |ê · p|ψ vk | 2 separately using the Yambo code, and use the Wannier interpolated eigenenergies from DFT and GW , respectively, to obtain the imaginary part of the dielectric function. The optical absorption coefficient is calculated as α(ω) = ω 2 (ω)/(nc), where we considered the refractive index to be a constant, n = 2.4 21 and c is the speed of light in vacuum. In all our calculations the optical matrix elements are averaged over the three polarization directions.
In order to reduce the computational cost, we calculate the optical matrix elements for a 40×40×40 k-point mesh, using norm-conserving LDA pseudopotentials for Pb and I without the semicore d-states, and a plane wave cutoff of 50 Ry. We have checked the difference between optical matrix elements calculated with and without semicore states for 7 k-points along the Γ-X direction, and obtained that in the latter case the optical matrix elements are underestimated by up to 16%. In addition, all Gaussian smearings are rigidly blue-shifted by 15 meV in order to account for inhomogeneous broadening, as discussed in Ref. 18.
Effective masses. We calculate the DFT and GW effective masses as the inverse of the second derivative of the conduction and valence band edges with respect to the wave vector k, as described in Ref. 9. We calculate the second derivative numerically using the finite difference method as in Ref. 9, by calculating the eigenvalues at a distance of 0.01 2π a from Γ in reciprocal space, where a is the smallest of the three lattice parameters. In Supplementary   Table 1 we show the effective masses calculated using DFT/LDA and using the GW band structure and the self-consistent scissor correction, respectively.
Supplementary Table 1: Comparison between the hole, electron and reduced effective masses calculated using DFT/LDA and the using the GW band structure and the self-consistent scissor correction. m e is the rest electron mass.  Proportionality constant µ Reduced effective mass of electron-hole system N (u, σ 2 ) Normal distribution with mean u and variance σ 2 ln(N (u, σ 2 )) Log-normal distribution arg E max () Location of maximum w Vector of local free parameters (E G , E X , σ T ) and global parameters (u s , σ s ). y Measured data Our fitting procedure follows the standard least squares minimisation method: where S is the sum of the squared differences between the experimental data (y i ) and the model for the absorption coefficient (α). Our goal is to minimise S with respect to the free parameters w of the model. In this work we have used Elliott's model for the intensity of optical absorption by excitons. 22 The absorption model is a linear combination of the absorption from bound excitons and electron-hole continuum states and for a direct semiconductor is written as, with, The joint density of states (JDoS) is given by, where, The first term in brackets in Supplementary Equation 8 is the contribution from bound excitonic states and has the form of a line series at energies −E X /n 2 below the band gap (E G ) where n is a positive integer and the magnitude is proportional to 1/n 3 . A pseudo-continuum is produced as n tends to infinity, which is continuous with the absorption of continuum states at the band gap. The absorption associated with continuum states deviates from the It represents the probability of an electron and a hole existing in the same space and is proportional to the overlap of the electron and hole wavefunctions. 23 However, we also need to be in a position to evaluate the absorption coefficient in a system where the Coulombic attraction is not felt between electrons and holes. Physically, this can occur when the attraction is screened owing to a high background charge-carrier density in the system (i.e. above the Mott transition, as discussed further below). Turning off the Coulombic attraction is exactly equivalent to evaluating the limit of the absorption coefficient as the exciton binding energy tends towards zero, since: Taking this limit for Supplementary Equation 8 defines the absorption coefficient α Free (E) as the value expected when Coulomb interactions are fully screened: Hence, we obtain from Supplementary Equation 8 with Therefore, in summary, we may write our expression for the total absorption coefficient α as the sum of bound excitonic (below-gap) absorption and the free (screened) electron-hole absorption multiplied by the Coulombic enhancement factor: Supplementary Figure 4 shows the contributions of bound excitonic states (α X ), Coulombcorrelated continuum states (α C ) and free (screened, α Free ) electrons and holes to the total absorption coefficient spectrum α(E) for a range of different temperatures, obtained from fits of Elliott Theory to the experimentally obtained absorption coefficient spectra. With decreasing temperature, excitonic effects are enhanced, leading to a larger discrepancy between α C and α Free .  Electron-phonon interactions broaden the absorption, with larger effects at higher temperatures. 18 An additional energetic disorder may lead to some element of inhomogeneous broadening. We found that a log-normal distribution was useful to parameterise this additional disorder brought into the system. Previous studies have found a distribution of band gap positions varying with crystallite size and even within a single crystal. 24,25 Our model can be described as the sum of the contributions of excitonic and continuum states convolved with a broadening function. The broadening function is a normalised function centred at the maximum and is the convolution of a normal distribution, caused by electron-phonon coupling, and a log-normal distribution, caused by disorder and local fluctuations of the stoichiometry of the material. In mathematical notation, where ⊗ represents a convolution, u s and σ s are fitted globally and σ T are temperature dependent parameters fitted for each temperature representing the electron-phonon coupling which starts from the principle of detailed balance, the radiative recombination rate (R Rad ) at thermal equilibrium for an infinitesimal energy interval dE at energy E is equal to the generation rate of electron-hole pairs by thermal radiation at that energy. This rate is P (E, T )ρ(E, T ) per unit volume and energy interval, where ρ(E, T )d(E) is the density of photons in the material in the interval dE, and P (E, T ) is the probability per unit time that a photon of energy E is absorbed. The total radiative recombination rate per unit volume may be written as the integral over E, Supplementary Equation 14 of van Roosbroeck and Shockley's article 26 gives P (E)ρ(E) in terms of frequency ν and absorption index κ, Replacing the expression for the absorption index κ with the absorption coefficient α = 4πnκν/c, and using a change of variables hν = E gives, where β is the reciprocal of the thermodynamic temperature. We make the approximation 27 that the refractive index in the medium is constant relative to large variations in the photon distribution ρ with temperature; n r = 2.4. 21 We consider this assumption to be reasonable since the refractive index is directly linked with absorption through the Kramers-Kronig relations. 23 We find that while the magnitude of the absorption coefficient may theoretically vary with temperature, such changes are very minor for CH 3 NH 3 PbI 3 near the band edge (see Supplementary Figure 3) suggesting that changes in refractive index are similarly negligible.
Since e βE 1 near the band gap, the Bose-Einstein distribution term in Supplementary   Equation 21 can be written as e −βE . Using this approximation we obtain the following expression for the total radiative recombination rate per unit volume: Intrinsic charge-carrier concentration. The electron density n 0 in the conduction band is given by: where m * e is the effective mass for the electrons as determined from GW and shown in Supplementary which is equivalent to: Using a change of variables (x = (E − E C )β; dx = βdE): The integral evaluates to √ π/2, and so the density of electrons in the conduction band is: Similarly, the density of holes p 0 can be expressed as, where m * h is the effective mass for the holes as determined from GW and shown in Supplementary Table 1 Bimolecular recombination rate. The intrinsic charge-carrier density may be used to convert the radiative recombination rate to the bimolecular recombination constant k 2 , according to where X and n are the exciton and free charge carrier densities with generation rates G C and G X respectively. k XC , k CX , and k XG are rate constants for transitions between the ground, exciton and continuum states as detailed in Supplementary Figure where X 0 and n i are the equilibrium densities of excitons and free electrons and holes respectively, with n i given by Supplementary Equation 30. The total radiative recombination will be the sum of the excitonic recombination and free electron-hole recombination to the ground state, however, only the latter is a second order process under equilibrium conditions.
Hence, the absorption coefficient that should be used in the van Roosbroeck and Shockley relation in order to determine the bimolecular rate constant (k 2 ) is the absorption coefficient of the continuum states, with or without Coulombic enhancement, and here we discuss these two possibilities: where a B is the Bohr radius of the exciton, given by: with a H the Bohr radius for hydrogen (5.29×10 −11 m) and µ the reduced mass of the electronhole system (here assumed to be 0.11 as determined from GW and shown in Supplementary Table 1).
Evaluation of Equations 39-43 and 44 also requires knowledge of the value of the dielectric function at the frequency corresponding to the exciton binding energy. We are able to determine this value of X r as a function of temperature from the exciton binding energy E X (as extracted from fits based on Elliott's Theory to the absorption spectra) using:   Supplementary Figure 8 shows that in the tetragonal phase (> 160 K) there is one peak resulting from band-to-band recombination. At higher temperatures, the PL peak energy blueshifts monotonically with increasing temperature, which is consistent with the known temperature dependence of the bandgap energy in hybrid lead halide perovskites. 46,47 In the orthorhombic phase (< 160 K) the PL spectra exhibit additional broad peaks at lower photon energies much lower than the band-gap, which has been observed many times for solution-processed CH 3 NH 3 PbI 3 thin films 18,40 and single crystals. 48 The additional emission peaks result from impurity states and small inclusions of the higher-temperature tetragonal phase that do not feature in the absorption because they are present as a minority species. 1,49 Immediately after photoexcitation, these states are therefore not populated and only become so as a result of charge-carrier diffusion 1 over the long life times of charge-carriers at the low fluences employed for these PL measurements. Hence the population of these states is only completed at times much greater than the time range probed by the transient spectroscopy measurements (2 ns) and the PL from these states will not contribute to the emitted radiation during the first 2 ns. Therefore, self-absorption over the first 2 ns after excitation is governed mostly by the PL and absorption from the majority orthorhombic phase. To account for this, we reconstruct the band-to-band recombination peak of the majority orthorhombic phase by mirroring the PL spectrum about the central peak position from the high-energy half onto the low-energy half, as displayed in Supplementary Figure 8. In addition, the PL is corrected for self-absorption effects which results in a small blue-shift in comparison with the uncorrected PL.
Using our model for PL re-absorption, we link the effective bimolecular recombination rate constants determined from transient spectroscopy 40  as a function of temperature. This correction factor allows for transfer of bimolecular recombination rate constants derived from a simple rate equation model to values that also take into account re-absorption of photons and charge-carrier diffusion.