Supporting Information Zero-Dimensional Hybrid Organic–Inorganic Halide Perovskite Modeling: Insights from First Principles Giacomo Giorgi1,2,* and Koichi Yamashita1,3,* 1
Department of Chemical System Engineering, School of Engineering, The University of Tokyo, 7-3-1, Hongo, Bunkyo-ku, Tokyo, Japan
2
Dipartimento di Ingegneria Civile e Ambientale, Università degli Studi di Perugia, Via G. Duranti, 06125 Perugia, Italy. 3
CREST-JST, 7 Gobancho, Chiyoda-ku, Tokyo 102-0076, Japan.
By means of a combination of ab initio molecular dynamics (AIMD) and scalar relativistic density-functional-based calculations, we theoretically characterized ‘bulk cut’ nanoclusters of increasing sizes, focusing on their structural and electronic properties. We have performed several AIMD simulations using the pseudoatomic orbitals (PAO) based SIESTA code.1 The generalized gradient approximation for the exchange and correlation functional of Perdew, Burke, and Ernzerhof (GGA/PBE) was used2 along with the norm-conserving pseudopotentials of the Troullier–Martins type3 for the description of the core electrons. For the Pb atoms, we used a scalar relativistic pseudopotential with 14 valence electrons,4 ten as semicore states (5d10) associated with a single-ζ plus polarization basis set, and four in the valence state (6s26p2) with a double-ζ plus polarization basis set. For the remaining species, nonrelativistic pseudopotentials with a standard double-ζ plus polarization basis set were employed. To test the validity of our theoretical setup, we initially fully
optimized both the tetragonal and the cubic polymorphs of MAPbI3, starting from previously reported experimental5,
6
and theoretical7,
8
structures, and relaxing the
positions of all the atoms until the maximum force was < 0.02 eV/Å. In this case, a Γ centered k-point sampling of 8 × 8 × 6 (8 × 8 × 8) [12 × 12 × 10 (14 × 14 × 14) for the electronic property calculations] was employed in the optimization of the tetragonal (cubic) polymorph. The electronic properties were in very good agreement with the experimentally and theoretically determined values. As widely documented in the literature,7-10 the agreement between Density Functional Theory (DFT) formalism and experimental ones is merely fortuitous because the inclusion of fully relativistic effects (SOC) has a great effect on the band edges;10 the calculated bandgap including SOC is indeed reduced by > 50% compared with the experimentally determined value. We are furthermore aware that to be extremely rigorous the true experimental bandgap would be recovered by adding many-body corrections to the fully relativistic calculated systems (SOC+GW): the very large computational cost associated with such setup excludes its use in our systems. Other recent studies11 have adopted a computational setup that takes into account the on-site Coulomb interaction along with the inclusion of relativistic effects (SOC-U) in order to calculate the electronic properties of Cd-doped MAPbI3. Concerning the present study the agreement between ‘bare’ DFT and experimental data makes our computational choice adequate for achieving our purpose. The main optimized geometrical parameters for tetragonal and cubic cells are reported in Table 1s. Figure 1s shows the band structure, DOS, and band edge wave functions of the tetragonal polymorph.
Table 1s. Main geometrical and electronic properties of the tetragonal and cubic PAO/PBE optimized structures. (Lattice parameters and bond length in Å, volume in Å3, and bandgap in eV.) Tetragonal
Cubic
a=8.77; b=8.77; c=12.93 a
[exp (8.85; 8.85; 12.66)
(8.80; 8.80; 12.68)b
a=6.35; b=6.36; c=6.35 [exp: (6.33;6.33; 6.33)a [theo: (6.45; 6.48; 6.45)g]
(8.849; 8.849; 12.642)c] [theo (8.86; 8.86; 12.66)d
V
EG
a
995.05
256.28
(exp: 992.6a; 990.0c)
(exp: 253.5a)
1.57 (on Γ)
1.48 (on R)
1.66d
1.55f
1.57e
1.64g
Ref. 12; bRef. 5; cRef. 6; d Ref. 7; eRef. 13; fRef. 14; gRef 8.
Figure 1s. (a) PAO/PBE calculated band structure and (b) DOS/PDOS for tMAPbI3. (c) Lateral and top view of Valence Band Maximum (VBM), and (d) lateral and top view of the Conduction Band Minimum (CBM)of t-MAPbI3.(Gray, Pb; purple, I; brown, C; light blue, N; white, H atoms. Isosurface level 0.08 eV/Å.)
Following such encouraging agreement, we considered a 3 × 2 × 2 supercell of the tetragonal polymorph (576 atoms) and performed an AIMD of Parrinello–Rahman [15] with variable shape and volume of the unit cell (NPT ensemble). Following previous works,16 we first performed an initial run of 7.5 ps with a time step of 1 fs at T = 300 K without imposing any constraints on the atomic positions.17, 18 Different values of the mass in the Nose Parrinello–Rahman mass barostat (NPT) have been tested: a reduced (10 Ry*fs**2) and a larger (150 Ry*fs**2) mass. To assemble our clusters, we extracted trajectories after about 2.5 ps for the former calculation. Figure 2s shows the time evolution of volume and temperature for the full AIMD simulation of the bulk supercell in the case of the larger mass. From the extracted trajectories, we obtained three clusters of increasing size, characterized by different terminations: both MAI-rich and PbI2-rich. We therefore adjusted their lateral terminations ‘by hand’ and performed a minimal number of AIMD steps (in this case NVT) to refine the atomic coordination of the clusters so obtained. In all the clusters investigated, additional MA molecular cations were randomly included in the empty semiconductor cavities19 to ensure charge neutrality of the final systems. On the OIHP 0D models we performed short AIMD simulations (∼3 ps for the first three clusters and ∼2.3 ps for the last one) with a time step of 0.25 fs in a microcanonical ensemble (NVE), extracting trajectories every 0.1 ps in the range 0.6– 1.5 ps, and every 0.05 ps for the remaining time of the simulation. On the extracted trajectories we performed spin-polarized single-point DFT calculations, also by means of the GGA/PBE exchange-correlation functional, along with a mesh cutoff of 400 Ry. Figure 3s reports both the Energy/simulation time profile for the whole AIMD run for the 1s cluster and the correlation between energy gap variation and reduced effective mass for the same cluster. Figure 4s shows the relationship between averaged volume of the three clusters here considered and their averaged bandgap.
Figure 2s. Time evolution of volume and temperature for the full AIMD simulation on the 3 × 2 × 2 supercell of MAPbI3 (NPT).
Figure 3s. (a) Energy vs time for the whole NVE simulation for the 1s cluster, and (b) correlation between bandgap variation (∆Eg, eV) and reduced effective mass (µ).
Figure 4s. (a) Relationship between averaged volume (D, Å3) of the three clusters (1s, 2i, and 3l) and the averaged bandgap (Eg, eV).
For the sake of comparison, Figure 5s shows the localization of the wave function at the band edges for a cluster of same size of 1s of CsPbI3 (after ~0.3 ps).
Figure 5s. Top view of (a) VBM and (b) CBM of the Cs54Pb27I108 cluster. (Gray, Pb; purple, I; red, Cs. Isosurface level 0.04 eV/Å.)
Figure 6s. (a) Lateral and (b) top view of the VBM wave function of the trajectory of the 2i cluster extracted after 2.65 ps, and (c) and (d) the same views for the CBM. Same ((e) lateral-(f) top view for the VBM, (g) lateral-(h) top view for the CBM wave function) for the trajectory of the 2i cluster extracted after 0.90 ps (Gray, Pb; purple, I; brown, C; light blue, N; white, H atoms. Isosurface level 0.04 eV/Å.)
To calculate the reduced effective mass according to Brus we have used the following equation: 20, 21 ∆Eg ≈ (2h2π2/µD2) ‒ (3.6e2/εD),
(1s)
where ∆Eg represents the difference between the gap of the cluster and that of the bulk, D is the average diameter of the cluster [(6V/π)1/3]. ε is the dielectric constant of the bulk whose choice clearly influences the final value of the reduced effective mass. First, we have considered the experimentally reported dielectric constant of Hirasawa, 22
which is actually only a few tenths of an eV higher than the static dielectric
constant (∼6), calculated by means of density functional perturbation theory.23 Adoption of the dielectric constant calculated by Brivio et al. (∼25) clearly leads to an increase in µ of ~ 30%.24 Volumes of the cavities are calculated with Gaussian0925 at the b3lyp/genECP level by means of the polarizable continuum model (PCM). In particular, 6-311G++ basis set was used for C, N, and H atoms, MWB46 and MWB78 pseudopotentials for I and Pb atoms respectively. To estimate the angular velocity in our systems we have compared the variation of the z-coordinate of the same C atoms between two selected trajectories in the AIMD simulation using the following equation: arcsin (zCI–zCII) = dφ.
(1s)
Calculating the angle (dφ) swept by the same C atom between two trajectories along the z-coordinate (zCI and zCII) allows indeed to obtain the angular velocity (ω) as the ratio dφ/dt, where dt is the time period between the two considered trajectories. To calculate the angular velocity of MA cations in the clusters, we considered three trajectories: the first extracted after ∼2 ps, followed by the other two, one extracted after 2.0625 ps (~250 trajectories after the first at 2 ps) and the second after 2.125 ps (~500 trajectories after the first at 2 ps).
References [1] Soler, J. M.; Artacho, E.; Gale, J. D.; García, A.; Junquera, J.; Ordejón, P.; Sánchez-Portal, D. The SIESTA method for ab initio order-N materials simulation. J. Phys.: Condens. Matter 2002, 14, 2745-2779. [2] Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865-3868. [3] Troullier, N.; Martins, J.L. Efficient pseudopotentials for plane-wave calculations. Phys. Rev. B 1991, 43, 1993-2006. [4]http://departments.icmab.es/dmmis/leem/siesta/Databases/BasisSets/PbTiO3/GG A-Wu-Cohen/. [5] Kawamura, Y.; Mashiyama, H.; Hasebe, K. Structural Study on Cubic– Tetragonal Transition of CH3NH3PbI3. J. Phys. Soc. Jpn. 2002, 71, 1694– 1697. [6] Stoumpos, C. C.; Malliakas, C. D.; Kanatzidis, M. G. Semiconducting Tin and Lead Iodide Perovskites with Organic Cations: Phase Transitions, High Mobilities, and Near-Infrared Photoluminescent Properties. Inorg. Chem. 2013, 52, 9019– 9038. [7] Mosconi, E.; Amat, A.; Nazeeruddin, Md. K.; Grätzel, M.; De Angelis, F. FirstPrinciples Modeling of Mixed Halide Organometal Perovskites for Photovoltaic Applications. J. Phys. Chem. C 2013, 117, 13902–13913. [8] Giorgi, G.; Fujisawa, J.-I.; Segawa, H.; Yamashita, K. Cation Role in Structural and Electronic Properties of 3D Organic–Inorganic Halide Perovskites: A DFT Analysis. J. Phys. Chem. C 2014, 118, 12176–12183. [9] Giorgi, G.; Fujisawa, J.-I.; Segawa, H.; Yamashita, K. Small Photocarrier Effective Masses Featuring Ambipolar Transport in Methylammonium Lead Iodide Perovskite: A Density Functional Analysis. J. Phys. Chem. Lett. 2013, 4, 4213–4216.
[10] Even, J.; Pedesseau, L.; Jancu, J.-M.; Katan, C. Importance of spin–orbit coupling in hybrid organic/inorganic perovskites for photovoltaic applications. J. Phys. Chem. Lett. 2013, 4, 2999– 3005. [11] Navas, J.; Sánchez-Coronilla, A.; Gallardo, J. J.; Martín, E. I.; Hernández, N. C.; Alcántara, R.; Fernández-Lorenzo, C.; Martín-Calleja, J. Revealing the role of Pb2+ in the stability of organic–inorganic hybrid perovskite CH3NH3Pb1−xCdxI3: an experimental and theoretical study Phys. Chem. Chem. Phys. 2015, 17, 23886-23896 [12]
Poglitsch,
A.;
Weber;
methylammoniumtrihalogenoplumbates
(II)
D.
Dynamic
observed
by
disorder
in
millimeter ‐ wave
spectroscopy. J. Chem. Phys. 1987, 87, 6373-6378. [13] Qiu, J.; Qiu, Y.; Yan, K.; Zhong, M.; Mu, C.; Yan, H.; Yang, S. All-solid-state hybrid solar cells based on a new organometal halide perovskite sensitizer and onedimensional TiO2 nanowire arrays. Nanoscale 2013, 5, 3245– 3248. [14] Kojima, A.; Teshima, K.; Shirai, Y.; Miyasaka, T. Organometal Halide Perovskites as Visible-Light Sensitizers for Photovoltaic Cells. J. Am. Chem. Soc. 2009, 131, 6050– 6051. [15] Parrinello, M.; Rahman, A. Crystal Structure and Pair Potentials: A MolecularDynamics Study. Phys. Rev. Lett. 1980, 45, 1196-1199. [16] Mosconi, E.; Quarti, C.; Ivanovska, T.; Ruani, G.; De Angelis, F. Structural and electronic properties of organo-halide lead perovskites: a combined ir-spectroscopy and ab initio molecular dynamics investigation. Phys. Chem. Chem. Phys. 2014, 16, 16137– 1614. [17] Frost, J. M.; Butler, K. T.; Walsh, A. Molecular ferroelectric contributions to anomalous hysteresis in hybrid perovskite solar cells. APL Mater. 2014, 2, 081506-10. [18] Carignano, M. A.; Kachmar, A.; Hutter, J. Thermal Effects on CH3NH3PbI3 Perovskite from Ab Initio Molecular Dynamics Simulations. J. Phys. Chem. C 2015, 119, 8991–8997.
[19] Papavassiliou, G. C.Three- and Low-Dimensional Inorganic Semiconductors. Prog. Solid State Chem. 1997, 25, 125– 270. [20] Li, Y.-F.; Liu, Z.-P. Particle Size, Shape and Activity for Photocatalysis on Titania Anatase Nanoparticles in Aqueous Surroundings. J. Am. Chem. Soc. 2011, 133, 15743–15752. [21] (a) Brus, L.E. Electron–electron and electron ‐ hole interactions in small semiconductor crystallites: The size dependence of the lowest excited electronic state. J. Chem. Phys. 1984, 80, 4403-4409. (b) Brus, L. Electronic wave functions in semiconductor clusters: experiment and theory. J. Phys. Chem. 1986, 90, 2555–2560. [22] Hirasawa, M.; Ishihara, T.; Goto, T.; Uchida, K.; Miura, N. Magnetoabsorption of the lowest exciton in perovskite-type compound (CH3NH3)PbI3, Phys. B: Condens. Matter 1994, 201, 427– 430. [23] Zhou, Y.; Huang, F.; Cheng, Y.-B.; Gray-Weale, A. Photovoltaic performance and the energy landscape of CH3NH3PbI3, Phys. Chem. Chem. Phys. 2015, 17, 2260422615. [24] Brivio, F.; Walker, A. B.; Walsh, A. Structural and electronic properties of hybrid perovskites for high-efficiency thin-film photovoltaics from first-principles. APL Mater. 2013, 1, 042111-5. [25] Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Gaussian 09, Revision A.02; Gaussian, Inc.: Wallingford, CT, 2009.