Exploring the guest–host relationship for zeolite Y: a synergistic structural and theoretical investigation

Abstract

Understanding the mobility and distribution of cations in the presence of additional guest species within the zeolite framework is inherently very complex because several thermal and chemical processes occur simultaneously. In this study, the dehydration of zeolite Na-Y was studied using in situ powder X-ray diffraction (XRD). The crystal structure parameters, both the evolution of lattice parameters and the occupancies of sodium at specific positions in the model, were monitored upon heating. A complementary computational study was performed to understand interactions between the framework, cations, and water at the molecular level. Grand Canonical Monte Carlo (GCMC) simulations provided initial insights into water adsorption. Machine Learning Potentials (MLPs) were then trained to the ab initio Potential Energy Surface (PES) using a deep neural network, modeling the dynamics of the framework and cations at various water loadings. Strong agreement between computational results and experimental data reveal that upon dehydration, zeolite Na-Y initially contracts due to water removal, but subsequently expands as sodium cations migrate to the double 6-membered rings (site I). This study demonstrates significant benefits of integrating parametric Rietveld refinement and Machine Learning assisted Molecular Dynamics simulations in understanding dynamic behavior of guest molecules in nanoporous materials at operating conditions and interpreting complex and convoluted experimental data.

---

Introduction

Zeolites are crystalline aluminosilicate materials characterized by their three-dimensional porous structures, enabling their use in a wide range of industrial applications including: catalysis, adsorption, and ion exchange.¹ The zeolite framework is defined as tetrahedral units of silicon or aluminum atoms linked by oxygen atoms (TO₂). The presence of aluminum (Al³⁺) within the structure introduces a net negative charge, which is counterbalanced by cations such as protons, alkali metal ions, or transition metal ions. These counter ions are situated in the pores and voids and are usually mobile. Zeolites can be organized into numerous frameworks and there are currently 258 known topologies.²

The faujasite type zeolite is commonly used in catalytic applications, particularly in fluid catalytic cracking. Faujasite belongs to the cubic crystal system, with a unit cell parameter of approximately 24.7 Å. A low-silica faujasite is referred to as zeolite X, while a high-silica variant (Si/Al > 2) is referred to as zeolite Y. It contains three distinct types of cages through which counter ions can migrate: double six-membered ring units, sodalite cages, and 12-membered ring units often called the supercage (see Fig. 1). Several cationic sites within the faujasite framework have been identified. Site I is found in the double six-membered rings, while site II is located in the supercage adjacent to the hexagonal face of a sodalite cage. Additionally, site I′ and site II′ are situated within the sodalite cage, where site I′ is adjacent to the double six-membered ring and site II′ is adjacent to the hexagonal ring shared with the supercage.³

Fig. 1 Graphical representation of a faujasite-type zeolite structure. Colored circles represent different sites between which the counter ions migrate depending on the hydration state of the zeolite.

Sodium faujasite is a precursor for various faujasite-type zeolites and has been extensively studied experimentally and with molecular mechanics to understand cation mobility in relation to zeolite hydration.^4–10 Mortier et al.⁴ employed powder X-ray diffraction to determine the location of sodium cations in zeolite Na-Y upon water adsorption at various temperatures. Their refinement procedures identified preferred sodium occupied sites in hydrated and dehydrated samples, concluding that the presence of water enhances Na occupancy at site I′ while reducing occupancy at sites I and II. These results were reinforced in a computational study by Beauvais et al.,⁵ who used replica-exchanged canonical Monte Carlo simulations to model water adsorption into Na-Y zeolite, maintaining a rigid framework while permitting sodium cation mobility. Results from their work provided numerical evidence of cation redistribution dependent on water content and indicated that increased hydration results in the migration of cations from site I to site I′. The group led by Boutin⁷ further explored cation redistribution upon dehydration of Na₅₈Y zeolite through neutron scattering and Monte Carlo simulations, observing consistent experimental and computational results regarding changes in cation locations with decreasing water content. Discrepancies between experimental and computational results at high water content were primarily linked to limitations in force fields and the need to consider framework flexibility.

Previous studies, in addition to the improved understanding of zeolite framework flexibility,¹¹ demonstrate the need for continued investigation attributed to the influence of guest species on the zeolite framework with methodology beyond a force field approach. Further, a combined experimental and computational approach will enable the most thorough examination of these phenomena. Experimentally, in situ X-ray diffraction (XRD) is a particularly sensitive method for probing these interactions. It enables real-time monitoring of crystal structure parameters and identification of sodium occupancies through Rietveld refinement. However, achieving quantum-level computational accuracy in studying these interactions is computationally heavy due to the complexity and size of the zeolite Y model. Therefore, we employ a novel, cost-effective approach to overcome these challenges. By training an accurate Machine Learning Potential (MLP) with a water adsorption dataset generated using Grand Canonical Monte Carlo (GCMC) simulations, we were able to accurately model the intermolecular interactions at an accelerated computational speed. With the MLP we studied the dynamics of the dehydration process and cation migration at various temperatures. The combination of the in situ XRD technique with these complementary simulation methods allows us to examine both the mobility of cations and water within the zeolite, as well as the dynamic response of the zeolite framework. Our findings reveal an unusual structural behavior of the zeolite framework, associated with the mobility of cations to specific sites within the structure. The insights gained from our combined experimental and computational approach contribute to a deeper understanding of the nature of sodium zeolite Y and its relationship with water molecules. Understanding of the non-monotonous behavior of Na-Y upon dehydration allows to interpret experimental data with greater precision.

Methodology

Experimental

Catalyst characterization. The catalyst was provided by Zeolyst International (CBV 100) with silicon to aluminum ratio of 2.55, corresponding to 54 Al sites per unit cell. The amount of water in the catalyst was quantified by TGA on a Netzsch STA 449 TG. Approximately 20 mg of catalyst was heated under 25 ml min⁻¹ of synthetic air to 800 °C at a rate of 10 °C min⁻¹. Total mass loss until 300 °C is assigned to water release.

In situ synchrotron radiation powder X-ray diffraction (SR-PXRD). Powder X-ray diffraction (PXRD) data were collected at the DanMAX beamline of the MAX IV Laboratory in Lund, Sweden, with a wavelength of 0.496632 Å, using a DECTRIS PILATUS3 X 2M CdTe area detector. The wavelength was determined from individual runs with Si standard. The capillary reactor (0.7 mm OD borosilicate glass) was mounted on a Swagelok connected clamp on top of the hexapod stage capable of translations for alignment.¹² The data were azimuthally integrated to one-dimensional powder patterns with in-house software with a 2θ range from 0 to 33°. For in situ dehydration experiments, the capillaries were mounted with one end open to atmospheric air and heated from 50 °C up to 400 °C at 5 °C min⁻¹ while flowing 10 ml min⁻¹ of helium. PXRD patterns were continuously recorded on a single point in the capillary with an extremely fast time resolution (6 scans per s). The selected scans were analyzed by parametric Rietveld refinement¹³ using TOPAS software¹⁴ and refined unit cell parameters were extracted. While the T-atom of the framework, sodium and water (dummy oxygen atoms) atomic positions were fixed, the oxygen atomic positions of the framework were allowed to move by using bond length and angle restraints and constraints. Occupancies of sodium and dummy oxygen were calculated in the refinement. The water content was estimated from the dummy oxygen occupancies considering that a correction factor of 0.89 is needed since the model does not include hydrogen atoms. A full list of refined parameters and results can be found in SI S2.

Computational

Training dataset generation (GCMC methods). To estimate the degree of water loading at various temperatures and to create an extensive dataset necessary for training machine learning potentials (MLPs), we employed the Grand Canonical Monte Carlo (GCMC) methodology, wherein the temperature, volume, and chemical potential of the system remain constant. GCMC simulations were performed using the open-source RASPA code.¹⁵ The zeolite model has Si/Al = 2.55 with Al randomly distributed, following Lowenstein's rule.¹⁶ Lennard-Jones potentials described the van der Waals interactions between guest–host and guest–guest molecules using previously reported parameters,¹⁷ the all-atom TIP4P model was used for water¹⁸ and atomic charges were taken from the work by Jaramillo et al.¹⁹ This approach effectively captures non-bonded interactions, providing accurate qualitative insights into the adsorption behavior of guest molecules, as supported by previous studies.^20,21 All details of GCMC simulations can be found in SI S3.

The GCMC simulations provided an initial guess of number of water molecules at various temperatures. Additionally, from GCMC 475 structural snapshots with different water loadings were extracted to construct the MLP training set.

Machine learning potentials training. An active learning approach was used to train MLPs using the psiflow simulation engine.²² The neural network (NN) interatomic potential model architecture used was MACE,²³ a type of equivariant message-passing graph NN. Random perturbations to the atomic positions and cell parameters were applied twice to each of the training structures obtained with the GCMC method, resulting in 950 new snapshots. An initial MACE MLP was trained using these new structures, labeled with the energy and force calculated using Density Functional Theory (DFT).

This initial MLP was used to start four sets of active learning loops, where MLPs are iteratively trained and applied to perform MLP MD simulations at increasing temperatures to explore the phase space. The first three iterations of MLP NPT-MD were conducted at 200 K and run for 0.1 ps. Then, five more iterations of MLP NPT-MD were run at 500 K for 1.5 ps, followed by three iterations at 700 K for 1.0 ps. Finally, three iterations at 900 K for 1.0 ps were performed. After each iteration, the final MD snapshot was recomputed with DFT and added to the total training set. More details of the MLPs training with the schematic explanation of active loop are provided in the SI S4 and Fig. S1.

All DFT evaluations were performed using version 8.2 of the CP2K simulation package.²⁴ The Gaussian-type pseudopotentials Goedecker–Teter–Hutter (GTH)²⁵ with the PBE exchange–correlation functional,²⁶ Grimme D3 dispersion corrections,²⁷ and TZVP MOLOPT basis set²⁸ were applied. The plane-wave cutoff energy was set to 1400 Ry to ensure energy and force convergence.

MLP NPT-MD simulations. MD simulations with trained MLPs were performed using OpenMM software in the NPT ensemble.^23,29 Temperature and pressure were held constant, with a timestep of 1.0 fs. Each simulation was run for 400 ps (400 000 MD steps). The temperature was controlled using a Langevin thermostat with a time constant of 100 fs, and the pressure was controlled using a Monte Carlo barostat with a time constant of 1 ps. Information about energy, temperature and unit cell volume were recorded every 100 steps. Trajectories were recorded every 1000 steps.

Identification of sodium occupancies and unit cell parameter. Sodium occupancies were calculated as an average from the last 100 trajectory snapshots of MLP NPT-MD simulations. To classify the position of each cation, a distance-based method with a specified cutoff radius for each reference position is applied. The cutoff distances were chosen dependent on the cation site: site I – 1.8 Å, site I′ – 2.0 Å, site II′ – 2.0 Å, site II – 2.4 Å. These cutoff values allow the classification of cation positions into known sites based on their proximity, or they leave the positions unclassified if they lie outside these defined boundaries (SI S5). The coordinates for each site position are derived from Rietveld refinement analysis.

The unit cell parameter was calculated as an average of the data recorded in the last 100 ps of the MD simulations. It was assumed that the symmetry of the framework remained unchanged, and that the zeolite Na-Y maintained a cubic structure throughout each MD simulation.

Results and discussion

Assessment of MLP validity

To validate the computational predictions for the Na-Y zeolite generated by the MLP model, we executed a validation protocol based on the methodology outlined in SI S12. In short, MLP energies and forces of the snapshots extracted from the MLP MD simulations were compared with those recomputed at the DFT level of theory. The calculated root mean square error (RMSE) values for energy and forces, are presented in Table 1. The RMSE values of 12 meV per atom for energy and 59 meV Å⁻¹ for forces were derived. While these values are comparable to those reported by Erlebach et al.³⁰ we recognize that the energy error is relatively high, which could be attributed to the limited size of the initial training dataset. A small dataset may not have adequately captured the full range of configurations and interactions present in the system. Despite these limitations, the comparison with experimental data, as explained below, suggests that the trained MLP model effectively captures the complex dynamics of zeolite Na-Y, including cation and water mobility. Specific error values for each trajectory comparison are documented in the SI, Fig. S10.

Table 1 Root mean square errors in energy [meV per atom] and forces [meV Å⁻¹] from MD simulations trajectories obtained using MLP and with respect to the PBE+D3 reference

Energy	Forces
12	59	This work
5.3	186	Ref. 30

---

Water quantification during zeolite Na-Y dehydration

Prior to discussing the water content of the as-received zeolite, it is essential to note that the elemental analysis conducted via energy-dispersive X-ray spectroscopy (EDX) revealed a composition that closely aligns with the expected values (SI S12). This is further corroborated by solid-state ²⁷Al-NMR measurements, which indicate that the crystalline structure of the material is nearly ideal (SI S13). The initial water content of the as-received zeolite Na-Y, as determined by TGA, was approximately 21.75 wt%, which corresponds to roughly 196 water molecules per unit cell. At 50 °C the same zeolite contained around 187 water molecules per unit cell (SI S6, Fig. S3). Water molecules were also identified and quantified through parametric Rietveld refinement of data obtained from in situ X-ray diffraction, and using three distinct positions predominantly occupied by adsorbed water molecules previously introduced by Agostini et al.³¹ (SI S9). The XRD experiment identified 245 water molecules per unit cell at the initial 50 °C temperature. The discrepancy between TGA (187) and XRD refinement (245) water quantification results is reasonably small and can be attributed to experimental errors of both techniques, as well as the variations in external conditions, such as humidity. The known hydrophilicity of zeolite Y, in combination with differing conditions during TGA and in situ XRD, and relatively large uncertainty of the calculated water content from XRD, accounts for varying numbers of adsorbed water molecules.

To have a proper estimate of the number of adsorbed water molecules in terms of temperature, multiple GCMC simulations were conducted, and the average water content was recorded following system equilibration. Apart from directly comparing the changes in water content both from theory and experiments, a proper assessment of the number of adsorbed water molecules and their position is necessary to initiate the MLP training. At an initial temperature of 50 °C, the H₂O content was determined to be 183 molecules per unit cell, closely aligning with TGA calculations and similar in magnitude to the data obtained from Rietveld refinement.

When comparing the refined and simulated changes in water content as temperatures increase, both methods exhibit an exponential decay pattern, with a steeper slope seen in the Rietveld refinement results compared to the GCMC simulation outcomes (Fig. 2). Thus, the dehydration process as simulated using our force-field method is slower compared to the experimental analysis. These differences in decay rate could be attributed to variations in the partial pressure of water vapor between the simulations and the XRD dehydration experiment, during which the vapor partial pressure was not controlled or measured. It is also possible that the force field description used in the GCMC simulations might overestimate number of adsorbed water molecules at higher temperatures. Overall, the qualitative match is excellent, and the quantitative match is very good given the approximations made in the simulations and the limitations of the experimental techniques.

Fig. 2 Changes in number of adsorbed H₂O molecules into zeolite Na-Y as calculated from Rietveld refinement and from GCMC simulations.

Sodium cations and water distribution

Rietveld refinement was also applied to estimate the relative occupancies of each of the established Na sites within the faujasite-type framework.³² This approach enabled the calculation of sodium occupancy at each of the identified sites throughout the duration of the zeolite dehydration. The evolution of the relative occupancy during the dehydration process is plotted in Fig. 3A. As dehydration progresses, migration of sodium cations is evident. The initially vacant site Na I increasingly accumulates cations within the double six-rings and surpasses 0.5 relative occupancy when zeolite is fully dehydrated. This migration is balanced by the decrease of occupancies at sites Na I′ and Na II′. The relative occupancy of Na I′ declines from over 0.75 to around 0.4, while Na II′ transitions from half-occupied to empty as dehydration proceeds. Based on MLP NPT-MD simulations, further insights into the relative occupancies of sodium sites and water molecules were assessed. The results are shown in Fig. 3B. In the fully hydrated state, sites I′ and II are nearly half-occupied, site II′ has a low occupancy of around 0.13, and site I is empty of cations. Approximately half of the cations in the unit cell are unclassified and can be found in the supercage surrounded by water molecules. As water content in the zeolite Y system decreases, the number of unclassified sodium cations drops significantly, site II′ becomes unoccupied, and sodium is found more frequently at sites II and I with relative occupancies of 0.75 and 0.33 respectively upon complete dehydration. Both experimental and computational results reveal preferences of sodium cations to occupy specific conventional sites at different hydration states of the zeolite Y.

Fig. 3 Relative occupancy of each type of Na⁺ cation per unit cell obtained from (A) the Rietveld refinement and (B) from MLP NPT-MD simulations, as a function of temperature during dehydration, with the grey line indicating the number of unclassified sodium cations. (C) Schematic representation of the changes in the cationic distribution and water content at low temperature and at high temperature.

An additional analysis of the cations and water molecules mobility was performed by analyzing their mean square displacement (MSD). Self-diffusion coefficients of sodium cations and water were obtained from the MLP MD trajectories using the Einstein relation, as implemented in the MDAnalysis Python library.^33–35 The MSD curves and corresponding self-diffusion coefficients are reported in SI S7, Fig. S5 and S6. The analysis revealed that water molecules exhibit greater mobility compared to sodium ions. This observation is attributed to the presence of preferred cationic sites within the zeolite material, which inhibit sodium mobility by positioning it either surrounded by water molecules or close to the zeolite framework. No increase in mobility was observed at elevated temperatures. Two main factors can be associated with this thermodynamically counterintuitive phenomenon. First, the reduced number of water molecules will very likely lead to a reduction in water assisted cationic mobility. Second, the increased sodium ions population of sites I and II at high temperatures: these sites are preferred as they minimize cation–cation repulsion and are stabilized by strong interaction with the zeolite framework. Once sodium is located there, the mobility remains unchanged despite rising temperatures. More extensive discussion can be found in SI S7.

In order to understand the discrepancies between experimental and computational initial relative occupancies, the total number of sodium cations per unit cell was examined. Rietveld refinement indicated a decrease in sodium cations from 60 to 45 during dehydration (SI S8). This reduction, unlikely to result from actual cation removal, is attributed to the difficulty of capturing all the real sodium positions with a simple average model containing only 4 sodium sites (I, I′, II′ and II). Thus, some of the electron density that was identified as sodium cations may actually correspond to water molecules, hence the observed reduction in sodium count. In this study, atomic coordinates for dummy oxygen atoms, representing H₂O molecules, were employed at three distinct supercage sites (SI S9) and the model proposed by Agostini et al.³¹ was applied. However, the MLP NPT-MD trajectories at high water content indicated water presence also in the sodalite cage (Fig. 4), a site not considered in the Rietveld refinement model. This indicates that the discrepancies in relative occupancies at low temperatures between experimental data and simulations can be attributed to two factors. First, the misidentification of water molecules in the sodalite cage as sodium cations during XRD data postprocessing. Second, the absence of sodium cations located in the supercage in the Rietveld model, which were detected during MLP MD simulations. This demonstrates the necessity of a synergistic experimental–theoretical approach to obtain a comprehensive understanding of the system.

Fig. 4 Snapshots focusing on the same sodalite cage and D6R, capturing the presence of water within the sodalite cage, even at high temperature (300 °C). Each snapshot was taken at 385 ps during the MLP MD simulation. The presence of water in the sodalite cage was not accounted for in the Rietveld refinement model, resulting in discrepancies at low temperatures between the sodium occupancies obtained from experimental and computational data.

While the initial relative sodium occupancies differ between experimental and computational results, evidence of sodium cation migration was detected in both the refined XRD experimental data and MLP NPT-MD simulations. As dehydration progressed, more cations occupied the D6R (site I) and the supercage close to the hexagonal window to the sodalite cage (site II), with fewer remaining in the sodalite cage (sites I′ and II′). This trend was apparent in both experimental and computational results. Based upon the agreement of both analyses, as water molecules begin to exit the framework, sodium cation migration is triggered. Sodium cations prefer to be tightly coordinated; thus, cations that were previously stabilized by ion–dipole interactions with water in the supercage migrate to more confined environments near the zeolite oxygen atoms, occupying site II.³⁶ In general, the transport of sodium cations through six-membered rings is challenging, as demonstrated by the established practice of partially breaking apart the as-synthesized zeolite structure through mild steaming to facilitate sodium removal by ion exchange.^37,38 However, the migration from site I′ to I upon water removal is facilitated by the increase in coordination number afforded by the tight space in the double 6-ring. It is therefore essential to examine the volumetric changes of this material to fully comprehend the behavior of the Na-Y zeolites during dehydration.

Unit cell parameter changes upon dehydration

The evolution of the unit cell parameter of zeolite Na-Y as a function of temperature is plotted as grey points in Fig. 5, with data obtained from both XRD refinement and MLP NPT-MD simulations. Both methods reveal identical trends and the agreement between experiment and theory is excellent. Initially, as the number of water molecules decreases, a reduction in the unit cell parameter is observed, which is consistent with the intuitive contraction upon the removal of molecules from the system. However, when the number of water molecules is roughly halved, an unexpected behavior is observed wherein the unit cell parameter begins to increase.

Fig. 5 Evolution of the unit cell parameter (grey) and number of adsorbed water molecules (blue) as a function of temperature. Obtained using (A) parametric Rietveld refinement, (B) MLP NPT-MD.

In the experimentally obtained values, this increase starts at temperatures above 100 °C, with the cell parameter of the fully dehydrated sample surpassing the values observed when the zeolite was hydrated at 50 °C. In the computational system, the cell parameter increases at temperatures above 200 °C; although, the magnitude of this increase is smaller than the growth derived from experimental data. Changes in the unit cell parameter in faujasite zeolite due to increasing temperatures have been reported previously. However, to our knowledge, there have been no documented reports detailing an increase in the cell parameter. In early work by Catlow's group,³⁹ changes in the lattice parameter upon heating of purely siliceous faujasite and Na-X zeolites were assessed using high-resolution powder XRD. Their study demonstrated shrinkage in the zeolite cell parameter upon heating regardless of the aluminum or sodium content but did not extend to include the response of the Na-Y zeolite.

To explore whether cations actively influence volumetric changes in zeolite Na-Y, several MLP NPT-MD simulations of the fully siliceous zeolite Y framework were conducted, with simulation details provided in the SI S10. Changes in the unit cell parameter as a function of temperature for both Na-Y and siliceous faujasite are illustrated in Fig. 6. In contrast to the increase observed in Na-Y, the siliceous model (beige points) demonstrates continuous decrease in the unit cell volume throughout the temperature increase. This behavior is consistent with previously reported XRD studies on siliceous faujasite,^39,40 thereby validating the MLPs' capability to accurately describe atomic interactions of siliceous faujasite with water at increasing temperatures. These computational results indicate that the observed changes in the unit cell parameter for Na-Y are primarily driven by the reduction in the number of water molecules and are significantly influenced by the presence of sodium cations. As water molecules exit the zeolite framework, there is an associated contraction in volume. This water removal induces cationic mobility, particularly their migration towards the previously unoccupied site I.

Fig. 6 Comparison of changes in the unit cell parameter as a function of temperature for Na-Y (grey) and siliceous faujasite (beige) with the same water loading.

The relationship between changes in the unit cell parameter and the relative occupancy of site I by sodium cations with respect to temperature is depicted in Fig. 7. As the occupancy of Na I increases, transitioning from an initially empty D6R state to almost half-occupied, the cell parameter of the zeolite framework simultaneously expands. This phenomenon is observed through both XRD refinement and MLP NPT-MD simulations. We therefore suggest that, as sodium cations migrate to occupy the double six-ring (site I), the structural elements of the zeolite framework adjust to accommodate them. The attraction towards oxygen atoms and repulsion from silicon or aluminum atoms alters the T–O–T angle within the double six-membered ring. The resulting strain in the T–O–T angle is reflected in the expansion of the unit cell parameter of zeolite Na-Y. These findings emphasize the critical role sodium cations play in the structural dynamics of zeolite Na-Y during dehydration.

Fig. 7 Evolution of the unit cell parameter (grey) and relative occupancy of site I by sodium (pink) as a function of temperature. Obtained using (A) parametric Rietveld refinement, (B) MLP NPT-MD.

Conclusions

The dehydration process and consequent cation migration dynamics in zeolite Na₅₄Y were investigated by parametric Rietveld refinement of in situ XRD data and effectively combined with Grand Canonical Monte Carlo simulations and Machine Learning Potential Molecular Dynamics simulations. Detailed structural insights were provided by the parametric Rietveld analysis, allowing us to monitor the evolution of the lattice parameters and occupancies of guest species as the temperature increased. Additionally, the role of cation migration in driving the observed structural changes during dehydration was successfully predicted and resolved through the integration of GCMC and MLP MD simulations. The predictive capabilities of MLP MD complemented the experimental findings.

Experimental and computational results reveal a sequential pattern of cell parameter changes during dehydration. Initially, water removal causes unit cell contraction. Sodium cations, which prefer coordination, exhibit increased mobility as water exits the zeolite. Cations that were initially surrounded by water in the supercage migrate closer to the hexagonal window of the sodalite cage (site II), driven by the need to coordinate with the zeolite atoms. Continued water removal triggers sodium migration from the sodalite cage (sites I′ and II′) towards the double six-ring (site I). As cations occupy site I, the zeolite framework begins to expand, with further site I occupancy leading to additional expansion.

The nuanced image of zeolite behavior in response to dehydration demonstrates the value of integrating parametric Rietveld refinement and MLP MD simulations in uncovering such intricate dynamics in complex environments. This lays the foundation for further work which could involve an exploration, though computational means, whether the observed behavior during dehydration is typical of faujasite materials. Such investigation would require the expansion of training dataset to include a broader range of Si/Al ratios and retraining a new MLP, tailored for different faujasite compositions.

As a future perspective, it is noteworthy that new methods have been developed to train an MLP capable not only of performing molecular dynamics but also Monte Carlo simulations, considering both flexible frameworks and the quantum accuracy of intermolecular interactions.⁴¹ Additionally, the novel computational approach applied here can be extended to analyze interactions between complex materials, like metal–organic frameworks, and adsorbed guest molecules with high accuracy, offering enhanced computational efficiency relative to DFT MD simulations.

Author contributions

Agnieszka Seremak: methodology, investigation, formal analysis, data curation, writing – original draft, writing – review & editing, visualization. Ruben Goeminne: methodology, software, resources, writing – review & editing. Izar Capel Berdiell: investigation, data curation, formal analysis, writing – review & editing. Lars F. Lundegaard: conceptualization, investigation, project administrator, writing – review & editing. Veronique Van Speybroeck: conceptualization, resources, supervision, writing – review & editing. Stian Svelle: conceptualization, investigation, funding acquisition, supervision, writing – review & editing.

Conflicts of interest

There are no conflicts of interest to declare.

Data availability

Input, output files and necessary scripts are publicly available via Zenodo at https://doi.org/10.5281/zenodo.15594037.

Supplementary information: details of Rietveld refinement model, EDX elemental analysis and ²⁷Al solid-state NMR, computational details of GCMC and MLPs training and computational details of MD simulations of siliceous Y zeolite. See DOI: https://doi.org/10.1039/d5ta04948b.

Acknowledgements

This work was supported by the European Union's Horizon 2020 research and innovation program under the Marie Skłodowska-Curie grant agreement no. 945371 and by The Research Council of Norway via the TomoCAT researcher project (project no. 301619). GCMC simulations were performed on resources provided by Sigma2 – the National Infrastructure for High-Performance Computing and Data Storage in Norway (Project nn4683k). We acknowledge the EuroHPC Joint Undertaking for awarding this project access to the EuroHPC supercomputer LUMI, hosted by CSC (Finland) and the LUMI consortium through a EuroHPC Regular Access call which allowed us to perform the Machine Learning Potential Training (Grant 465000259 and 465000239). MLP MD simulations were performed using computational resourced provided by the VSC (Flemish Supercomputer Center), funded by Ghent University. VVS acknowledges the Research Board of Ghent University (BOF). The NMR component of this work was supported by the Research Council of Norway through the Norwegian NMR Platform 1 & 2 (226244, 322373). Finally, we would like to thank Mads Ry Jørgensen, David Wragg, Nicolai Haaber Junge and Robin Christensen for their help at the DanMAX beamline, Bjørnar Arstad for NMR analysis, Sander Vandenhaute and Massimo Bocus for their assistance with the MLP simulations, Torstein Fjermestad and Pablo Beato for useful discussions, and Nathan Campbell for proofreading and providing constructive criticism of the manuscript.

References

1. V. Van Speybroeck, K. Hemelsoet, L. Joos, M. Waroquier, R. G. Bell and C. R. A. Catlow, Advances in Theory and Their Application within the Field of Zeolite Chemistry, Chem. Soc. Rev., 2015, 44(20), 7044–7111,  10.1039/C5CS00029G.
2. Database of Zeolite Structures, https://europe.iza-structure.org/IZA-SC/ftc_table.php, accessed 12 September 2024.
3. W. H. Baur, On the Cation and Water Positions in Faujasite, Am. Mineral., 1964, 49(5–6), 697–704.
4. W. J. Mortier, E. Van den Bossche and J. B. Uytterhoeven, Influence of the Temperature and Water Adsorption on the Cation Location in NaY Zeolites, Zeolites, 1984, 4(1), 41–44,  DOI:10.1016/0144-2449(84)90071-X.
5. C. Beauvais, A. Boutin and A. H. Fuchs, A Numerical Evidence for Nonframework Cation Redistribution upon Water Adsorption in Faujasite Zeolite, ChemPhysChem, 2004, 5(11), 1791–1793,  DOI:10.1002/cphc.200400195.
6. C. Abrioux, B. Coasne, G. Maurin, F. Henn, M. Jeffroy and A. Boutin, Cation Behavior in Faujasite Zeolites upon Water Adsorption: A Combination of Monte Carlo and Molecular Dynamics Simulations, J. Phys. Chem. C, 2009, 113(24), 10696–10705,  DOI:10.1021/jp902274t.
7. W. Louisfrema, B. Rotenberg, F. Porcher, J. L. Paillaud, P. Massiani and A. Boutin, Cation Redistribution upon Dehydration of Na₅₈Y Faujasite Zeolite: A Joint Neutron Diffraction and Molecular Simulation Study, Mol. Simul., 2015, 41(16–17), 1371–1378,  DOI:10.1080/08927022.2015.1027889.
8. A. D. Lella, N. Desbiens, A. Boutin, I. Demachy, P. Ungerer, J.-P. Bellat and A. H. Fuchs, Molecular Simulation Studies of Water Physisorption in Zeolites, Phys. Chem. Chem. Phys., 2006, 8(46), 5396–5406,  10.1039/B610621H.
9. S. Buttefey, A. Boutin and A. H. Fuchs, Cation Distribution in Faujasite-Type Zeolites: A Test of Semi-Empirical Force Fields for Na Cations, Mol. Simul., 2002, 28(12), 1049–1062,  DOI:10.1080/0892702021000011070.
10. W. Louisfrema, J.-L. Paillaud, F. Porcher, E. Perrin, T. Onfroy, P. Massiani, A. Boutin and B. Rotenberg, Cation Migration and Structural Deformations upon Dehydration of Nickel-Exchanged NaY Zeolite: A Combined Neutron Diffraction and Monte Carlo Study, J. Phys. Chem. C, 2016, 120(32), 18115–18125,  DOI:10.1021/acs.jpcc.6b05657.
11. S. Ghojavand, E. Dib and S. Mintova, Flexibility in Zeolites: Origin, Limits, and Evaluation, Chem. Sci., 2023, 14(44), 12430–12446,  10.1039/D3SC03934J.
12. P. Beato, L. F. Lundegaard, S. Svelle and D. S. Wragg, Case Studies: Crystallography as a Tool for Studying Methanol Conversion in Zeolites, in Springer Handbook of Advanced Catalyst Characterization, ed. I. E. Wachs and M. A. Bañares, Springer International Publishing, Cham, 2023, pp. 541–563,  DOI:10.1007/978-3-031-07125-6_26.
13. G. W. Stinton and J. S. O. Evans, Parametric Rietveld Refinement, J. Appl. Crystallogr., 2007, 40(1), 87–95,  DOI:10.1107/S0021889806043275.
14. A. A. Coelho, TOPAS and TOPAS-Academic: An Optimization Program Integrating Computer Algebra and Crystallographic Objects Written in C++, J. Appl. Crystallogr., 2018, 51(1), 210–218,  DOI:10.1107/S1600576718000183.
15. D. Dubbeldam, S. Calero, D. E. Ellis and R. Q. Snurr, RASPA: Molecular Simulation Software for Adsorption and Diffusion in Flexible Nanoporous Materials, Mol. Simul., 2016, 42(2), 81–101,  DOI:10.1080/08927022.2015.1010082.
16. W. Loewenstein, The Distribution of Aluminum in the Tetrahedra of Silicates and Aluminates, Am. Mineral., 1954, 39(1–2), 92–96.
17. M. Erdős, D. F. Geerdink, A. Martin-Calvo, E. A. Pidko, L. J. P. van den Broeke, S. Calero, T. J. H. Vlugt and O. A. Moultos, In Silico Screening of Zeolites for High-Pressure Hydrogen Drying, ACS Appl. Mater. Interfaces, 2021, 13(7), 8383–8394,  DOI:10.1021/acsami.0c20892.
18. W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey and M. L. Klein, Comparison of Simple Potential Functions for Simulating Liquid Water, J. Chem. Phys., 1983, 79(2), 926–935,  DOI:10.1063/1.445869.
19. E. Jaramillo and S. M. Auerbach, New Force Field for Na Cations in Faujasite-Type Zeolites, J. Phys. Chem. B, 1999, 103(44), 9589–9594,  DOI:10.1021/jp991387z.
20. J. M. Castillo, J. Silvestre-Albero, F. Rodriguez-Reinoso, T. J. H. Vlugt and S. Calero, Water Adsorption in Hydrophilic Zeolites: Experiment and Simulation, Phys. Chem. Chem. Phys., 2013, 15(40), 17374–17382,  10.1039/C3CP52910J.
21. N. Jiang, M. Erdős, O. A. Moultos, R. Shang, T. J. H. Vlugt, S. G. J. Heijman and L. C. Rietveld, The Adsorption Mechanisms of Organic Micropollutants on High-Silica Zeolites Causing S-Shaped Adsorption Isotherms: An Experimental and Monte Carlo Simulation Study, Chem. Eng. J., 2020, 389, 123968,  DOI:10.1016/j.cej.2019.123968.
22. S. Vandenhaute, M. Cools-Ceuppens, S. DeKeyser, T. Verstraelen and V. Van Speybroeck, Machine Learning Potentials for Metal-Organic Frameworks Using an Incremental Learning Approach, npj Comput. Mater., 2023, 9(1), 1–8,  DOI:10.1038/s41524-023-00969-x.
23. I. Batatia, D. P. Kovács, G. N. C. Simm, C. Ortner and G. Csányi, Higher Order Equivariant Message Passing Neural Networks for Fast and Accurate Force Fields, arxiv, 2022 preprint, arXiv:2206.07697,  DOI:10.48550/arXiv.2206.07697.
24. T. D. Kühne, M. Iannuzzi, M. Del Ben, V. V. Rybkin, P. Seewald, F. Stein, T. Laino, R. Z. Khaliullin, O. Schütt, F. Schiffmann, D. Golze, J. Wilhelm, S. Chulkov, M. H. Bani-Hashemian, V. Weber, U. Borštnik, M. Taillefumier, A. S. Jakobovits, A. Lazzaro, H. Pabst, T. Müller, R. Schade, M. Guidon, S. Andermatt, N. Holmberg, G. K. Schenter, A. Hehn, A. Bussy, F. Belleflamme, G. Tabacchi, A. Glöß, M. Lass, I. Bethune, C. J. Mundy, C. Plessl, M. Watkins, J. VandeVondele, M. Krack and J. Hutter, CP2K: An Electronic Structure and Molecular Dynamics Software Package – Quickstep: Efficient and Accurate Electronic Structure Calculations, J. Chem. Phys., 2020, 152(19), 194103,  DOI:10.1063/5.0007045.
25. S. Goedecker, M. Teter and J. Hutter, Separable Dual-Space Gaussian Pseudopotentials, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54(3), 1703–1710,  DOI:10.1103/PhysRevB.54.1703.
26. J. P. Perdew, K. Burke and M. Ernzerhof, Generalized Gradient Approximation Made Simple, Phys. Rev. Lett., 1996, 77(18), 3865–3868,  DOI:10.1103/PhysRevLett.77.3865.
27. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu, J. Chem. Phys., 2010, 132(15), 154104,  DOI:10.1063/1.3382344.
28. G. Lippert, J. Hutter and M. Parrinello, A Hybrid Gaussian and Plane Wave Density Functional Scheme, Mol. Phys., 1997, 92(3), 477–488,  DOI:10.1080/002689797170220.
29. P. Eastman, R. Galvelis, R. P. Peláez, C. R. A. Abreu, S. E. Farr, E. Gallicchio, A. Gorenko, M. M. Henry, F. Hu, J. Huang, A. Krämer, J. Michel, J. A. Mitchell, V. S. Pande, J. P. Rodrigues, J. Rodriguez-Guerra, A. C. Simmonett, S. Singh, J. Swails, P. Turner, Y. Wang, I. Zhang, J. D. Chodera, G. De Fabritiis and T. E. Markland, OpenMM 8: Molecular Dynamics Simulation with Machine Learning Potentials, J. Phys. Chem. B, 2024, 128(1), 109–116,  DOI:10.1021/acs.jpcb.3c06662.
30. A. Erlebach, M. Šípka, I. Saha, P. Nachtigall, C. J. Heard and L. Grajciar, A Reactive Neural Network Framework for Water-Loaded Acidic Zeolites, Nat. Commun., 2024, 15(1), 4215,  DOI:10.1038/s41467-024-48609-2.
31. G. Agostini, C. Lamberti, L. Palin, M. Milanesio, N. Danilina, B. Xu, M. Janousch and J. A. Van Bokhoven, In Situ XAS and XRPD Parametric Rietveld Refinement to Understand Dealumination of Y Zeolite Catalyst, J. Am. Chem. Soc., 2010, 132(2), 667–678,  DOI:10.1021/ja907696h.
32. T. Frising and P. Leflaive, Extraframework Cation Distributions in X and Y Faujasite Zeolites: A Review, Microporous Mesoporous Mater., 2008, 114(1–3), 27–63,  DOI:10.1016/j.micromeso.2007.12.024.
33. A. Ghysels, S. L. C. Moors, K. Hemelsoet, K. De Wispelaere, M. Waroquier, G. Sastre and V. Van Speybroeck, Shape-Selective Diffusion of Olefins in 8-Ring Solid Acid Microporous Zeolites, J. Phys. Chem. C, 2015, 119(41), 23721–23734,  DOI:10.1021/acs.jpcc.5b06010.
34. R. J. Gowers; M. Linke; J. Barnoud; T. J. E. Reddy; M. N. Melo; S. L. Seyler; J. Domański; D. L. Dotson; S. Buchoux; I. M. Kenney and O. Beckstein, MDAnalysis: A Python Package for the Rapid Analysis of Molecular Dynamics Simulations, SciPy, 2016,  DOI:10.25080/Majora-629e541a-00e.
35. N. Michaud-Agrawal, E. J. Denning, T. B. Woolf and O. Beckstein, MDAnalysis: A Toolkit for the Analysis of Molecular Dynamics Simulations, J. Comput. Chem., 2011, 32(10), 2319–2327,  DOI:10.1002/jcc.21787.
36. E. Dooryhee, C. R. A. Catlow, J. W. Couves, P. J. Maddox, J. M. Thomas, G. N. Greaves, A. T. Steel and R. P. Townsend, A Study of Cation Environment and Movement during Dehydration and Reduction of Nickel-Exchanged Zeolite Y by x-Ray Absorption and Diffraction, J. Phys. Chem., 1991, 95(11), 4514–4521,  DOI:10.1021/j100164a062.
37. W. Lutz, Zeolite Y: Synthesis, Modification, and Properties—A Case Revisited, Adv. Mater. Sci. Eng., 2014, 2014(1), 724248,  DOI:10.1155/2014/724248.
38. K. Sato, Y. Nishimura, N. Matsubayashi, M. Imamura and H. Shimada, Structural Changes of Y Zeolites during Ion Exchange Treatment: Effects of Si/Al Ratio of the Starting NaY, Microporous Mesoporous Mater., 2003, 59(2), 133–146,  DOI:10.1016/S1387-1811(03)00305-6.
39. J. W. Couves, R. H. Jones, S. C. Parker, P. Tschaufeser and C. R. A. Catlow, Experimental Verification of a Predicted Negative Thermal Expansivity of Crystalline Zeolites, J. Phys.: Condens. Matter, 1993, 5(27), L329,  DOI:10.1088/0953-8984/5/27/001.
40. M. P. Attfield, Strong Negative Thermal Expansion in Siliceous Faujasite, Chem. Commun., 1998, 601–602,  10.1039/A707141H.
41. R. Goeminne and V. Van Speybroeck, Ab Initio Predictions of Adsorption in Flexible Metal–Organic Frameworks for Water Harvesting Applications, J. Am. Chem. Soc., 2025, 147(4), 3615–3630,  DOI:10.1021/jacs.4c15287.

---