Subscriber access provided by UNIV OF CALIFORNIA SAN DIEGO LIBRARIES
Article
On the Development of a Classical Interatomic Potential for MAPbBr3 Tomoyuki Hata, Giacomo Giorgi, Koichi Yamashita, Claudia Caddeo, and Alessandro Mattoni J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.6b11298 • Publication Date (Web): 30 Jan 2017 Downloaded from http://pubs.acs.org on January 31, 2017
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
The Journal of Physical Chemistry C is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 36
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
On the Development of a Classical Interatomic Potential for MAPbBr3 Tomoyuki Hata1,2*, Giacomo Giorgi3, Koichi Yamashita1,2*, Claudia Caddeo4, Alessandro Mattoni4*
1
Department of Chemical System Engineering, School of Engineering, The University
of Tokyo, 7-3-1, Hongo, Bunkyo-ku, Tokyo 113-8656, Japan, 2
CREST-JST, 7 Gobancho, Chiyoda-ku, Tokyo 102-0076, Japan
3
Dipartimento di Ingegneria Civile e Ambientale, Università degli Studi di Perugia,
06125 Perugia, Italy 4
Istituto Officina dei Materiali, CNR-IOM SLACS Cagliari, Cittadella Universitaria,
Monserrato (CA) 09042-I, Italy * corresponding authors E-mail:
[email protected] (Tomoyuki Hata)
[email protected] (Koichi Yamashita)
[email protected] (Alessandro Mattoni)
1 ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 2 of 36
ABSTRACT: We develop a classical interatomic potential for MAPbBr3. The model belongs to the class of MYP force-fields for hybrid perovskites based on two-body Buckhingam-Coulomb and dispersive terms to describe organic-inorganic interactions and already successfully applied to MAPbI3. The model calibration is based on a simplified procedure able to extend one existing parameterization to a different halide by suitable scaling of selected subgroups of parameters. The main static and dynamical properties of MAPbBr3 are well reproduced by the developed model: the lattice constant, cohesive energy curve, bulk modulus, energy barriers for cation rotations (both static and dynamic), the phase transition temperatures and structural parameters evolution with temperature. The model also provides valid relationship between MAPbBr3 and MAPbI3: MAPbBr3 has shorter lattice constant, higher cohesive energy, lower phase transition temperatures, and lager anisotropy in orthorhombic phase. The good comparison extends also to the vibrational properties at finite temperatures that have been benchmarked on experimental and DFT results. The developed MAPbBr3 model is further used to calculate the MA dynamics in MAPbBr3 at room temperature finding a reorientation time of ~3ps in good agreement with experimental data. Present work represents an important step towards the large-scale atomistic modeling of MAPbBr3 and the development of a general class of force fields for hybrid perovskites.
2 ACS Paragon Plus Environment
Page 3 of 36
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
MAIN TEXT: I. INTRODUCTION Hybrid organo-metal halides with formula ABX3 (A = organic molecular cation; B = Ge, Sn, Pb; X = halide such as Cl, Br, I) have received enormous attention in recent years because of their importance for third generation solar cells.1–4 Power conversion efficiency of hybrid perovskite solar cells has rapidly increased from 3.8% in 20095 to 22.1% in early 2016,6 making hybrid perovskites one of the fastest advancing solar technology to date.7 The hybrid nature of such perovskites makes it possible to combine the easy processing from solutions, typical of organic materials, to the unique semiconducting properties of inorganic lead halides. The result is a flexible class of materials of great importance not only for photovoltaics8 but as optoelectronic,9,10 piezoelectric,11–13 or thermoelectric materials.14,15 An intense research effort has tried to clarify the fundamental properties of hybrid perovskites focusing in particular on methylammonium lead iodide (MAPbI3) that is the archetypal hybrid perovskite material. It is nowadays well established that the photovoltaic properties of such materials ultimately derive from the lead halide inorganic sublattice, providing high optical absorption, optimal band gap and electronic structure,16–18 low exciton binding energies,19,20 and high carrier mobilities based on small effective masses21–23 associated also to long diffusion lengths.24,25 The role of molecules is nevertheless important; though not directly involved in the optoelectronic properties, molecules are necessary to achieve structural stability (replacing inorganic cations that even for the largest I-group elements such as Cs are not large enough to satisfy the optimal tolerance for lead halides). However, molecules also introduce sizable molecular entropy (e.g. due to molecular orientations) and give rise to a complex dynamical interplay between the possible molecular orientations and the tilting of lead halide octahedra (and in turn 3 ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 4 of 36
affect the electronic properties related to Pb-I overlap).8,18 This interplay shows up, for example, in structural, vibrational and thermal properties. For example, the equilibrium lattice spacing can vary sizably depending on the actual molecular pattern.26,27 Furthermore, though the heat insulating properties (i.e. the low thermal conductivity) are mainly due to the soft metal halide inorganic lattice, nevertheless molecular rotations and molecular substructure can tune sizably the actual thermal conductivity.14,28–30 Besides MAPbI3, other hybrid perovskites are important such as formamidinium lead iodide or methylammonium lead bromide (CH3NH3PbBr3: MAPbBr3). Such systems have been proposed as alternatives to stabilize hybrid structures4,31 or modify the electronic properties.32–34 Hybrid perovskites offer impressive design possibilities extending in a wide compositional region of metals and halides and molecular ions.35,36 Over 700 compounds can be in principle made with a large number of degrees of freedom for materials design.36 Atomistic theoretical methods are playing an important role in the study and the design of such kind of materials. Static ab initio calculations have been extensively used to study MAPbI326 but also the effect of halide substitutions on band structures,37 electron-phonon coupling38 or sampling the vibrational properties on meta-stable structures.39 However many relevant properties of hybrid perovskites (e.g. the stability of structures, molecular ordering, phase transitions, vibrational and thermal properties) are dynamic in nature involving large structural fluctuations40 with a rich variety of octahedral distortions and molecular disorder. Molecular dynamics naturally goes beyond the static picture since it can sample the configurational space of hybrid crystals under controlled thermodynamic conditions. An increasing number of ab initio molecular dynamics studies, mainly focusing on MAPbI3, have recently appeared in 4 ACS Paragon Plus Environment
Page 5 of 36
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
literature to clarify the effect of structural fluctuations on the electronic gap,4142 the ability of the molecules to rotate,43 the phase transitions,44 the vibrational properties, degradation phenomena45 and single events of defects diffusion.46 Combined strategies in which density functional theory (DFT) calculations are used with effective hamiltonians have been also used to study molecular dynamics in MAPbBr3.47 However ab initio methods cannot fulfill the need of atomistic modeling of hybrid perovskites.27 Large atomistic models (beyond the possibility of ab initio) and long simulation times are necessary to study, for example, the thermal properties,28 calculate defects diffusivity48 or to take into account the compositional disorder of mixed-perovskites49 or the structural complexity of realistic samples containing defects. Classical molecular dynamics (MD) can be a valuable tool in this perspective, provided that model potentials are available to accurately calculate interatomic forces. Only a few model potentials have been proposed for hybrid perovskites and all of them focus on MAPbI3. The first model potential for hybrid perovskites (MYP) was reported by Mattoni et al.50. The MYP model has been applied successfully to study thermally activated properties of MAPbI3, including molecular rotations,50 vibrations,51 defects diffusivity,48 thermal transport29,30 and
electrocaloric effects52 providing results that,
despite the numerical simplicity of the model, compare well with ab initio data27. Models other than MYP have been also proposed. For example, the one developed by Gutierrez-Sevillano et al.53 is able to study organohalide perovskite precursors in solvents; unfortunately it is not able to reproduce the crystalline structure of MAPbI3. Models for crystalline MAPbI3 have been reported by Hata et al.28 and by Qian et al54 to study thermal conductivity. These models involve many-body54 or anharmonic terms28 that are fitted by an automatic and straightforward parameterization on ab initio forces calculated during finite temperature ab initio MD. Though accurate, such models require a fixed bond topology and are hardly transferrable to phases different from the 5 ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 6 of 36
targeted one54 or to disordered/defected systems. According to the above considerations, in this work we focus on the MYP potential considered as a good compromise between accuracy, simplicity and physical foundation and we show that the same MYP functional form can be successfully applied to different perovskites with halides substitutions. In particular we describe a calibration procedure consisting in a step-by-step and physically-sound detuning of parameters based on a limited number of scaling factors. The strategy is demonstrated by developing a force-field for MAPbBr3, intended to reproduce typical results of DFT calculations and experiments including phase transitions. The model accuracy is further assessed by the comparison between MAPbI3 and MAPbBr3 force fields. Notably the Pb-Pb, Pb-MA, and MA internal interactions already set in the MAPbI3 are not modified in MAPbBr3. It means that the developed potential functions are suitable for the large class of mixed systems MAPbIxBr1-x as it is, once mixed I-Br interactions are calibrated. We show that the MAPbBr3 model is able to reproduce the main structural and thermodynamic properties of the material. In particular, the evolution of molecular disorder with temperature is consistent with that of MAPbI3 exhibiting a transition from low temperature in which molecules are not able to rotate to a high temperature regime in which molecular order is lost and the molecules can rotate within a few ps. The force field is also used to study the reorientation time of MA in the crystal and the vibrational density of states of MAPbBr3 as a function of temperature.
II. MYP MODEL POTENTIAL AND CALIBRATION PROCEDURE In this work we adopt the classical MYP interatomic model50 that is extensively described elsewhere27. The total energy is described by the sum of organic-organic, 𝑈!! , organic-inorganic, 𝑈!" , and inorganic-inorganic interactions, 6 ACS Paragon Plus Environment
Page 7 of 36
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
𝑈!! . 𝑈!! is based on the AMBER functional form55 and standard GAFF parameters;56 𝑈!! consists of the Buckingham-Coulomb (BC) potential which includes electrostatic interactions.57 𝑈!" is described as the sum of three terms,50 namely, the Buckingham potentials for (Pb, I)-(C, N) interactions, electrostatic terms, and Lennard-Jones for (Pb, I)-H, as following: 𝑈!" =
𝐴!" 𝑒
!!!" /!!"
!,!
𝜎!" − ! + 𝑟!"
!,!
1 𝑍! 𝑍! + 4𝜋𝜀! 𝑟!"
𝜖!" !,!
𝑠!" 𝑟!"
!"
𝑠!" − 𝑟!"
!
,
(1)
where 𝐴, 𝜌, and 𝜎 are the parameters of Buckingham potential, 𝜀! is the vacuum dielectric constant, 𝑍 is the partial atomic charge, 𝜖 and 𝑠 are the parameters of Lennard-Jones potential, and 𝑟 is the interatomic distance. Considering all the possible interactions in a ternary system as MAPbI3, there are 46 parameters altogether. The original parameters of MAPbI3 were obtained by fitting on a dataset obtained by ab initio DFT calculations.50 Here we detune these parameters to develop a novel force-field for the MAPbBr3 hybrid perovskite by following the concept of MAPbI3 model potential, viz., placing much importance to manage both quantitative accuracy and simplicity of the model. The procedure consists of the following three steps: (i) The first step is aimed at reproducing the experimental lattice constant and DFT cohesive energy curve of cubic MAPbBr3. In the MYP scheme, the lattice constant is mainly controlled by the inorganic BC interactions 𝑈!! , while it is moderately affected by the electrostatic hybrid interactions and only slightly by the Buckingham hybrid interactions. It is discussed in Ref. 27 that it is possible (at fixed charges and taking into account that the dispersive terms of the Buckingham have a minor effect with respect to ionic interactions) to change the equilibrium distance of the BC model !
!
(i.e. 𝑟! → 𝑟! /𝛼) by rescaling its parameters as 𝐴 = 𝛼𝐴′, 𝜌 = ! 𝜌′, and 𝜎 = ! 𝜎′, where primed letters indicate the iodide case. Accordingly it is natural to use 𝛼 as an effective 7 ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 8 of 36
fitting coefficient to adapt the force-field parameters of MAPbI3 to different halides such as MAPbBr3; furthermore, we include the possibility to change the ionic charges 𝑍 by a factor 𝛽. In summary, 𝐴!! = 𝛼𝐴′!! , 𝜌!! =
1 𝜌′ , 𝛼 !!
𝜎!! =
1 𝜎′ , 𝑍 = 𝛽𝑍′. 𝛼 !!
(2)
Here 𝛼 > 1 and 𝛽 > 1 are associated to hybrid crystals with smaller halogens (i.e. Br or Cl). Instead of calibrating both parameters, we
fixed 𝛽 by using ab initio
calculated atomic Bader charges of MAPbBr3 and MAPbI358 (see Supporting Information for details). (ii) The second step is aimed at reproducing the energy barriers for molecular rotations, that is mainly controlled by organic-inorganic 𝑈!" interactions, i.e., Br-(C, N, H). Both Buckingham and Lennard-Jones terms of the new hybrid interactions involving Br are obtained by calibrating a factor 𝛾 scaling the corresponding parameters of the iodine case as 1 1 1 1 𝐴!" = 𝐴′!" , 𝜌!" = 𝜌′!" , 𝜖!" = 𝜖′!" , 𝑠!" = 𝑠′!" . 𝛾 𝛾 𝛾 𝛾
(3)
The other interactions (Pb-(C, N) Buckingham and Pb-H Lennard-Jones) remain unchanged with respect to the iodine case and in the BC hybrid interactions 𝜎!" = 𝜎′!" = 0. 𝛾 > 1 corresponds to shorter inorganic-organic interactions, and at variance with Equation (2), the scaling factors of A coefficients are inverse. There is no fundamental physical motivation for this choice but practical arguments related to a better sensitivity of the rotational barriers on 𝛾 variations. Different choices result in moderate modifications, which allow us to preserve the lattice constant and cohesive energy curve obtained at the previous step during the calibration. In addition, occasionally during the 𝛾 fitting, we refine 𝛼 of step (i) to keep the experimental lattice constants. 8 ACS Paragon Plus Environment
Page 9 of 36
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
(iii) The third step of the procedure consists in refining the dynamical properties and the transition temperatures of MAPbBr3; we observe that after the calibration of static barriers at steps (i-ii) the transition temperatures for the orthorhombic-to-tetragonal phases are not necessarily accurately reproduced. This is well understood because of the importance of dynamical effects that are not necessarily caught by static analysis. In the MAPbI3 case, the improvement in the dynamical properties was obtained by refining the I-N interactions (the size of the Buckingham term, in particular) finding that weaker hybrid interactions (longer repulsive Buckingham) result in lower phase transition temperature. While the organic-inorganic interactions are also important for MAPbBr3 force fields, we observe that the static properties of MAPbBr3 (likely due to shorter lattice spacing) are more sensitive to parameters with respect to the MAPbI3 case and global refinements are necessary. Accordingly, during this step, we modified 𝐴!" and 𝜌!" for Br-(C, N) interactions but we also refined the previous parameters to preserve the static properties obtained at previous steps. The interactions between lead atoms and cations are unchanged with respect to the MAPbI3 case. The above procedure is designed to reduce the number of fitting parameters and simplify the search by associating group of parameters to different physical properties. Overall, the fitting procedure only uses five parameters: three for scaling (α,β,γ) and two for phase refinements (AOI, ρOI). The parameters are optimized based on the Powell conjugate direction method, which does not require any derivatives of functions, followed by manual refinement. As it will be shown below, despite the relative simplicity of the approach, the obtained MAPbBr3 model reproduces the main static and thermodynamic properties of the material at an acceptable level of accuracy.
III. STATIC ENERGY PROFILES 9 ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 10 of 36
In this section, we report the results for the cohesive energy curve and rotational barriers versus the DFT data. MAPbI3 results are also shown for comparison.50 The cohesive energy curves are evaluated as a function of lattice spacing, reported in the left panel of Figure 1. We performed classical calculations by using the DL_POLY code59 on a cubic 8 × 8 × 8 MAPbBr3 periodic crystal
keeping the
(Pm 3 m) symmetry fixed during lattice expansion. Reference DFT results with GGA-PBE exchange-correlation functional60 were obtained by the VASP code.61,62 The reciprocal space was sampled by a 8 × 8 × 8 Monkhorst-Pack grid.63 The energy curves reported in Figure 1 are evaluated by fitting atomistic data to the universal energy relation.50,64,65 The DFT curves data are shifted to the experimental lattice constants.66,67 MYP model reproduces the experimental lattice constant (5.91 Å versus 5.9 Å), and the cohesive energy curvature. The fitted universal energy relation gives the cohesive energy of MAPbBr3 (8.60 eV) and bulk modulus (0.235 Mbar) that can be compared with the corresponding
DFT and experimental
values66–68 (see Supporting Information for details). The MAPbBr3 crystal in cubic phase has a shorter lattice constant than MAPbI3 (5.9 versus 6.27 Å) and higher cohesive energy (8.89 eV versus 7.8 eV). Bulk modulus of MAPbBr3 is also higher (0.235 versus 0.18 Mbar), consistently with experiments and DFT,6768 though the MYP values are slightly overestimated. This corresponds to the fact that, in Figure 1, the MYP energy curves grow more rapidly than DFT with moving from the minimum. This is expected since the electrostatic interactions are due to fixed-point charges and no redistribution is possible. Overall the comparison between MAPbBr3 and MAPbI3 clearly shows that the smaller ionic radius and the larger electron affinity of Br are the reasons for the smaller lattice constant and the stiffer mechanical properties with respect to MAPbI3 (in agreement with the scaling of BC model27). The static rotational barriers are reported in the right panel of Figure 1. Also 10 ACS Paragon Plus Environment
Page 11 of 36
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
in this case we focus on the cubic crystal and rotate cations in (001) crystallographic plane keeping the symmetry and volume fixed. Methodological details correspond to the case of energy-strain calculations described above. We can see that the rotational barrier of MYP is as high as 0.63 eV and almost identical to the one calculated by DFT. The barrier of MAPbI3 is lower that of MAPbBr3 for the whole angle range. The most stable configuration (𝜃 = 0°) is stabilized by the attractive (I, Br)-N interaction and the stronger one for Br results in the higher rotational barrier; however, we point out that at finite temperature, without the rigid symmetry, the actual rotational barrier can be much smaller (as already found for the MAPbI3 case and as discussed in the finite temperature section). The higher static barrier of MAPbBr3 with respect to MAPbI3 is in seeming contrast to the observed transition from orthorhombic to tetragonal phase66,69 occurring at a lower temperature for MAPbBr3; however we anticipate that this is reproduced also by the present model. This means that not all properties related to rotations can be explained by only considering the organic-inorganic interactions inferred from the static rotational barriers. We finally observe that the MYP energy profile shows convexoconcaves at around 90 and 270° that are not present in the DFT curves. We attribute the difference to the point charge approximations. Yet at finite temperature, such features are smoothed by thermal averages.
IV. FINITE TEMPERATURE ANALYSIS In this section we study the structural properties of MAPbBr3 at finite temperature and the phase transitions. There are still several ambiguities of the phase transitions in organic-inorganic halides, especially for MAPbBr3. According to X-ray diffractions,66 MAPbBr3 exhibits four different phases: the orthorhombic phase (Pnma) is found at low temperature up to 140-150 K; then the MAPbBr3 crystal turns into the 11 ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 12 of 36
first tetragonal phase (P4/mmm), but it is easily deformed into a second tetragonal phase (I4/mcm) at around 155 K; by further increasing the temperature beyond 237 K the crystal finally transforms into the cubic phase (Pm3m). In contrast to XRD findings, other experimental studies including neutron diffraction measurements70 have not clearly confirmed the occurrence of two different tetragonal structures (in particular the first P4/mmm tetragonal structure) suggesting rather a metastable phase with dynamic disorder. Yet at least, we can expect that there are roughly three phases, namely orthorhombic, tetragonal, and cubic, whose transition temperatures are located at around 150 K and 240 K. Here we adopt the MYP potential to study the MAPbBr3 crystal at different temperatures and check the model accuracy in reproducing experimental volumes, symmetry, and dynamics of cations as well as temperature boundaries for phase transitions. To this aim we consider a perfect MAPbBr3 orthorhombic 4 × 4 × 4 crystal and we thermalize it in a dense grid of temperatures within 0-350 K under the 𝑁𝜎𝑇 thermobarostat conditions at zero pressure. Each system is equilibrated for at least 0.15 ns, and the physical properties are averaged over successive 0.15 ns. First, we pick up three typical structures at 50, 170, and 300 K to check the crystalline phase of MAPbBr3, (see Figure 2).
The structure at 50 K corresponds to the orthorhombic
phase with characteristic octahedral tilt, 𝑎! 𝑎! 𝑏 ! ; conversely, at 170 K the crystal structure is distorted only along 𝑐-axis, with pattern 𝑎! 𝑎! 𝑐 ! , corresponding to the tetragonal phase.71,72 The 300 K snapshot is described as the averaged structures around 0.30 ns and it clearly shows cubic perovkskite crystal with no tilting, 𝑎! 𝑎! 𝑎! . We also show the results at 155 K, close to the expected transition temperature from orthorhombic-to-tetragonal. At such temperature the crystalline phase can fluctuate, with the possibility to transform from one phase to the other and switch the crystallographic axes during the simulation time. For example, we found that the crystal 12 ACS Paragon Plus Environment
Page 13 of 36
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
is first transformed from orthorhombic to tetragonal-like as shown by the snapshots at 0.15 ns with partial 𝑎! 𝑎! 𝑐 ! tilts and rotations of cations. Such dynamics stands during the subsequent 0.20 ns, and eventually the structure turns back into the orthorhombic phase but with a different tilting direction. The above qualitative observation is consistent to the experimental orthorhombic-to-tetragonal phase transition temperatures. The top panel of Figure 3 shows the evolution with temperature of the volume of MAP(Br/I) expressed as pseudo cubic lattice constant, 𝑉 !/! . The MD confirms the static analysis of cohesive energy; the volume of MAPbBr3 is smaller than that of MAPbI3 throughout the whole temperature range investigated. The volume increases with temperature with a change of slope at 140-170 K for MAPbBr3 and 150-200 K for MAPbI3. In the original work on MAPbI3, this feature was related to the orthorhombic-to-tetragonal phase transition: this is also the case of MAPbBr3. The anisotropy of lattice is evaluated as the ratio between the largest and minimum lattice axis, 𝑐/𝑎, shown in the bottom panel of Figure 3 as a function of temperature. At low temperature, 𝑐/𝑎 > 1 , corresponding to the orthorhombic Pnma structure. Passing over 140-170 K, 𝑐/𝑎 abruptly decreases to ~1, indicating the phase transition from orthorhombic to tetragonal. Comparing with MAPbI3, the experimental change of anisotropy in MAPbBr3 has almost the same trend, but 𝑐/𝑎 in orthorhombic phase is slightly higher than that of MAPbI3. Table 1 shows the detailed comparison between the experimental results.66 The equilibrium structural parameters calculated by the MYP potentials are in reasonable agreement with the experiments in all crystallographic phases. Concerning the orthorhombic phase it is found that 𝑉 !/! and 𝑐/𝑎 of MAPbBr3 are slightly underestimated while the opposite occurs for MAPbI3. In contrast, for the tetragonal phase of both MAPbBr3 and MAPbI3 the calculated anisotropy is close to 1 with no notable changes observed at the tetragonal-to-cubic phase transition. This was already reported for MAPbI350 and it is a limit of the simple two-body potentials for 13 ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 14 of 36
inorganic lattices that are not able to reproduce the covalent character of the bonds and the uniaxial tetragonal distortion of octahedra. We finally focus on the rotational dynamics of MA cations in MAPbBr3 at finite temperature. Similarly to previous analysis,50 in order to analyze the orientation of each molecule, the normalized direction vector of N-C backbone, 𝒏, is described in polar coordinates, as shown in the top panels of Figure 4. The xy-plane corresponds to (001) crystallographic plane and 𝒛 is parallel to tetragonal tilting axis. 𝒏 can be written as 𝒏 = (sin𝜃cos𝜑, sin𝜃sin𝜑, cos𝜃), where 𝜑, 𝜃 are the azimuthal and polar angles, respectively,. The quantity 𝑑𝜑𝑑(cos𝜃) corresponds to the surface element d𝑆 of the molecule rotational space, whose cumulative frequency gives the distribution of cation orientations. In the bottom panels of Figure 4, we report the orientation map at different temperatures. 𝜑 and 𝜃 are represented in the 𝜑, cos𝜃 = 𝑛𝒛 plane for 50
