Phase transition model of ferroelectric 2D hybrid perovskites

Abstract

Understanding structural phase transitions and cation dynamics in low-dimensional hybrid lead halide perovskites is essential for optimizing their photovoltaic and ferroelectric properties. Here, we propose a model describing the ordering of MA and BA cations in the two-dimensional ferroelectric hybrid perovskites BA₂MAPb₂X₇ (X = I, Br). Our model is based on the available structural data and involves short-range framework-mediated interactions between molecular cations. The effective model Hamiltonian is solved using Monte Carlo simulations, which allow us to explore the interactions between the molecular cations reproducing the structural phase transitions and long-range molecular cation order in both perovskites. Our results reveal the appearance of spontaneous electric polarization at the phase transition confirming the experimentally observed ferroelectricity of BA₂MAPb₂Br₇. We discuss the broader applicability of our model for investigating other two-dimensional Ruddlesden–Popper perovskites incorporating different molecular cations.

---

1 Introduction

Hybrid perovskites have attracted significant attention from the scientific community owing to their exceptional optoelectronic properties and potential for use in efficient, solution-processable photovoltaic devices.^1–5 The majority of research in this field has focused on three-dimensional (3D) lead halide hybrid perovskites, which feature corner-sharing PbX₆ octahedra forming extended inorganic frameworks that accommodate molecular cations such as methylammonium (MA), formamidinium (FA), methylhydrazinium (MHy), or aziridinium (AZR) within the lattice cavities.^6–10

Alongside the 3D lead halide perovskites, a diverse family of lower-dimensional compounds has emerged, in which the inorganic slabs are separated by bulky organic spacer cations.^11–13 This structural modulation can give rise to the so-called Ruddlesden–Popper (RP) phases, with the general formula (n ≥ 2), where A and A′ denote a small ‘perovskitizer’ cation and a large spacer monovalent cation, respectively.^14–17 RP-type hybrid perovskites have attracted growing interest because of their highly anisotropic diverse functionalities and potential applications ranging from photovoltaics to ferroics.^18–20 Notably, proper ferroelectric phase transitions are sometimes clearly observed in these low-dimensional systems, in contrast to their 3D counterparts where such effects are typically absent.¹⁰ In general, understanding the mechanisms of structural phase transitions in hybrid perovskites is of high importance because they can strongly influence dielectric response, lattice stability, and charge-transport properties,^10,21 all of which are crucial for device design and performance.

Among the RP hybrid perovskites, the BA₂MA_n−1Pb_nX_3n+1 (BA = butylammonium, X = I, Br) family has emerged as a representative model system exhibiting rich structural and functional behavior.^{15,17,22–25} The n = 2 bromide analogue BA₂MAPb₂Br₇ exhibits ferroelectric properties that emerge in conjunction with a first-order structural phase transition near 350 K. The crystal structure of this compound comprises perovskite-like slabs incorporating orientationally active MA cations. These slabs are separated by bilayers of bulky BA spacer cations, which disrupt the three-dimensional connectivity of the inorganic framework giving rise to a quasi-2D layered architecture (Fig. 1). In the high-temperature (HT) phase, both MA and BA cations display substantial orientational disorder, whereas upon cooling through the transition, long-range order develops. This ordering breaks inversion symmetry and results in the emergence of spontaneous electric polarization along the crystallographic c-axis in the low-temperature (LT) phase.¹⁵ The iodide analogue BA₂MAPb₂I₇ displays two structural phase transitions near 280 and 190 K.¹⁷ Recent works showed that these transitions are closely coupled to organic cation dynamics, as below 280 K the dynamics of BA cations slow down markedly, while the rotational freedom of MA freezes upon entering the LT phase.^17,25

Fig. 1 Structural motifs of BA₂MAPb₂Br₇ in the (a and b) LT and (c) HT phases as viewed along the (a) (011) and (b and c) (100) directions. Structural data taken from ref. 15.

The dynamics and ordering of molecular cations in hybrid perovskites have been extensively investigated using theoretical approaches such as density functional theory (DFT)^26–31 and molecular dynamics (MD) simulations.^29,31–34 DFT calculations provide detailed atomistic insights into cation orientations and interaction energies but are inherently limited to relatively small system sizes and are less suited for capturing entropic effects associated with phase transitions. Temperature-dependent phenomena can be addressed using MD simulations. However, these methods often suffer from pronounced finite-size effects and high computational cost, particularly in the case of ab initio MD. Recently, machine-learning-based approaches have emerged as a promising alternative, enabling high accuracy for larger systems, although their application to hybrid perovskites remains computationally demanding and methodologically challenging.^9,35–42

Another powerful method for studying phase transitions and ordering phenomena is Monte Carlo (MC) simulations of effective model Hamiltonians. This method enables the exploration of many-particle systems and thermal effects over large length scales, however, the atomistic complexity must typically be reduced to a coarse-grained representation. Lahnsteiner et al. demonstrated that, with an appropriately chosen model Hamiltonian, such coarse-grained MC simulations can achieve predictive accuracy comparable to that of machine-learning approaches.³⁵ Effective Hamiltonian models have been successfully employed to reproduce structural phase transitions and related phenomena in 3D lead halide^43–47 and other hybrid perovskites.^48–50 Nevertheless, corresponding modeling approaches for low-dimensional hybrid perovskites, where multiple types of molecular cations coexist within distinct sublattices, remain unexplored.

In this work, we propose a model to describe the molecular cation ordering in low-dimensional hybrid perovskites. As a representative system, we focus on the 2D BA₂MAPb₂Br₇ compound. The model is constructed based on the available structural data and solved using MC simulations. Our simulations successfully reproduce the structural phase transition and the associated ordering of molecular cations observed in both the HT and LT phases. We complement the simulations with experimental heat-capacity data, which we use to estimate the effective interaction energies. Our calculations further reveal the emergence of spontaneous electric polarization indicating the ferroelectric nature of the transition. Finally, we demonstrate that the model can also be readily extended to describe the two phase transitions in BA₂MAPb₂I₇.

2 Sample preparation and experimental details

In order to grow crystals of BA₂MAPb₂Br₇, 4.5 mmol of PbBr₂ was dissolved in 20 ml of HBr under stirring. Then 4 mmol of butylamine and 3 mmol of methylamine was added and the solution was heated to 90°. The solution was cooled to RT with a rate of 5°/h and left at room temperature. Light-yellow plate-like crystals, which grew at the bottom of the vial, were separated from the mother liquid after 3 days and dried at room temperature.

The heat capacity of BA₂MAPb₂Br₇ was measured using the thermal relaxation method implemented in the heat capacity option of a commercial Physical Property Measurement System (PPMS). The measurements were performed in two separate experimental runs, as a single continuous measurement was not possible due to insufficient thermal coupling between the sample and the measurement platform, which was affected by the properties of the used thermal grease. Consequently, two different greases were employed: one suitable for the low-temperature range up to 275 K, and another for the high-temperature range from 250 to 375 K. Details of the measurement method are described in ref. 51. In both runs, the measured sample consisted of several small single crystals to ensure efficient thermal response, with total masses of 1.55 mg (first run) and 1.76 mg (second run).

3 Model and simulation details

The development of a microscopic model for MA and BA cation ordering in BA₂MAPb₂Br₇ relies on a detailed analysis of the available structural data. Our model builds on the X-ray diffraction (XRD) results of Li et al., who resolved the structure of this compound in both its HT (disordered) and LT (ordered) orthorhombic phases.¹⁵ The experimentally determined structural motif of this compound is presented in Fig. 1a, where inorganic slabs containing MA cations are separated along the a-axis by a layer containing two BA cations. Note that the neighbouring inorganic slabs are shifted with respect to each other by half a lattice periodicity along the b-axis, such that the stacking pattern repeats every second slab.

According to this structural model, the LT phase exhibits complete long-range order, with both the framework-confined MA and the spacer BA cations showing no disorder. In this phase, the MA cations within each slab are oriented approximately perpendicular to one another forming a checkerboard pattern in the bc-plane (Fig. 1a and b), similar to that observed in the related 3D MAPbX₃ (X = I, Br, Cl) hybrid perovskites.⁴⁵ To accommodate this cation arrangement, the lead bromide octahedra are tilted forming a network of alternating rhombic cavities, as shown in Fig. 1b. In contrast, the HT phase features a significantly less distorted lead bromide framework and disordered MA cations with four possible orientations within the cavities (Fig. 1c).

In the LT phase, the nitrogen termini of BA cations are directed toward the ammonium groups of the MA cations thereby imposing a similar checkerboard-like arrangement on the BA cations situated in the same plane adjacent to the inorganic slab (Fig. 1a). In the neighboring plane, the BA cations are similarly oriented toward the corresponding inorganic layer mirroring the MA cation arrangement. Note that the displacement of the neighbouring inorganic slabs relative to each other imposes an identical shift on the two BA cation planes within the same BA layer. In the HT phase, each BA cation site exhibits two-fold disorder arising from mirror symmetry normal to the c-axis.

We model the system by mapping the experimental structure on a cubic lattice (xyz-coordinate frame), where each lattice point corresponds to a single MA or BA cation (see Fig. 2a and b). Consistent with the previously noted displacement of neighbouring slabs, this shift is incorporated by offsetting the corresponding model sites by (0, 1/2, 1/2) (Fig. 2b). To mimic the arrangement of the MA cations in the LT phase, we assume that the C–N axis (electric dipole moment) of each MA cation is oriented perpendicular to a cube face resulting in six distinct MA states (S^M = 1, …, 6) (Fig. 2c) consistent with the model of the 3D perovskites.⁴⁵ We also consider eight BA cation orientational states (S^B = 1, …, 8) corresponding to configurations where the BA dipole moment, defined by the axis connecting the terminal N and C atoms, points toward the cube edges (Fig. 2c). Note that in the real structure, BA cations point toward the inorganic layers, and there are no orientations parallel to these layers. As a result, our model considers only eight of the twelve possible BA cation states (Fig. 2c).

Fig. 2 Mapping of (a) the experimental crystal structure of BA₂MAPb₂Br₇ onto (b) a cubic model lattice. (c) The six MA and eight BA cation orientational states included in our model. (d) MA–MA (ε_MM1 and ε_MM2), MA–BA (ε_MB) and BA–BA (ε_BB1 and ε_BB2) interactions between the molecular cations considered in the model.

We take into account interactions between the neighbouring MA and BA molecular cations at specific orientations, while we do not explicitly consider the inorganic framework, though it serves as a medium for interaction transfer. Based on our previous work on the 3D MAPbX₃ system,⁴⁵ MA–MA interactions are considered only when the cations are oriented perpendicular to each other. Depending on whether the dipole of one cation points toward or away from the center of its neighbour, the interaction is assigned an energy of ε_MM1 or ε_MM2, respectively. Examples of these configurations are shown in Fig. 2d, while all possible interacting MA–MA arrangements considered in the model are summarized in Fig. S1.

Based on the cation arrangement in the LT phase, we assume that the heteromolecular interaction between MA and BA cations ε_MB occurs, when the BA cation is directed toward the N-terminus of the MA cation. An example of this interaction is shown in Fig. 2d, while other equivalent configurations obtained by symmetry are presented in Fig. S2. Note that the MA cation states S^M = 3 and 6 do not interact with BA cations, as there are no corresponding BA states that satisfy the required interaction geometry. Consequently, these states are largely suppressed in our simulations, which is consistent with the four MA orientations observed in the XRD data of the HT phase (Fig. 1c).¹⁵

In our model, the inorganic slabs are coupled through the effective BA–BA interactions. Based on the cation distances in the XRD structure, we assume that interactions occur only between BA cations in adjacent planes within the spacer region, while cations within the same plane do not interact. As illustrated in the example presented in Fig. 2d, such interactions take place when the terminal carbon atoms of both BA cations are oriented approximately toward each other. As a result, a single BA cation interacts with four neighbouring BA cations as shown in Fig. S4. The BA configurations in our model give rise to two slightly different BA–BA interactions denoted as ε_BB1 and ε_BB2. The interaction ε_BB1 applies when the two cation orientations form a right angle, whereas ε_BB2 corresponds to a relative angle of 120°. Fig. S3 shows all BA–BA interaction arrangements considered in our model.

The effective Hamiltonian of the model can be written as

where the three contributions correspond to MA–MA, MA–BA and BA–BA interaction terms, respectively. As demonstrated in our previous works,^45,47–49 the explicit form of these Hamiltonians can be conveniently expressed using the Kronecker delta formalism. The full expressions are provided in the SI.

The ground state of the Hamiltonian corresponds to MA cation planes arranged in a checkerboard pattern, similar to the classical antiferromagnetic Ising model.^52,53 This ordering is transferred to the adjacent BA cation plane through the ε_MB interaction. We were not able to map the ordering of the BA cation sublattice onto any conventional lattice model.

Based on our previous work on the 3D MAPbX₃ system,⁴⁵ we adopt the interaction energies ε_MM1 = −43 and ε_MM2 = −13 meV. We expect ε_MB, ε_BB1 and ε_BB2 to be of a similar magnitude, and we explore these parameters in the results and discussion section. It is worth noting that several previous studies of hybrid perovskites and related compounds have also considered dipolar interactions between molecular cations.^{26,43,45–48,54} For instance, we previously demonstrated that such interactions are essential in 3D perovskites for coupling the planes of MA cations arranged in a checkerboard pattern, thereby reproducing the complete phase transition sequence in MAPbX₃.⁴⁵ In the present work, for simplicity, we neglect dipolar interactions, as the MA cations in a slab are oriented perpendicular to each other (yielding zero dipolar interaction), while the interactions involving the much larger BA cations are primarily steric in origin due to their restricted rotational freedom.

We employed a single-flip Metropolis MC algorithm⁵² to solve the model Hamiltonian given by eqn (1). All simulations were initialized from a randomly generated lattice configuration. At each MC step, a molecular cation was randomly selected, and its initial interaction energy E_i was evaluated using the full Hamiltonian of the system. The cation state was then randomly changed to one of the remaining possible orientations, and the corresponding final energy E_f was computed. The new state was accepted with a probability min(1,e^−ΔE/k_BT), where ΔE = E_f − E_i, T is the temperature, and k_B is the Boltzmann constant. Subsequently, another lattice site was randomly chosen, and the Metropolis procedure was repeated. The simulations were carried out on a cubic lattice with periodic boundary conditions containing a total number of V = 3L³ molecular cations with L up to 16. For each temperature, the initial 5 × 10⁵ MC steps were discarded for equilibration, and up to 2.5 × 10⁶ steps were used to calculate the thermodynamic averages.

To investigate the phase transition properties of the model, we calculated the heat capacity at constant volume using the standard fluctuation equation:⁵²

where angle brackets denote MC averages.

The electric polarization was evaluated by averaging the total electric dipole moment of the MA and BA cation lattices:

where is the dipole moment corresponding to the orientation of the cation at lattice site i. We performed DFT calculations (B3LYP/6-31G*) for isolated molecular cations to obtain dipole moment magnitudes of |p_MA| = 2.26 and |p_BA| = 8.94 D.¹⁰

4 Results and discussions

4.1 BA₂MAPb₂Br₇, ε_BB1 = ε_BB2

We first investigated the phase-transition behavior of our model assuming the same magnitude ε_BB for both BA–BA interactions (ε_BB1 = ε_BB2 = ε_BB). We then systematically varied the interaction parameters ε_MB and ε_BB and constructed the corresponding phase diagram from the calculated heat capacity C_V data (Fig. 3). For small ε_BB values relative to ε_MM1, two distinct anomalies appear in the temperature dependence of C_V indicating two phase transitions, where the higher temperature transition is related to the ordering of the MA sublattice (Fig. 3a). Both phase transitions are of a continuous character as evident by the broad anomalies and the unimodal energy histograms close to the transition points (see inset of Fig. 3a). This behavior contrasts with experimental observations, which show a single first-order phase transition in BA₂MAPb₂Br₇.¹⁵ As ε_BB increases, the two anomalies merge into a single transition (Fig. 3b and c), consistent with the experiment. The first-order nature of this transition is confirmed by the appearance of two distinct peaks in the energy histogram calculated at the transition point (inset in Fig. 3b).

Fig. 3 Simulated heat capacity of our model for (a) ε_MB/ε_MM1 = 0.6, ε_BB/ε_MM1 = 0.4, and (b) ε_MB/ε_MM1 = 0.6, ε_BB/ε_MM1 = 0.8. Insets show energy histograms obtained at transition points indicated by the arrows showing the (a) continuous and (b) first-order phase transitions. The higher-temperature transition in (a) originates from the MA–MA interactions as evident by the calculated heat capacity of the MA–MA interactions only (gray curve). The dots in (a and b) indicate temperature points at which the snapshots of cation arrangement were obtained. (c) Simulated heat capacity data for different values of ε_BB/ε_MM1. White curves follow the peaks of the phase transition anomalies showing a crossover from the single-to double-transition region. (d) Phase diagram of our model in the ε_BB and ε_MB parameter space. Crosses indicate interaction sets which reproduce the experimentally observed phase transition temperature of 350 K. Filled squares show parameters at which the heat capacity simulations presented in (a and b) were performed. The region at small ε_MB values marks the approximate boundary of the glassy phase.

The transition behavior of our model also depends on the ε_MB interaction energy, as illustrated in the phase diagram plotted in the (ε_MB, ε_BB) parameter space (Fig. 3d). This diagram defines the range of interaction strengths that reproduce the experimentally observed single first-order phase transition. Notably, the effect of ε_BB is partially compensated by ε_MB as sufficiently strong MA–BA interactions stabilize the cooperative ordering of the cations. In contrast, for weak MA–BA interaction (ε_MB/ε_MM1 ≲ 0.2), the system no longer develops long-range order and instead forms a glassy phase.

Having mapped the phase diagram of our model, we next examine the MA and BA cation arrangement by analyzing simulation snapshots at different temperatures. For the single-transition region, representative configurations of the HT and LT phases are shown in Fig. 4. In the HT phase, both MA and BA cations exhibit orientational disorder with all possible cation states present. Upon cooling into the LT phase, long-range order emerges, with MA cations forming a checkerboard arrangement within the slabs, which are coupled by the ordered BA cations, in agreement with the experimental structural data. Snapshots obtained in the two-transition region of the phase diagram are presented in Fig. S5 and discussed in the SI.

Fig. 4 Snapshots of MA (red) and BA (blue) cation arrangements (a and b) above (k_BT/ε_MM1 = 0.6) and (c and d) below (k_BT/ε_MM1 = 0.1) the phase transition temperature (k_BT₀/ε_MM1 = 0.57) as obtained by MC simulations of our model for ε_MB/ε_MM1 = 0.6 and ε_BB/ε_MM1 = 0.8. View along the (a and c) z-axis (multiple MA and BA layers) and (b and d) x-axis (single MA slab).

We further constrained the magnitudes of the ε_MB and ε_BB interaction energies by fitting our simulations to the experimentally measured heat capacity data (Fig. 5). To isolate the contribution associated with the phase transition, we subtracted the heat capacity arising from lattice and molecular vibrations unrelated to the structural transformation. This background subtraction was performed by fitting a combination of Einstein and Debye models^55–57 to the low-temperature (<150 K) experimental data (inset in Fig. 5). The fits yielded Debye and Einstein temperatures of θ_D = 305.3 ± 6.7 K and θ_E = 66.6 ± 1.1 K, respectively, which are comparable to those reported for other hybrid perovskites and related systems.^55–57 The remaining excess heat capacity was attributed primarily to the phase transition and used for comparison with the simulated results (Fig. 5). The shape of the simulated heat capacity anomaly closely reproduces the experimental data, particularly in the vicinity of the transition, thereby validating our model. We found that multiple combinations of ε_MB and ε_BB parameters reproduce the experimental phase transition temperature of 350 K and provide a satisfactory fit to the measured data (Fig. 5). This allowed us to isolate the region of parameter space consistent with the phase transition behavior of BA₂MAPb₂Br₇ as shown in Fig. 3d. Further narrowing of the identified range would require either atomistic calculations or direct experimental determination (using e.g. NMR⁵⁸ or neutron scattering⁴³ techniques) of the interaction strengths.

Fig. 5 Temperature dependence of the experimental heat capacity of BA₂MAPb₂Br₇ after subtracting the lattice vibrational contribution unrelated to the phase transition (see inset). Simulated heat capacity curves for two markedly different combinations of ε_MB and ε_BB parameters are shown for comparison. A weak anomaly associated with an additional phase transition is observed at 105 K in the experimental data (inset) and can be related to subtle changes in the system.

We selected ε_MB = −43 and ε_BB = −35.9 meV as representative values from the parameter region identified above for the BA₂MAPb₂Br₇ system and used them to calculate the spontaneous electric polarization. Note that other combinations of interaction parameters within this region yield practically identical polarization results. The temperature dependence of the polarization recalculated to match the experimental crystallographic c-axis is shown in Fig. 6. The total polarization was obtained by summing the contributions from the MA and BA sublattices, with molecular dipole moments estimated from DFT calculations (|p_MA| = 2.26 and |p_BA| = 8.94 D). The dipole orientations were assumed to coincide with the model states defined in Fig. 2. Upon cooling, the calculated polarization increases sharply at the phase transition point reflecting the onset of cooperative orientational ordering of both molecular sublattices. The calculated polarization curve agrees well with the experimentally reported non-linear optical (NLO) response of BA₂MAPb₂Br₇ (Fig. 6).¹⁵ Moreover, the calculated saturation polarization of 10.5 µC cm⁻² is of the same order of magnitude as the experimentally measured value of ∼3.6 µC cm⁻². These results demonstrate that the polarization in BA₂MAPb₂Br₇ arises primarily from the molecular cations and confirm that the proposed model successfully captures the key features of the phase transition and molecular dipole ordering in this system.

Fig. 6 Temperature dependence of the simulated total and partial (MA and BA) electric polarizations along the c-axis obtained using ε_MB = −43 and ε_BB = −35.9 meV. The experimental NLO data taken from ref. 15 is presented for comparison.

4.2 BA₂MAPb₂I₇, |ε_BB1| > |ε_BB2|

The iodide analogue BA₂MAPb₂I₇ exhibits two structural phase transitions rather than a single one.¹⁷ The phase diagram in Fig. 3 shows two transitions for small ε_BB values, where the higher-temperature transition corresponds to MA ordering, followed by BA ordering at lower temperatures. This, however, contradicts recent experimental findings, which report the opposite cation ordering sequence.^17,25 To reproduce the correct sequence, we modified our model by assigning different magnitudes to the two BA–BA interactions. We assume |ε_BB1| > |ε_BB2|, since the BA cations in the ε_BB1 configuration are slightly closer to each other than in the ε_BB2 arrangement (see Fig. S4). In the subsequent simulations, we used a ratio of ε_BB1/ε_BB2 = 2.

The calculated heat capacity data for varying ε_BB1 values are shown in Fig. 7a. For weak and moderate BA–BA interactions, the resulting phase diagram closely resembles the ε_BB1 = ε_BB2 case (Fig. 3c) displaying a crossover from two transitions to a single transition region with increasing BA–BA interaction strength. For strong BA–BA interactions (ε_BB1/ε_MM1 ≳ 0.7) the single transition splits again into two anomalies. In this new regime, the higher temperature transition is mostly dominated by the partial BA cation ordering, whereas the lower temperature transition reflects the complete ordering of the system, as confirmed by the C_V calculations for the MA sublattice only (see Fig. 7b). This phase transition sequence is in good agreement with experimental observations.^17,25

Fig. 7 (a) Simulated heat-capacity data for different values of ε_BB1/ε_MM1 assuming ε_MB/ε_MM1 = 0.95 and ε_BB1/ε_BB2 = 2. The white curves trace the phase-transition anomalies and reveal a crossover from a double-transition region to a single-transition regime, followed by entry into a second double-transition region. (b) Simulated temperature dependence of the heat capacity for ε_MB/ε_MM1 = 0.95, ε_BB1/ε_MM1 = 0.8, and ε_BB1/ε_BB2 = 2. The higher-temperature transition originates from the partial BA cation ordering, as confirmed by the heat capacity calculated for the MA–MA interactions only (gray curve). The transition temperatures are close to the BA₂MAPb₂I₇ system. The dot indicates temperature at which the snapshot of cation arrangement in the intermediate phase was obtained.

By exploring the interaction parameter space, we found that the phase transition temperatures extracted from the calculated heat-capacity data for ε_BB1/ε_MM1 = 0.8 closely reproduce those reported for BA₂MAPb₂I₇ (ref. 17) (Fig. 7b). To achieve this agreement, ε_MB/ε_MM1 had to be increased from 0.6 to 0.95 relative to the bromide analogue, while keeping ε_MM1 fixed at −43 meV. Our preliminary simulations indicate that allowing ε_MM1 to vary would broaden the region of parameter space that reproduces the experimentally observed phase transition temperatures, while leaving the qualitative ordering behavior unchanged.

We also examined the ordering of MA and BA cations for this case. In the HT phase, the MA and BA cations exhibit full disorder, whereas they become fully ordered in the LT phase, analogous to the behavior observed in BA₂MAPb₂Br₇ (see Fig. 4). In the intermediate phase (see snapshots in Fig. S6), the BA cations exhibit partial ordering in agreement with a recent NMR study.⁵⁹ This ordering of BA also confines most of the MA states to the yz plane, although the MA sublattice does not develop long-range order in the intermediate phase.

5 Summary and conclusions

In this study, guided by the available structural data, we constructed and numerically solved a phase transition model describing the ordering of MA and BA cations in BA₂MAPb₂X₇ (X = I, Br) 2D hybrid perovskites. The model incorporates six MA and eight BA cation orientations, which interact through short-range strain-mediated interactions. By exploring the intricate phase diagram in the interaction parameter space and fitting the experimental heat capacity data of BA₂MAPb₂Br₇, we determined the magnitude range of effective interactions between the molecular cations. Using these parameters, we demonstrated that the model accurately reproduces the structural phase transitions and cation order in different phases of BA₂MAPb₂Br₇. Our simulations also reproduced the experimentally observed emergence of macroscopic electric polarization confirming the ferroelectric nature of the phase transition in BA₂MAPb₂Br₇. By relaxing the assumption of identical BA–BA interaction magnitudes, our model also captures two phase transitions in BA₂MAPb₂I₇ and reproduces the corresponding sequence of molecular cation ordering highlighting the importance of interaction anisotropy within the BA sublattice.

The model is easily modifiable and thus can be applied to study other low-dimensional hybrid perovskite systems. By adjusting the interaction parameters and cation orientations, it can capture the ordering phenomena in materials with different organic cations and halide compositions. In particular, it can be extended to layered perovskites with different thicknesses of the inorganic slabs (n > 2) allowing investigation of how dimensionality influences the structural phase transitions and ferroelectric properties.

Conflicts of interest

There are no conflicts to declare.

Data availability

Data for this article are available at https://doi.org/10.18279/MIDAS.312864.

Supplementary information (SI) is available. See DOI: https://doi.org/10.1039/d6ta00025h.

Acknowledgements

This research has been carried out in the framework of the “Universities' Excellence Initiative” program by the Ministry of Education, Science and Sports of the Republic of Lithuania under the agreement with the Research Council of Lithuania (project no. S-A-UEI-23-6).

References

1. H. J. Snaith, J. Phys. Chem. Lett., 2013, 4, 3623–3630.
2. M. Grätzel, Nat. Mater., 2014, 13, 838–842.
3. Z. Li, T. R. Klein, D. H. Kim, M. Yang, J. J. Berry, M. F. A. M. van Hest and K. Zhu, Nat. Rev. Mater., 2018, 3, 18017–EP.
4. L. N. Quan, F. P. García de Arquer, R. P. Sabatini and E. H. Sargent, Adv. Mater., 2018, 30, 1801996.
5. W. Li, Y. Ma, Y. Liu, H. Xu, W. Guo, Q. Fan, L. Tang, H. Rong, J. Luo and Z. Sun, Nat. Commun., 2025, 16, 7535.
6. J. Jeong, M. Kim, J. Seo, H. Lu, P. Ahlawat, A. Mishra, Y. Yang, M. A. Hope, F. T. Eickemeyer, M. Kim, Y. J. Yoon, I. W. Choi, B. P. Darwich, S. J. Choi, Y. Jo, J. H. Lee, B. Walker, S. M. Zakeeruddin, L. Emsley, U. Rothlisberger, A. Hagfeldt, D. S. Kim, M. Grätzel and J. Y. Kim, Nature, 2021, 592, 381–385.
7. J. J. Yoo, G. Seo, M. R. Chua, T. G. Park, Y. Lu, F. Rotermund, Y.-K. Kim, C. S. Moon, N. J. Jeon, J.-P. Correa-Baena, V. Bulović, S. S. Shin, M. G. Bawendi and J. Seo, Nature, 2021, 590, 587–593.
8. M. Ma̧czka, A. Gagor, J. K. Zareba, D. Stefanska, M. Drozd, S. Balciunas, M. Simenas, J. Banys and A. Sieradzki, Chem. Mater., 2020, 32, 4072–4082.
9. M. Ma̧czka, M. Ptak, A. Ga̧gor, J. K. Zarȩba, X. Liang, S. Balčūnas, O. A. Semenikhin, O. I. Kucheriv, I. A. Gural’skiy, S. Shova, A. Walsh, J. Banys and M. Šimėnas, Chem. Mater., 2023, 35, 9725–9738.
10. M. Simenas, A. Gagor, J. Banys and M. Ma̧czka, Chem. Rev., 2024, 124, 2281–2326.
11. H. Lin, C. Zhou, Y. Tian, T. Siegrist and B. Ma, ACS Energy Lett., 2018, 3, 54–62.
12. K. Hong, Q. V. Le, S. Y. Kim and H. W. Jang, J. Mater. Chem. C, 2018, 6, 2189–2209.
13. Y. Shang, Y. Liao, Q. Wei, Z. Wang, B. Xiang, Y. Ke, W. Liu and Z. Ning, Sci. Adv., 2019, 5, eaaw8072.
14. C. C. Stoumpos, D. H. Cao, D. J. Clark, J. Young, J. M. Rondinelli, J. I. Jang, J. T. Hupp and M. G. Kanatzidis, Chem. Mater., 2016, 28, 2852–2867.
15. L. Li, X. Liu, Y. Li, Z. Xu, Z. Wu, S. Han, K. Tao, M. Hong, J. Luo and Z. Sun, J. Am. Chem. Soc., 2019, 141, 2623–2629.
16. X. Li, J. M. Hoffman and M. G. Kanatzidis, Chem. Rev., 2021, 121, 2230–2291.
17. M. Ma̧czka, S. Smółka and M. Ptak, Materials, 2024, 17, 2503.
18. Q. Zhang, A. Solanki, K. Parida, D. Giovanni, M. Li, T. L. C. Jansen, M. S. Pshenichnikov and T. C. Sum, ACS Appl. Mater. Interfaces, 2019, 11, 13523–13532.
19. Z. Wu, W. Zhang, H. Ye, Y. Yao, X. Liu, L. Li, C. Ji and J. Luo, J. Am. Chem. Soc., 2021, 143, 7593–7598.
20. P. Siwach, P. Sikarwar, J. S. Halpati and A. K. Chandiran, J. Mater. Chem. A, 2022, 10, 8719–8738.
21. I. Anusca, S. Balčiūnas, P. Gemeiner, S. Svirskas, M. Sanlialp, G. Lackner, C. Fettkenhauer, J. Belovickis, V. Samulionis, M. Ivanov, B. Dkhil, J. Banys, V. V. Shvartsman and D. C. Lupascu, Adv. Energy Mater., 2017, 7, 1700600.
22. L. Li, Z. Sun, P. Wang, W. Hu, S. Wang, C. Ji, M. Hong and J. Luo, Angew. Chem., 2017, 56, 12150–12154.
23. M. C. Gélvez-Rueda, E. M. Hutter, D. H. Cao, N. Renaud, C. C. Stoumpos, J. T. Hupp, T. J. Savenije, M. G. Kanatzidis and F. C. Grozema, J. Phys. Chem. C, 2017, 121, 26566–26574.
24. W. Paritmongkol, N. S. Dahod, A. Stollmann, N. Mao, C. Settens, S.-L. Zheng and W. A. Tisdale, Chem. Mater., 2019, 31, 5592–5607.
25. A. A. Koegel, I. W. H. Oswald, C. Rivera, S. L. Miller, M. J. Fallon, T. R. Prisk, C. M. Brown and J. R. Neilson, Chem. Mater., 2022, 34, 8316–8323.
26. C. Motta, F. El-Mellouhi and S. Sanvito, Phys. Rev. B, 2016, 93, 235412.
27. D. H. Fabini, T. A. Siaw, C. C. Stoumpos, G. Laurita, D. Olds, K. Page, J. G. Hu, M. G. Kanatzidis, S. Han and R. Seshadri, J. Am. Chem. Soc., 2017, 139, 16875–16884.
28. W. Sukmas, U. Pinsook, P. Tsuppayakorn-aek, T. Pakornchote, A. Sukserm and T. Bovornratanaraks, J. Phys. Chem. C, 2019, 123, 16508–16515.
29. M. T. Weller, O. J. Weber, P. F. Henry, A. M. Di Pumpo and T. C. Hansen, Chem. Commun., 2015, 51, 4180–4183.
30. M. Simenas, S. Balciunas, J. N. Wilson, S. Svirskas, M. Kinka, A. Garbaras, V. Kalendra, A. Gagor, D. Szewczyk, A. Sieradzki, M. Ma̧czka, V. Samulionis, A. Walsh, R. Grigalaitis and J. Banys, Nat. Commun., 2020, 11, 5103.
31. K. Druzbicki, R. Laven, J. Armstrong, L. Malavasi, F. Fernandez-Alonso and M. Karlsson, J. Phys. Chem. Lett., 2021, 12, 3503–3508.
32. M. A. Carignano, Y. Saeed, S. A. Aravindh, I. S. Roqan, J. Even and C. Katan, Phys. Chem. Chem. Phys., 2016, 18, 27109–27118.
33. O. J. Weber, D. Ghosh, S. Gaines, P. F. Henry, A. B. Walker, M. S. Islam and M. T. Weller, Chem. Mater., 2018, 30, 3768–3778.
34. R. Namakian, M. A. Garzon, Q. Tu, A. Erdemir and W. Gao, ACS Nano, 2024, 18, 22926–22937.
35. J. Lahnsteiner, R. Jinnouchi and M. Bokdam, Phys. Rev. B, 2019, 100, 094106.
36. H. Grüninger, M. Bokdam, N. Leupold, P. Tinnemans, R. Moos, G. A. De Wijs, F. Panzer and A. P. M. Kentgens, J. Phys. Chem. C, 2021, 125, 1742–1753.
37. M. Ahmadi, M. Ziatdinov, Y. Zhou, E. A. Lass and S. V. Kalinin, Joule, 2021, 5, 2797–2822.
38. Q. Tao, P. Xu, M. Li and W. Lu, npj Comput. Mater., 2021, 7, 23.
39. K. Fykouras, J. Lahnsteiner, N. Leupold, P. Tinnemans, R. Moos, F. Panzer, G. A. de Wijs, M. Bokdam, H. Grüninger and A. P. M. Kentgens, J. Mater. Chem. A, 2023, 11, 4587–4597.
40. J. Wiktor, E. Fransson, D. Kubicki and P. Erhart, Chem. Mater., 2023, 35, 6737–6744.
41. S. Dutta, E. Fransson, T. Hainer, B. M. Gallant, D. J. Kubicki, P. Erhart and J. Wiktor, J. Am. Chem. Soc., 2025, 147, 37019–37029.
42. T. Hainer, E. Fransson, S. Dutta, J. Wiktor and P. Erhart, Nat. Commun., 2025, 16, 8775.
43. A. M. A. Leguy, J. M. Frost, A. P. McMahon, V. G. Sakai, W. Kochelmann, C. Law, X. Li, F. Foglia, A. Walsh, B. C. O'Regan, J. Nelson, J. T. Cabral and P. R. F. Barnes, Nat. Commun., 2015, 6, 7124.
44. L. Z. Tan, F. Zheng and A. M. Rappe, ACS Energy Lett., 2017, 2, 937–942.
45. M. Simenas, S. Balciunas, M. Ma̧czka, J. Banys and E. E. Tornau, J. Phys. Chem. Lett., 2017, 8, 4906–4911.
46. M. Simenas, J. Banys and E. E. Tornau, J. Mater. Chem. C, 2018, 6, 1487–1494.
47. M. Simenas, S. Balciunas, M. Ma̧czka and J. Banys, J. Mater. Chem. C, 2022, 10, 5210–5217.
48. M. Simenas, S. Balciunas, M. Ma̧czka, J. Banys and E. E. Tornau, Phys. Chem. Chem. Phys., 2016, 18, 18528–18535.
49. M. Simenas, A. Ibenskas, A. Stroppa, A. Gagor, M. Ma̧czka, J. Banys and E. E. Tornau, J. Phys. Chem. C, 2019, 123, 19912–19919.
50. C. Coates, H. Gray, J. Bulled, H. Boström, A. Simonov and A. Goodwin, Philos. Trans. R. Soc., A, 2019, 377, 20180219.
51. J. S. Hwang, K. J. Lin and C. Tien, Rev. Sci. Instrum., 1997, 68, 94–101.
52. D. P. Landau and K. Binder, A Guide to Monte Carlo Simulations in Statistical Physics, Cambridge University Press, 3rd edn, 2009.
53. F. Y. Wu, Rev. Mod. Phys., 1982, 54, 235–268.
54. J. M. Frost, K. T. Butler and A. Walsh, APL Mater., 2014, 2, 081506.
55. S. Kawachi, M. Atsumi, N. Saito, N. Ohashi, Y. Murakami and J.-I. Yamaura, J. Phys. Chem. Lett., 2019, 10, 6967–6972.
56. A. van Hattem, J.-C. Griveau, E. Colineau, A. J. E. Lefering, R. J. M. Konings and A. L. Smith, J. Phys. Chem. C, 2023, 127, 22808–22816.
57. C. Tower, F. S. Razavi, J. Dion, J. Nuss, R. K. Kremer and M. Reedyk, Low-temperature structural instabilities of the halide double perovskite Cs₂AgBiBr₆ investigated via x-ray diffraction and infrared phonons, arXiv, 2025, preprint arXiv:2505.10563,  DOI:10.48550/arXiv.2505.10563.
58. L. Piveteau, V. Morad and M. V. Kovalenko, J. Am. Chem. Soc., 2020, 142, 19413–19437.
59. F. Lyu, X. Zheng, Z. Li, Z. Chen, R. Shi, Z. Wang, H. Liu and B.-L. Lin, Chem. Mater., 2021, 33, 3524–3533.

---