Role of hydrogen-bonding and its interplay with octahedral tilting in CH₃NH₃PbI₃

Abstract

First principles calculations on the hybrid perovskite CH₃NH₃PbI₃ predict strong hydrogen-bonding which influences the structure and dynamics of the methylammonium cation and reveal its interaction with the tilting of the PbI₆ octahedra. The calculated atomic coordinates are in excellent agreement with neutron diffraction results.

---

CH₃NH₃(MA)PbI₃ is one of the most extensively studied hybrid halide perovskites due to its impressive power conversion efficiency for solar cell applications.¹ In spite of numerous studies on its structural properties,^2–8 relatively little progress has been made in our understanding on the role of hydrogen-bonding in the orthorhombic phase (o-MAPbI₃). In this study we have used density functional theory (DFT) calculations to examine the thermodynamic consequences of hydrogen-bonding and its interplay with the octahedral rotations that are ubiquitous in perovskites. The computed atomic positions, including those of the hydrogen atoms, are compared to a recent powder neutron diffraction study.⁹

The DFT calculations were performed using the generalized gradient approximation (GGA) and projector augmented wave (PAW)¹⁰ pseudopotentials as implemented in the Vienna ab initio Simulation Package (VASP).^11–13 The effects of spin–orbit coupling and van der Waals (vdW) interactions¹⁴ were included during structural and electronic relaxation. We adopted (i) a 4 × 3 × 4 Monkhorst–Pack k-point mesh centered at Γ,¹⁵ (ii) a 500 eV plane-wave kinetic energy cutoff, and (iii) the tetrahedron method with the Blöchl corrections for the Brillouin zone integrations.¹⁶ We explicitly treated 14 valence electrons for Pb (5d¹⁰6s²6p²), 7 for I (5s²5p⁵), 4 for C (2s²2p²), 5 for N (2s²2p³), and 1 for H (1s¹). The ions were relaxed until the forces on them were less than 0.01 eV Å⁻¹. The experimentally observed ratio of the orthorhombic unit cell vectors² was imposed when we optimized the lattice parameters by calculating the Kohn–Sham (K–S) energy as a function of volume.¹⁷

Fig. 1 shows the optimized crystal structure of o-MAPbI₃ with space group Pnma. The computed lattice parameters are a = 8.844 Å, b = 12.592 Å, and c = 8.563 Å, in good agreement with the experimental values at 100 K (a = 8.8362 Å, b = 12.5804 Å, and c = 8.5551 Å for the X-ray study² and a = 8.8657 Å, b = 12.6293 Å, and c = 8.5769 Å for the neutron study⁹). Note that the calculated values lie between the experimental ones. The perovskite structure exhibits corner-linked PbI₆ octahedra forming a three-dimensional network. The Pnma phase displays a particular pattern of rotations of the octahedra involving either in-phase or out-of-phase rotations around the three I–Pb–I octahedral axes, denoted by a⁻b⁺a⁻ in Glazer's notation.¹⁸ Regarding the organic MA ion, the C–N bonds are aligned approximately along the [101] and [101̄] directions with a head-to-tail configuration as shown in Fig. 1. Our optimized atomic positions for this structure are listed in Table 1, where H atoms that are connected to N atoms are denoted by H_N, and H atoms that are connected to C atoms are denoted by H_C. The subscripts A and E on the iodine atoms refer to axial and equatorial configurations, respectively, as shown in Fig. 1(c). Our optimized atomic positions agree very well with the previous experimental values. In particular, the coordinates of the H atoms are in very good agreement with the recent powder neutron diffraction study.⁹

Fig. 1 (a) Optimized crystal structure of o-MAPbI₃, viewed in the a-b and, (b) a–c planes. PbI₆ octahedra are highlighted in purple. (c) Local cage structure of o-MAPbI₃ is characterized by hydrogen-bonding between the three H atoms in the NH₃ group and one axial I [I_A(1)] atom and two equatorial I [I_E(2) and I_E(2′)] atoms. Hydrogen-bonding interactions are indicated by dotted lines.

Table 1 Calculated fractional coordinates for o-MAPbI₃ compared with experimental values. The H_C(1) and H_C(2) atoms are connected to C while the H_N(1) and H_N(2) atoms are connected to N in the MA cation. ND and XD denote neutron powder diffraction data⁹ and X-ray diffraction data,² respectively

	Wyckoff site	This work				Experiment
		x	y	z		x	y	z
Pb	4b	0.500	0.000	0.000	ND	0.500	0.000	0.000
					XD	0.500	0.000	0.000
IA(1)	4c	0.476	0.250	0.938	ND	0.484	0.250	0.944
					XD	0.486	0.250	0.947
IE(2)	8d	0.176	0.017	0.172	ND	0.189	0.015	0.184
					XD	0.190	0.017	0.186
N	4c	0.953	0.750	0.023	ND	0.942	0.750	0.030
					XD	0.932	0.750	0.029
C	4c	0.908	0.250	0.075	ND	0.937	0.250	0.058
					XD	0.913	0.250	0.061
HC(1)	4c	0.941	0.250	0.198	ND	0.937	0.250	0.187
HC(2)	8d	0.842	0.178	0.048	ND	0.866	0.170	0.029
HN(1)	8d	0.111	0.183	0.003	ND	0.128	0.189	−0.009
HN(2)	4c	0.978	0.750	0.142	ND	0.954	0.750	0.146

---

Fig. 1(c) shows the ground-state configuration of the H atoms in the surrounding inorganic Pb–I cage. The positions of the H atoms conform with the Pnma space group. If we assume that hydrogen-bonding is important when H⋯I < 3 Å, there are three hydrogen-bonds in the Pb–I cage between H_N(1), H_N(2), and H_N(2′) and I_A(1), I_E(2), and I_E(2′), respectively.

These three H⋯I bonds per MA ion are controlled by both the particular inorganic a⁻b⁺a⁻ tilt pattern and the organic MA conformation. Among these, the H_N(1)⋯I_A(1) hydrogen bond is the shortest (2.565 Å) and presumably the strongest. Accordingly, the N–H_N(1)⋯I_A(1) angle is almost 180° (Table 2). These results are in good agreement with the recent powder neutron diffraction study, as well as the structure of ammonium iodide (NH₄I)¹⁹ where hydrogen-bonding is known to play a key role in structural stabilization.²⁰ Details of the computed bond lengths and bond angles for the principal interactions are given in Table 2 and compared with the experimental neutron diffraction values. The N⋯I distances are much shorter than those of C⋯I and also the angle of C–H_C(1)⋯I_A(3) is smaller than that of N–H_N(1)⋯I_A(1). Accordingly, the hydrogen-bond interactions mainly originate from H atoms on nitrogen.

Table 2 Calculated H⋯I, N⋯I, and C⋯I bond lengths and N–H⋯I, C–H⋯I and Pb–I–Pb bond angles compared with NH₄I and ND data for MAPbI₃

	This work	NH4I^19	Experiment^9
Bond type (Å)
HN(1)⋯IA(1)	2.565	2.454	2.613
HN(2)⋯IE(2)	2.611		2.808
HC(1)⋯IA(3)	3.137		3.190
HC(2)⋯IE(4)	3.094		3.006
N⋯IA(1)	3.606	3.473	3.611
N⋯IE(2)	3.566		3.681
C⋯IA(3)	4.221		4.298
C⋯IE(4)	4.032		4.090
Angle type (°)
N–HN(1)⋯IA(1)	177	178	174
N–HN(2)⋯IE(2)	152		146
C–HC(1)⋯IA(3)	170		172
C–HC(2)⋯IE(4)	144		148
Pb–IA(1)–Pb	159		162
Pb–IE(2)–Pb	145		151

---

To understand the strength of the hydrogen-bonding in o-MAPbI₃, we have examined the energetics as a function of the H positions. The most stable structure shows that the MA cation adopts a staggered conformation [Fig. 2(b) inset], as expected and in agreement with the neutron results; this is denoted by θ = 0° in the following calculations in which we rotate the NH₃ and CH₃ groups around the C–N axis for all four MA-ions in the unit cell. We have calculated the Kohn–Sham (K–S) energies as a function of the torsion angle θ, fixing all other atoms, for three rotational modes as shown in Fig. 2(a): (i) in the H_C mode, the NH₃ group is held rigid while the CH₃ group is rotated around the C–N bond axis. This mode reveals the stability of the staggered conformation at θ = 0° compared with the eclipsed conformation at θ = 60° [Fig. 2(b)]. The staggered conformation is ∼66 meV per MA-ion more stable than the eclipsed structure. This suggests that the methyl groups may be rotating at high temperatures, but becomes ordered when the material is cooled to lower temperatures, as is seen in [(CH₃)₂NH₂]Mn(HCOO)₃.²¹ (ii) In the H_N mode, the NH₃ group rotates around the C–N bond axis. This mode reflects the stability of both the staggered conformation and the hydrogen-bonding [Fig. 2(b)], and is the stiffest mode. (iii) In the third mode, H_C–N, the MA cation is held in the staggered conformation and rotated around the C–N axis. This mode is dominated by the strength of the hydrogen-bonding, which represents ∼102 meV per MA-ion. This value is quite close to that obtained by NMR for the corresponding rotation in the dimethylammonium perovskite formates (88 meV per DMA-ion).²¹ Note that the difference in the energy barriers between the H_N and H_C modes is similar to the energy barrier for the H_C–N mode. The three K–S energy curves plotted as a function of the torsion angle θ for the three rotational modes have 120° periodicity in common. This is because the most stable structure represents not only the staggered conformation but also the optimised hydrogen-bonding in the Pb–I cage.

Fig. 2 (a) Schematic view of three MA rotational modes (see text). (b) The computed Kohn–Sham (K–S) energies as a function of the torsion angle θ for the three rotational modes and a staggered conformation (inset) with the torsion angle θ = 0°.

We have also examined the relationship between hydrogen-bonding and the octahedral tilting in o-MAPbI₃. We began by fixing the H_C–N modes at θ = 0° (hydrogen bonding maximised) and 60° (hydrogen bonding minimised) but allowing the inorganic Pb–I cage to relax. We found that the resulting structure and octahedral rotations were almost unaltered. This suggests, unsurprisingly in perovskites with tolerance factors less than 1,^22,23 that the octahedral rotations are robust.

To understand how the robust octahedral rotations (which change amplitude and pattern with temperature), can affect the hydrogen bonding, we performed a second computational experiment. Here, we computed the K–S energies as a function of the torsion angle θ for the same three rotational modes but in a hypothetical high-symmetry structure. This high-symmetry structure is characterized by no octahedral tilting, but with the MA ions in the same configuration as the orthorhombic phase [Fig. 3(a)]. When we calculated the energy of the high-symmetry structure, we relaxed only the MA units within this configuration to optimise any hydrogen-bonding in the Pb–I cage.

Fig. 3 (a) A high-symmetry structure characterised by the removal of octahedral tilting in o-MAPbI₃. (b) The K–S energies plotted as a function of the torsion angle θ for three rotational modes in the high-symmetry structure. Inset compares the energy barriers for the H_C–N mode in the high-symmetry structure and the ground-state structure with octahedral tilting.

In the high-symmetry structure, interestingly, the strength of hydrogen-bonding is remarkably reduced. This is seen for the H_N and H_C–N modes where the energy barriers are significantly shallower in the high-symmetry structure [Fig. 3(b)]. The energy barrier of the H_N mode is ∼90 meV per MA-ion, while that of H_C–N mode is ∼12 meV per MA-ion. The latter is an order of magnitude smaller than the corresponding barrier in the orthorhombic (ground state) phase with octahedral tilting [Fig. 3(b) inset]. On the other hand, the energy barrier of the H_C mode (∼96 meV per MA-ion) of the high-symmetry structure slightly increases. This suggests that the MA cation is held in the staggered conformation and easily rotated around the C–N axis in the high-symmetry structure without invoking the octahedral tilting, in agreement with NMR results.²⁴ Thus, one can conclude that the octahedron tilting is strongly correlated with the hydrogen-bonding interaction. Our findings are consistent with a recent report of a stiffening of the lattice of a hybrid metal–organic formate [NH₄][Zn(HCOO)₃] when cooling through a hydrogen bond ordering transition, as revealed by Resonant Ultrasound Spectroscopy.²⁵ Our results could also explain why the MA-ion ordering does not appear above the orthorhombic transition, but is coincident with the phase transition driven by the robust octahedral rotations.

As a way of visualising the hydrogen-bonding in o-MAPbI₃, we present in Fig. 4 the non-covalent interaction isosurfaces calculated using the Critic2 code^26,27 for the ground-state and high-symmetry structures at θ = 0°. The ground-state case in Fig. 4(a) shows non-covalent interaction density located between three H atoms [H_N(1), H_N(2), and H_N(2′)] and three I atoms [I_A(1), I_E(2), and I_E(2′)], which is consistent with hydrogen bonding in the energy minimised structure. On the other hand, there is only one non-covalent interaction density contour between H_N(1) and I_A(1) in the high-symmetry structure. The octahedron tilting is clearly coupled to the strength of the hydrogen-bonding. To further quantify this strength we have calculated the hydrogen bond strength index²⁸ from the N–H stretching frequencies of the MA molecule in isolation and in the perovskite. For the H_N(1)⋯I_A(1) bond, for example, this index is 0.11 in the ground-state structure but only 0.03 in the high-symmetry structure. For the H_N(2)⋯I_E(2) bond the values are 0.05 and 0.02 respectively. The 11% softening of the N–H stretching frequency in o-MAPbI₃ is comparable to the softening of the O–H stretching frequency in ice,²⁹ suggesting similar strengths of hydrogen bonding in both cases.

Fig. 4 Non-covalent interaction density isosurfaces for (a) the ground-state (θ = 0°) and (b) high-symmetry structures (θ = 0°). The reduced gradient density cut-off is 1.0 au.

We have calculated precise positions for the H atoms in o-MAPbI₃ and found excellent agreement with a recent powder neutron diffraction study.⁹ We have also explored the role of hydrogen bonding in relation to the behaviour of the MA cations and the tilting of the PbI₆ octahedra. The MA cation adopts a staggered conformation in the optimised structure. It is shown that rotation of the staggered MA cation around its C–N bond is quite hindered, with an activation energy of ∼102 meV per MA-ion, due to the strong hydrogen-bonding. The strength of the hydrogen-bonding is found to be highly correlated with the octahedral rotations in MAPbI₃, providing an explanation for why the MA-ordering transition is coincident with the tilting phase transition temperature.

Funding from the Winton Programme for the Physics of Sustainability at the University of Cambridge is gratefully acknowledged. NCB acknowledges financial support from the Royal Commission for the Exhibition of 1851 for a fellowship at Imperial College London. The calculations were performed at the Cambridge HPCS and the UK National Supercomputing Service. Access to the latter was obtained via the UKCP consortium and funded by EPSRC Grant No. EP/K014560/1.

Notes and references

1. http://www.nrel.gov/ncpv/images/efficiency_chart.jpg .
2. T. Baikie, Y. Fang, J. M. Kadro, M. Schreyer, F. Wei, S. G. Mhaisalkar, M. Graetzel and T. J. White, J. Mater. Chem. A, 2013, 1, 5628–5641.
3. C. C. Stoumpos, C. D. Malliakas and M. G. Kanatzidis, Inorg. Chem., 2013, 52, 9019–9038.
4. Y. Wang, T. Gould, J. F. Dobson, H. Zhang, H. Yang, X. Yao and H. Zhao, Phys. Chem. Chem. Phys., 2014, 16, 1424–1429.
5. J. Feng and B. Xiao, J. Phys. Chem. Lett., 2014, 5, 1278–1282.
6. J. M. Frost, K. T. Butler and A. Walsh, APL Mater., 2014, 2, 081506.
7. E. Menéndez-Proupin, P. Palacios, P. Wahnón and J. Conesa, Phys. Rev. B: Condens. Matter Mater. Phys., 2014, 90, 045207.
8. W. Geng, L. Zhang, Y.-N. Zhang, W.-M. Lau and L.-M. Liu, J. Phys. Chem. C, 2014, 118, 19565–19571.
9. M. T. Weller, O. J. Weber, P. F. Henry, A. M. Di Pumpo and T. C. Hansen, Chem. Commun., 2015, 51, 4180–4183.
10. P. E. Blöchl, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953–17979.
11. G. Kresse and J. Hafner, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 558–561.
12. G. Kresse and J. Furthmüller, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169–11186.
13. G. Kresse and D. Joubert, Phys. Rev. B: Condens. Matter Mater. Phys., 1999, 59, 1758–1775.
14. S. Grimme, J. Comput. Chem., 2006, 27, 1787–1799.
15. H. J. Monkhorst and J. D. Pack, Phys. Rev. B: Solid State, 1976, 13, 5188–5192.
16. P. E. Blöchl, O. Jepsen and O. K. Andersen, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 49, 16223–16233.
17. D. A. Egger and L. Kronik, J. Phys. Chem. Lett., 2014, 5, 2728–2733.
18. A. Glazer, Acta Crystallogr., Sect. B: Struct. Crystallogr. Cryst. Chem., 1972, 28, 3384–3392.
19. D. P. Kozlenko, B. N. Savenko, V. P. Glazkov, V. A. Somenkov and S. Hull, Physica B, 2000, 276–278, 226–227.
20. A. M. Heyns, K. R. Hirsch and W. B. Holzapfel, J. Chem. Phys., 1980, 73, 105–119.
21. T. Besara, P. Jain, N. S. Dalal, P. L. Kuhns, A. P. Reyes, H. W. Kroto and A. K. Cheetham, Proc. Natl. Acad. Sci. U. S. A., 2011, 108, 6828–6832.
22. V. Goldschmidt, Naturwissenschaften, 1926, 14, 477–485.
23. G. Kieslich, S. Sun and A. K. Cheetham, Chem. Sci., 2014, 5, 4712–4715.
24. R. E. Wasylishen, O. Knop and J. B. Macdonald, Solid State Commun., 1985, 56, 581–582.
25. Z. Zhang, W. Li, M. A. Carpenter, C. J. Howard and A. K. Cheetham, CrystEngComm, 2015, 17, 370–374.
26. A. Otero-de-la-Roza, M. A. Blanco, A. M. Pendás and V. Luaña, Comput. Phys. Commun., 2009, 180, 157–166.
27. A. Otero-de-la-Roza, E. R. Johnson and V. Luaña, Comput. Phys. Commun., 2014, 185, 1007–1018.
28. X.-Z. Li, B. Walker and A. Michaelides, Proc. Natl. Acad. Sci. U. S. A., 2011, 108, 6369–6373.
29. B. Santra, J. Klines, D. Alfe, A. Tkatchenko, B. Slater, A. Michaelides, R. Car and M. Scheffler, Phys. Rev. Lett., 2011, 107, 185701.

---