The Effects of the Organic–Inorganic Interactions on the Thermal

Letter pubs.acs.org/NanoLett

The Effects of the Organic−Inorganic Interactions on the Thermal Transport Properties of CH3NH3PbI3 Tomoyuki Hata,*,†,‡ Giacomo Giorgi,*,†,§ and Koichi Yamashita*,†,‡ †

Department of Chemical System Engineering, School of Engineering, The University of Tokyo, 7-3-1, Hongo, Bunkyo-ku, Tokyo 113-8656, Japan ‡ CREST-JST, 7 Gobancho, Chiyoda-ku, Tokyo 102-0076, Japan § Dipartimento di Ingegneria Civile e Ambientale, Università degli Studi di Perugia, 06125 Perugia, Italy S Supporting Information *

ABSTRACT: Methylammonium lead iodide perovskite (CH3NH3PbI3), the most investigated hybrid organic−inorganic halide perovskite, is characterized by a quite low thermal conductivity. The rotational motion of methylammonium cations is considered responsible for phonon transport suppression; however, to date, the specific mechanism of the process has not been clarified. In this study, we elucidate the role of rotations in thermal properties based on molecular dynamics simulations. To do it, we developed an empirical potential for CH3NH3PbI3 by fitting to ab initio calculations and evaluated its thermal conductivity by means of nonequilibrium molecular dynamics. Results are compared with model systems that include different embedded cations, and this comparison shows a dominant suppression effect provided by rotational motions. We also checked the temperature dependence of the vibrational density of states and specified the energy range in which anharmonic couplings occur. By means of phonon dispersion analysis, we were able to fully elucidate the suppression mechanism: the rotations are coupled with translational motions of cations, via which inorganic lattice vibrations are coupled and scatter each other. KEYWORDS: Thermal conductivity, hybrid organic−inorganic halide perovskites, thermoelectrics, thermoelectric materials, molecular dynamics

I

parameters than fully inorganic ones and embedded cations move incoherently.13 These conditions are similar to the existing phonon−glass electron−crystal materials, namely host−guest structures like type-I clathrate and skutterudite compounds.14−16 In those systems, the rattling motions of the embedded atoms are coupled with the transverse acoustic mode and provide the phonon scattering.14 For MAPI, it is also suggested that the ionic interactions couple inorganic lattice (Pb−I) vibrations with methylammonium (MA) motions. The IR and Raman spectra show that the Pb−I vibration density of states (vDOS) mainly overlaps with the rotation of MA,17−19 so rotational motion is thought to play a critical role in the suppression.11 Even if this is the correct explanation, the specific mechanism remains unclear. In particular, the coupling mechanism between rotational motions and Pb−I vibrations, how it modifies the vibrational states, and the extent to which it suppresses the thermal conductivity are all issues that should be fully clarified. Such further understanding is mandatory not only to control the thermal transport properties of MAPI but

n addition to their well-known superior performance in photovoltaics,1−4 hybrid organic−inorganic halides of the formula ABX3 (A = organic molecular cation; B = Ge, Sn, Pb; X = halide) have recently attracted attention for possible thermoelectric applications.5 The thermoelectric figure of merit (ZT) is expressed as ZT = S2σT/κ, where S is the Seebeck coefficient, T is the temperature, σ is the electronic conductivity, and κ is the thermal conductivity.6 High carrier mobility, high Seebeck coefficient, and low thermal conductivity are desirable for efficient thermoelectric materials. Continuing efforts to understand the light conversion efficiency of ABX3 revealed that their high carrier mobilities are based on small effective masses and long diffusion lengths.7−10 Moreover, recent research shows that CH3NH3PbI3 (hereafter MAPI) has a high Seebeck coefficient, deriving from multidegenerated conduction and valence bands.5 Fortunately, MAPI is also characterized by a beneficial low thermal conductivity.11 Both single-crystal and polycrystalline MAPI are characterized by a thermal conductivity of ∼0.4 Wm−1 K−1 at room temperature, far below that of inorganic perovskites.12 On the one hand, these results show that MAPI meets all the optimal requirements of materials for thermoelectric conversion; on the other hand, the cause of the drastic suppression of phonon transport remains obscure. The hybrid organic−inorganic halide perovskites have larger lattice © XXXX American Chemical Society

Received: February 2, 2016 Revised: March 17, 2016

A

DOI: 10.1021/acs.nanolett.6b00457 Nano Lett. XXXX, XXX, XXX−XXX

Letter

Nano Letters

temperature is maintained at 280 K by the Nosé−Hoover thermostat for 10000 steps. In the fitting procedure, we use the partial least-squares regression and confirm it by 10-fold crossvalidation.31 With the developed potential, we evaluate the thermal conductivity of MAPI. Among all the possible approaches, NEMD methods are appropriate for systems with low κ.23 In this study, we employed the RNEMD method in the LAMMPS package, which attains steady-state heat flux and temperature gradient by exchanging velocities of atoms. After thermalization with the NVT thermostat, RNEMD calculations are performed for 107 steps with a 0.5 fs time step. In RNEMD simulations, we obtained the size-dependent thermal conductivity through which the bulk thermal conductivity is evaluated. For this estimation, we adopt the linear extrapolation method,32 where the bulk thermal conductivities are expressed as 1 1 1 = +a κ (L ) κ (∞ ) L (1)

also to provide guidelines for the design of hybrid thermoelectric materials. In this study, we analyze the thermal transport properties of MAPI. The suppression effects induced by the MA rotation are evaluated and explained through the vibrational coupling analysis.20,21 To perform thermal conductivity calculations including anharmonic effects, we have developed an empirical potential for MAPI based on ab initio molecular dynamics (MD) calculations. By using this force field, we calculate the thermal conductivity of MAPI using the reverse nonequilibrium molecular dynamics (RNEMD) method.22 The results are compared with model systems that include different embedded cations, and it shows that the suppression of heat conduction is mainly attributed to the MA rotations. By checking the phonon band structures, we find that coupling between translational and rotational motions of MA is the key to the suppression mechanism in MAPI. The molecular dynamics method is the best way to evaluate the thermal transport properties of materials with strong anharmonicity. Anharmonic potentials for thermal transport calculations can be developed via several approaches in which the interatomic force constants (IFCs) are fitted to force trajectories of ab initio MD by linear regression.23,24 Such procedures are quite accurate for crystals in which all atoms have rigid stable points for IFCs. However, in host−guest structures embedded components feel complex and shallow potentials, which are a combination of nonbonded interactions with host atoms. Because the method using IFCs is not adequate for such structures, we have improved it by using interatomic potential functions. In these functions, potential energies are attributed to structural components (bonds, angles, and dihedrals) and nonbonded interactions. We have expanded these structural component functions to higher orders and introduced cross terms that include anharmonicity (see Supporting Information for details). The nonbonded interaction terms include van der Waals (vdW) terms based on Lennard-Jones potentials and Coulomb interactions, which are nonlinear functions. To avoid a nonlinear fitting, we use the empirical vdW radii25 and adopt atomic charges from the Voronoi charge analysis.26 The force trajectories of ab initio MD are calculated using density functional theory as implemented in the atomic orbital based code SIESTA,27 employing the GGA-PBE exchange-correlation functional,28 along with the Troullier-Martins pseudopotentials29 for the description of the core electrons. Here we consider the tetragonal (I4/mcm, Z = 4) unit cell of MAPI with 48 atoms (Figure 1), using 3 × 3 × 2 Monkhorst−Pack grids30 for kpoint sampling. The time step is set to 0.5 fs and the

where L is the length in the conduction direction and a is the factor of proportionality. The thermal conductivities obtained are shown in Figure 2. We compared our results with previous

Figure 2. Temperature dependence of the thermal conductivity of MAPI evaluated with RNEMD (red diamonds) and experimentally11 (blue solid line). Experimental results of inorganic perovskites of TbMNO3 (orange dashed line) and YMnO3 (green dotted lines) fully inorganic perovskites.12

experimental ones for MAPI11 and with others for fully inorganic perovskite materials,12 observing very good agreement with the experimental thermal conductivity of MAPI single crystals; regardless, these results are quite lower than those of fully inorganic perovskites over a wide temperature range. The structural difference between MAPI and fully inorganic perovskites is represented by the A-site cation nature; thus, it is intuitive to ascribe the thermal conductivity suppression to the organic−inorganic interactions. To check whether motions of the MA cations are coupled with the inorganic host cages, we analyzed the thermal transport properties of the model structures shown in Figure 3. In Model A, the embedded component is a diatomic molecule and the twisting degrees of freedom are truncated. Model B includes a single atom without rotational motions. Model C is the bare inorganic frame. In all

Figure 1. (a) Lateral and (b) top view of the unit cell of MAPI in the tetragonal phase. (Gray, lead; purple, iodine; silver, carbon; blue, nitrogen; white, hydrogen atoms. Black solid line, unit cell). B

DOI: 10.1021/acs.nanolett.6b00457 Nano Lett. XXXX, XXX, XXX−XXX

Letter

Nano Letters 1 P(ω) = 2π

1 N

N

∫ ∑ mi⟨vi(τ)·vi(t + τ)⟩τ e−iωt dt i=1

(2)

where ω is the angular momentum, N is the number of atoms, m is the atomic mass, v is the velocity, τ is the time, and t is the time interval. To compare the vibration energies of organic and inorganic components, we calculate partial vDOS for two groups, Pb−I and MA. MAPI is thermalized at 150−250 K at constant volume before the velocity sampling. Peak attributions are based on the normal modes of MAPI at the Γ point, viz., q = (0 0 0) where q is the phonon wave vector. The results are shown in Figure 4 for two energy ranges, 400−3300 and 0−400

Figure 3. Crystal structures of MAPI and of Models A, B, and C. Model A has a diatomic molecule, Model B has a single atom as embedded cation, while Model C is the bare inorganic frame.

the three systems, the atomic charges are adjusted to give a final neutral system; in Model A and B the total charges of the embedded components are the same as those of MA cations, whereas Model C is neutralized by slightly lowering the charges of I atoms. Such reduction may induce smaller anharmonic effects but the subsequent overestimation of the thermal conductivity is expected to be rather small than the suppression effects caused by the embedded components. We have also assumed as reference the relaxed volume to modify and take into proper account the van der Waals interactions, avoiding the size effects of the embedded components (see Supporting Information for details). The thermal conductivities are evaluated with the finite calculation model but not extrapolated to bulk values; they are reported for each model in Table 1, Table 1. Thermal Conductivity of MAPI and of Models A, B, and C MAPI degrees of freedom thermal conductivity (W/ mK)

twist + trans + rot 0.185

A trans + rot 0.209

B

C

trans

none

0.302

0.326

Figure 4. Partial vDOS of MAPI: MA (blue solid line) and Pb−I (orange dashed line) components. The energy range is divided into two regions: (a) 400−3300 and (b) 0−400 cm−1. Their intensities are normalized considering the highest peak in each energy region and each group. Temperature is increased from 150 K (thickest line) to 250 K (thinnest line).

which shows that the thermal conductivity of Model C is larger than that of MAPI but still low. Heavy atoms and weak bonds in the inorganic lattice thus contribute to the low thermal transport properties of MAPI. The results also show that the main difference is between Model A and Model B, meaning that the rotational motion of MA efficiently suppresses the thermal conductivity of the inorganic lattice. For further understanding, we analyzed the vibrational properties of MAPI using MD simulations. On the basis of the Wiener−Khintchine theorem, the vDOS are evaluated as power spectra P(ω) of velocity autocorrelation functions.33 Velocities are mass weighted to normalize the intensity of spectra, so vDOS is expressed as

cm−1, and their intensities are normalized by the highest peak in each energy region and each group. The vibration energies are consistent with spectroscopic results in previous studies.17−19,34 In the higher energy region, the main contribution to the spectra is associated only with the MA motions, and these peaks are attributed as bond stretching and angle deformation of MA cations. Increasing the temperature, the intensities of the peaks are reduced but their energies remain unaffected, meaning that these vibrations are not ideally harmonic, but C

DOI: 10.1021/acs.nanolett.6b00457 Nano Lett. XXXX, XXX, XXX−XXX

Letter

Nano Letters

partial dispersion. These results are consistent with the fact that when two vibration systems with different mass are coupled, vibrational energies of the heavy system result almost unaffected.21 In the MA partial dispersion, the blue lines indicate the translational motions and the orange lines are the rotations. At the Γ point, these motions are distinct and symmetric; the N and C atoms have parallel vibrational vectors of the same length. When the energy increases along the wave vector and passes over the region A in Figure 5, some of these vibrations appear to be asymmetric (green lines). This means that the translational and rotational motions are mixed and branch off into two asymmetric modes, N atom dominant and C atom dominant motions. This coupling corresponds to the deviation of molecular motions from the crystal symmetry, which may come from long-range interactions or weak ordering between molecular dipoles. In MAPI, the transverse optic (TO) modes of Pb−I are mainly coupled with the rotational motions of MA. There is a small energy difference between such coupled modes and only limited energy shifts occur in MA branches. On the other hand, by comparing MAPI with MA, we can see that the translational branches are strongly shifted to higher energy because of the coupling with LO modes. These coupled modes are also described by the vDOS analysis. Similarly, also in MAPI translational and rotational branches are mixed into asymmetric vibrations in almost the same wave vector region. Moreover, the strong energy shift results in a crossover between MA branches located in region B. This crossing becomes clear in the phonon dispersion, where the dynamical matrix is projected onto the normal modes at the Γ point (Figure 5), revealing that the energy gap between the two phonon branches of MA in region B corresponds to an avoided crossing, where the translational motions switch to the rotations, and vice versa. Such mixing and crossing result in the coupling between LO and TO modes mediated by the interactions between the translational and rotational motions, respectively. The suppression effects of MA rotations and the anharmonic energy shift in the 70−80 cm−1 region suggest that these couplings limit the thermal conductivity of MAPI. That is, rotations do not directly scatter phonons but the quite low thermal conductivity of MAPI is attributed to the interaction between LO and TO modes via the mixed/crossed translational and rotational motions of MA. In summary, we have presented the mechanism of the low thermal transport properties of MAPI. After developing an empirical potential of MAPI based on ab initio MD and analyzing its thermal conductivities with RNEMD methods, we compared the results with those from model systems that include different embedded cations. We found that the suppression of thermal conductivity mainly derives from the rotational motions of MA cations. The vDOS analysis shows that the vibrational energy shifts as the temperature increases, where anharmonic couplings occur. From the analysis of the phonon band structures, coupled translational and rotational motions of MA are found to interact with the Pb−I cages and build couplings between the isolated lattice vibrations, which suppress the phonon transport in MAPI. We expect that this mechanism is applicable to all organic−inorganic hybrid materials with particular emphasis on the compounds including complex organic molecules, where more vibrations are reciprocally coupled. The symmetry of molecules also affects the coupling intensities and suppression, features that deserve further careful investigation.

that at the same time they do not couple with other vibrations via their anharmonic terms. Both Pb−I vibrations and MA motions are located in the lower energy region. The Pb−I spectra have two broad bands: a lower one located at 40 cm−1 with no energy shifts, and a higher band at 70−80 cm−1 that is clearly shifted to higher energy following a temperature increase, which is attributed to the longitudinal optic (LO) modes. MA components are also characterized by two bands; at 0−200 cm−1, we observe the band associated with the translational and rotational motions. The highest peak in this band indicates the translations and shifts to higher energy as the temperature increases. The second band is delocalized at around 300 cm−1, which is attributed to the twisting motions and has no peak shift. The peak shifts with temperature mainly derive from the anharmonic contributions or from the thermal expansion of crystals.35 Such results at constant volume suggest that anharmonic couplings occur between the translational motions of MA and LO modes of the Pb−I lattice, a feature that seems incompatible with the model analysis. The rotational motions suppress the thermal conductivity but the anharmonic couplings are mainly observed in the translations. A way to remedy this inconsistency stems from the phonon dispersion analysis. The phonon dispersions are directly measured from MD trajectories at 200 K, assuming an equipartition of vibrational energies.36 We also impose the cubic crystal symmetry (Pm3̅m) including a single MA cation in the unit cell, which allows us separate the intramolecular vibrations in the real space from the intermolecular interactions in the reciprocal lattice space. The results from Γ to X (q = (0 0 0) to (0 0 0.5)) are shown in Figure 5. The partial phonon

Figure 5. Phonon dispersion relation of MAPI in the cubic phase (center) with partial components of Pb−I and MA (right) and their projection onto the normal modes at the Γ point (left). We classify the phonon branches based on their major contributions to vibrations: Pb−I (black thin line), MA translational (blue thick line), rotational (orange thick line), and mixed motions (green thick line). The mixing and crossing among MA branches are observed in the wave vector regions A and B, respectively.

dispersions for Pb−I and MA are also evaluated from submatrices of the dynamical matrix. Figure 5 also shows how phonon branches are classified based on their major contributions to vibrations. Thin lines correspond to Pb−Idominant vibrations and we observe how the Pb−I branches in the MAPI dispersion are almost the same as those of the Pb−I D

DOI: 10.1021/acs.nanolett.6b00457 Nano Lett. XXXX, XXX, XXX−XXX

Letter

Nano Letters

■

(21) Hata, T.; Kawai, H.; Jono, R.; Yamashita, K. Nanotechnology 2014, 25, 245703. (22) Müller-Plathe, F. J. Chem. Phys. 1997, 106, 6082−6085. (23) Tadano, T.; Gohda, Y.; Tsuneyuki, S. J. Phys.: Condens. Matter 2014, 26, 225402. (24) Zhou, F. Phys. Rev. Lett. 2014, 113, 185501. (25) Batsanov, S. S. Inorg. Mater. 2001, 37, 871−885. (26) Fonseca Guerra, C.; Handgraaf, J.-W.; Baerends, E. J.; Bickelhaupt, F. M. J. Comput. Chem. 2004, 25, 189−210. (27) Soler, J. M.; Artacho, E.; Gale, J. D.; Garcia, A.; Junquera, J.; Ordejon, P.; Sanchez-Portal, D. J. Phys.: Condens. Matter 2002, 14, 2745−2779. (28) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1996, 77, 3865−3868. (29) Troullier, N.; Martins, J. L. Phys. Rev. B: Condens. Matter Mater. Phys. 1991, 43, 1993−2006. (30) Monkhorst, H. J.; Pack, J. D. Phys. Rev. B 1976, 13, 5188−5192. (31) Mevik, B.-H.; Wehrens, R. J. Stat. Soft. 2007, 18, 2. (32) Schelling, P. K.; Phillpot, S. P.; Keblinski, P. Phys. Rev. B: Condens. Matter Mater. Phys. 2002, 65, 144306. (33) Thomas, M.; Brehm, M.; Fligg, R.; Vöhringer, P.; Kirchner, B. Phys. Chem. Chem. Phys. 2013, 15, 6608−6622. (34) Mattoni, A.; Filippetti, A.; Saba, M. I.; Caddeo, C.; Delugas, P. J. Phys. Chem. Lett. 2016, 7, 529−535. (35) Postmus, C.; Ferraro, J. R.; Mitra, S. S. Phys. Rev. 1968, 174, 983. (36) Kong, L. T. Comput. Phys. Commun. 2011, 182, 2201−2207.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.nanolett.6b00457. Potential development setup, conditions of calculations and their validations, potential parameters of MAPI, and details of model structures. (PDF)

■

AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. *E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS T.H. is supported by JSPS Research Fellowships for Young Scientist (14J10561).

■

REFERENCES

(1) Kojima, A.; Teshima, K.; Shirai, Y.; Miyasaka, T. J. Am. Chem. Soc. 2009, 131, 6050−6051. (2) Kim, H.-S.; Lee, C.-R.; Im, J.-H.; Lee, K.-B.; Moehl, T.; Marchioro, A.; Moon, S.-J.; Humphry-Baker, R.; Yum, J.-H.; Moser, J. E.; Grätzel, M.; Parkb, N.-G. Sci. Rep. 2012, 2, 591. (3) Zhou, H.; Chen, Q.; Li, G.; Luo, S.; Song, T.; Duan, H.-S.; Hong, Z.; You, J.; Liu, Y.; Yang, Y. Science 2014, 345, 542−546. (4) Jeon, N. J.; Noh, J. H.; Yang, W. S.; Kim, Y. C.; Ryu, S.; Seo, J.; Seok, S. I. Nature 2015, 517, 476−480. (5) He, Y.; Galli, G. Chem. Mater. 2014, 26, 5394−5400. (6) Rowe, D. M. Thermoelectrics Handbook: Macro to Nano; CRC Press: Boca Raton, FL, 2006; Chapter 1. (7) Giorgi, G.; Fujisawa, J.-I.; Segawa, H.; Yamashita, K. J. Phys. Chem. Lett. 2013, 4, 4213−4216. (8) Xing, G.; Mathews, N.; Sun, S.; Lim, S. S.; Lam, Y. M.; Grätzel, M.; Mhaisalkar, S.; Sum, T. C. Science 2013, 342, 344−347. (9) Brivio, F.; Butler, K. T.; Walsh, A. Phys. Rev. B: Condens. Matter Mater. Phys. 2014, 89, 155204. (10) Kawai, H.; Giorgi, G.; Marini, A.; Yamashita, K. Nano Lett. 2015, 15, 3103−3108. (11) Pisoni, A.; Jaćimović, J.; Barišić, O. S.; Spina, M.; Gaál, R.; Forró, L.; Horváth, E. J. Phys. Chem. Lett. 2014, 5, 2488−2492. (12) Zhou, J.-S.; Goodenough, J. B. Phys. Rev. Lett. 2006, 96, 247202. (13) Leguy, A. M. A.; Frost, J. M.; McMahon, A. P.; Sakai, V. G.; Kockelmann, W.; Law, C. H.; Li, X.; Foglia, F.; Walsh, A.; O’Regan, B. C.; Nelson, J.; Cabral, J. T.; Barnes, P. R. F. Nat. Commun. 2015, 6, 7124. (14) Takabatake, T.; Suekuni, K. Rev. Mod. Phys. 2014, 86, 669. (15) Nolas, G. S.; Cohn, J. L.; Slack, G. A.; Schujman, S. B. Appl. Phys. Lett. 1998, 73, 178. (16) Nolas, G. S.; Morelli, D. T.; Tritt, T. M. Annu. Rev. Mater. Sci. 1999, 29, 89−116. (17) Quarti, C.; Grancini, G.; Mosconi, E.; Bruno, P.; Ball, J. M.; Lee, M. M.; Snaith, H. J.; Petrozza, A.; Angelis, F. D. J. Phys. Chem. Lett. 2014, 5, 279−284. (18) Quarti, C.; Mosconi, E.; De Angelis, F. D. Phys. Chem. Chem. Phys. 2015, 17, 9394−9409. (19) Brivio, F.; Frost, J. M.; Skelton, J. M.; Jackson, A. J.; Weber, O. J.; Weller, M. T.; Goñi, A. R.; Leguy, A. M. A.; Barnes, P. R. F.; Walsh, A. Phys. Rev. B: Condens. Matter Mater. Phys. 2015, 92, 144308. (20) Damnjanović, M.; Dobardžić, E.; Milošević, I.; Vuković, T.; Nikolić, B. New J. Phys. 2003, 5, 148. E

DOI: 10.1021/acs.nanolett.6b00457 Nano Lett. XXXX, XXX, XXX−XXX