Periodic modeling of zeolite Ti-LTA

Zeolite Linde Type A (LTA) was originally reported as an aluminosilicate material with formula $Na12[(AlO2)12(SiO2)12]$$⋅$$27H2O$⁠.^1,2 The periodic building unit (PerBU) of LTA comprises three-dimensionally interconnected $β$-cages [Fig. 1(a)] that generate an $α$-cage at the center of a cubic unit cell with Pm$3¯$m symmetry [Fig. 1(b)]. An alternative PerBU to the $β$-cage consists of double 4-membered rings [D4R, Fig. 1(c)].³ The pore system of LTA is characterized by eight-membered ring (8MR) apertures with a diameter smaller than 7 Å connecting $α$-cages along the three dimensions [Fig. 1(e)]. The application of LTA as a molecular sieve was identified early on following its discovery by Breck and collaborators, who reported how the ion exchange of $Na+$ by other cations with charges 1+ and 2+ was able to tailor the adsorption capacity and selectivity for a wide array of small molecules.¹ The relatively large size of the $α$-cages in LTA, combined with small apertures formed by 8MR, leads to attractive applications such as the formation and confinement of metal clusters, which are protected by the zeolite framework against poisoning, thereby allowing their selectivity towards the conversion of small molecules.⁴ 

Although zeolite LTA is mainly recognized by its high Al content (known as zeolite A), driving its applications as molecular sieve and sorbent,⁵ pure-silica LTA has been synthesized directly using self-assembled organic structure-directing agents (SDAs).⁶ Recently, high silica zeolite LTA, with Si/Al ratios of 11 and 16, has shown exceptional durability in the selective catalytic reduction of $NOx$ with $NH3$ (⁠$NH3$-SCR) under hydrothermal conditions at temperatures as high as 750 °C.⁷ Furthermore, aluminosilicate LTA with Si/Al ratios between 12 and 42 has shown promising results related to the selectivity of $C5$ and $C6$ products in the methanol-to-olefins (MTO) reaction.⁸ Boal et al. have also proposed the synthesis of titanosilicate LTA, suggesting that this material may be suitable for the oxidation of small molecules, taking advantage of the 8MR pore openings, the relatively large cavities, and the presence of intra-framework sites with Lewis acidity.⁸ Titanium-substituted zeolites have been used previously as catalysts for the epoxidation and oxidation of different organic molecules, such as phenol, ethylene, and bulkier olefins.^9–12 

A number of theoretical reports have analyzed Ti substitution in different zeolites, such as silicalite, using cluster models and periodic boundary conditions (PBC).^13–17 However, a systematic study of Ti substitution in zeolite LTA is still missing. The recent synthesis of high silica, Ti-substituted zeolite LTA⁸ provides the rationale for theoretical research to elucidate the local structure, electronic properties, and distribution of Ti(iv) within the framework of LTA, and this fundamental information may be of use in the design of future improved catalysts.

In the present work, we have examined single and double substitutions of Si by Ti in the T-sites of zeolite LTA. The combination of density functional theory (DFT) calculations with atomistic simulations allowed us to explore in detail the electronic structure of these systems whilst validating the use of low Si/Ti ratios against more realistic, lower Ti concentrations. We also considered the interaction of water with the Lewis acid site and its possible dissociation over the active site. In addition, the infrared (IR) spectra of relevant structures were simulated, which allowed us to propose a theoretical benchmark that may help us to interpret the experimental analysis of the zeolite structure.

We have performed the structural and electronic analysis of pure-silica and Ti-substituted zeolite LTA using DFT calculations as implemented in the Vienna Ab initio Simulation Package (VASP).^18–21 In addition, we have included atomistic simulations that employ the Born model of ionic solids,²² to investigate a more realistic Si/Ti ratio to complement the DFT results.

The DFT calculations were carried out under periodic boundary conditions, employing a basis set of plane waves to explicitly treat the valence electrons, whilst the interaction with the inner core of the atom and nodal features were described by the projector-augmented-wave (PAW)^23,24 method. The accuracy of the basis set was guaranteed by a kinetic energy cutoff of 550 eV. We have used the generalized gradient approximation (GGA), as proposed by Perdew, Burke, and Ernzerhof (PBE), to account for the exchange and correlation contributions to the electronic energy.²⁵ The effect of long-range dispersion forces was considered using Grimme’s atomic pairwise method that includes three-body terms to avoid the overbinding of only two-body methods.^26–28 Furthermore, we have chosen the Becke-Johnson damping function to avoid the double counting of the correlation effects at short distances.²⁹ The Brillouin zone was sampled with a k-point mesh of 3 $×$ 3 $×$ 3, generated with the Monkhorst–Pack scheme.³⁰ We have used the linear tetrahedron method, with Blöchl’s correction, to account for the occupation of the electronic bands and the integration over the k-point mesh during the structural optimization.³¹ This method was employed because the systems analyzed here are insulators, with all bands fully occupied or fully empty, thus avoiding the method’s drawback of not being variational with respect to partial occupancies. The thresholds for the electronic and ionic iterations were set at 10⁻⁶ eV and 0.01 eV/Å, respectively, ensuring proper convergence.

The GGA vibrational frequencies of the zeolite framework were obtained by the diagonalization of the Hessian matrix computed by density functional perturbation theory.³² The signal intensity of each normal mode was estimated using the following equation:

where $τi$ is the Born effective charge tensor of the ith atom, calculated by density functional perturbation theory and $νiα$ is the atomic displacement in the $α$th eigenvector divided by the square root of its mass. The IR spectra were simulated by widening each signal with a Gaussian function.

The interatomic potential (IP) calculations were performed using the General Utility Lattice Program (GULP),^33,34 where the energy of the system comprises Coulombic contributions, calculated according to the Ewald summation method,³⁵ and short range repulsions and dispersion forces described by Buckingham potentials.³⁶ For the potential model of the $SiO2$ units that form the zeolite framework, we have used the parametrization proposed by Sanders et al.,³⁷ which was initially developed for $α$-quartz but has shown satisfactory transferability to zeolites.^38–41 The presence of Ti within the zeolite framework was described by the potentials reported by Jentys and Catlow.¹⁴ The IP vibrational frequencies were computed after the diagonalization of the Hessian matrix, which was obtained by the finite difference method that keeps the system within the harmonic approximation. The signal intensity of each normal mode was evaluated by the following equation:

where q_i is the charge of the ith atom and $νiα$ is the atomic displacement in the $α$th eigenvector divided by the square root of its mass.

Within the GGA, we optimized the LTA unit cell by using a set of fixed-volume calculations, where only the atomic positions were allowed to relax. The equilibrium volume was then evaluated by fitting the Birch-Murnaghan equation of state to the correlation between the lattice energy versus the lattice volume.⁴² This procedure also gave the isotropic bulk modulus as a fitted parameter. The optimization of the cell and the atomic positions with the IP method took place at constant pressure, thus permitting the simultaneous relaxation of the cell volume and the atomic coordinates in a single calculation. The unit cell of zeolite LTA has 24 T-sites, all of them equivalent by an operation of symmetry within the Pm$3¯$m space group.

Two different cells were used for the calculations. The first one, the 1 $×$ 1 $×$ 1 unit cell containing 24 tetrahedral sites (or equivalently: 24 $TO2$ units), was used for the GGA study of pure-silica and Ti-substituted zeolite LTA. For this cell size, the single Ti substitution gave a Si/Ti ratio of 23 whilst the double Ti substitution decreased it to 11. These values are at least 4 to 5 times smaller than those reported experimentally.⁸ However, a bigger cell would have been computationally extremely demanding within the electronic structure calculations. Therefore, we have complemented the GGA study with IP simulations, which allow the examination of larger systems using much less computational resource. Both the GGA and the IP methods have shown to be suitably accurate to analyze the thermodynamic and physical properties of zeolites.^41,43–45

The second type of cell was formed by the 2 $×$ 2 $×$ 2 expansion of the 24-$TO2$ unit cell. As such, the number of $TO2$ units increased to 192, with a Si/Ti ratio of 191 and 95 for single and double substitutions, respectively, much closer to the experimental value of 100.⁸ All non-equivalent double substitutions of Ti in the zeolite framework were generated using the Site Occupancy Disorder (SOD) program, obtaining 12 structures for the 1 $×$ 1 $×$ 1 cell and 91 structures for the 2 $×$ 2 $×$ 2 supercell.⁴⁶ In order to simplify the calculations, the presence of the SDA and its effect on the Ti substitution were not considered.

The structural visualizations presented in this work were generated by the Visualization for Electronic and Structural Analysis (VESTA 3).⁴⁷ 

The GGA elastic constants were obtained following the protocol of deformations proposed by Ravindran et al., where the volume is held constant during the orthorhombic and monoclinic shear distortions.⁴⁸ Owing to the cubic symmetry of zeolite LTA, the elastic constants that we needed to evaluate reduced to C₁₁, C₁₂, and C₄₄ only. The elastic constants were obtained by fitting a second order polynomial equation to the correlation between the lattice energy and the applied strain,⁴⁸ 

where $EDn$ is the set of lattice energies obtained after applying the deformation D_n suggested by Ravindran and collaborators,⁴⁸ E₀ is the energy at the equilibrium volume V₀, $τi$ is an element in the stress tensor, and $δ$ is the applied strain (between −3% and 3% to minimize the contribution of higher order terms in the Taylor expansion).

The bulk modulus (B) and shear modulus (G) were obtained following the Voigt approximation, where the external strain is equated to the uniform strain in a polycrystalline ensemble.⁴⁹ The values of B and G correspond to the resistance of the material to fracture and plastic deformation, respectively. The Voigt equations use the elastic constants for the calculation of B and G,

In addition, Young’s modulus (Y) and Poisson’s ratio (⁠$σ$⁠), which give information regarding the rigidity of the crystal, and the shear anisotropic factor (A), which measures the level of anisotropy of the material, were calculated from the following equations:

In the case of the IP method, the elastic parameters were directly provided by the GULP following Voigt’s approximation.

After structural relaxation of pure-silica LTA within the GGA, we obtained a cubic unit cell parameter of 11.950 Å, which was only 0.8% longer than the reported value by Boal et al. (see Table I).⁸ The IP method decreased the difference even further, finishing with a value after optimization of 11.847 Å, only slightly shorter than the experimental parameter. The small discrepancy between the cell parameters evaluated by the GGA and IP was because of longer Si–O distances for the GGA compared to the IP method, by an average of 0.022 Å. A similar outcome has been reported for the optimization of the unit cell of zeolite MFI, where the GGA average value was 0.112 Å longer than the IP prediction.⁴¹ 

We calculated a bulk modulus B of 42 GPa by fitting the GGA-calculated energies to the Birch-Murnaghan equation,⁴² decreasing to 35 GPa when employing the GGA-calculated elastic constants within the Voigt equation (see Table I).⁴⁹ The IP method led to a value for B of 74 GPa. These results followed the same trend of zeolite MFI, where a B of 18 GPa was evaluated by GGA, rising to 54 GPa with IP,⁴¹ with the GGA value being much closer to the experimental 18.2 GPa.⁵¹ The validation of our results is more difficult in the case of pure-silica zeolite LTA because most of the experimental measurements provided in the literature come from zeolite-A, which has the same framework of LTA but with a Si/Al ratio of 1 and loaded with 6–12 extra-framework cations per unit cell to compensate for the negative charge of the framework. For instance, Arletti et al. have reported a value of 22.1 GPa for zeolite-A,⁵⁰ whilst Niwa et al. have measured values between 25.3 and 34.4 GPa depending on the pressure-transmitting medium.⁵² Our calculated B is larger than the experimentally measured values, but this may be attributed to the modifications in bonding strengths when the Si/Al ratio changes from 1 to infinity, which makes the structure more resistant to deformations; an experimental benchmark using pure-silica LTA is necessary to validate our results.

The values of the calculated shear moduli G, presented in Table I, of 16 and 23 GPa for GGA and IP, respectively, mean that zeolite LTA has a relatively low resistance to plastic deformation.^53,54 Huang and Havenga have actually reported the reversibility of the pressure-induced amorphization of zeolite-A in the range of 15–25 kbar, although they showed that the presence of extra-framework cations is essential for recovering the crystalline structure after the deformation.⁵⁵ G calculated in this study considered a pure-silica LTA without extra-framework cations, which should therefore make its G value larger than the one for zeolite-A.

At the same time, the ratio B/G measures the level of ductility (if high) or brittleness (if low) of a material, with the turning point between ductility and brittleness at approximately 1.75.⁵³ We obtained values of 2.2 and 3.2 for B/G, using the GGA and IP method, respectively. Both types of calculations depict zeolite LTA as ductile, although the larger value of IP should be due to the overestimation of B.

The Poisson ratio $σ$ is a measure of the ionic or covalent character of the bonds within the material. Values smaller than 0.1 are attributed to covalent bonding, whilst a ratio above 0.25 suggests ionic interactions.⁵⁶ The GGA and IP led to similar values of $σ$ of 0.30 and 0.37, respectively, with both methods giving values above 0.25. This result also validates the use of the Born model of ionic solids to describe the largely ionic zeolite structure.²² 

We have also calculated the shear anisotropic factor A that evaluates the level of anisotropy of the crystal for values smaller or larger than unity. In the present study, GGA had an A of 1.6, reducing to 1.2 in the case of IP. Accordingly, pure-silica LTA should show a significant degree of anisotropy according to the GGA. However, a previous report has shown that A may be overestimated by the GGA methodology, where a theoretical value of 1.34 is predicted for magnetite (⁠$Fe3O4$⁠) against an experimental measurement of 1.13.⁵⁷ 

We have studied the stability and distribution of single and double substitutions of Si by Ti using Boltzmann statistics, a method that has been used in the past to analyze single and double substitutions of aluminum in zeolite MFI.⁵⁸ The 1 $×$ 1 $×$ 1 unit cell was used for the GGA calculations, whilst an expanded 2 $×$ 2 $×$ 2 supercell was used with the IP method in order to complement the GGA results. Figure 2 shows the Boltzmann distributions at 125 °C, which is the synthesis temperature used by Boal et al. to obtain the titanosilicate.⁸ 

The GGA calculations identified five double Ti-substituted structures with significant probabilities out of the twelve non-equivalent substitutions in the 1 $×$ 1 $×$ 1 unit cell. Owing to its small size, only the arrangements $Ti–Ti$⁠, $Ti–Si–Ti$⁠, and $Ti–Si–Si–Ti$ were obtained for this cell, where configuration 2_a in Fig. 2(a) had the highest prevalence (the number in the label specifies the number of bridging Si atoms placed between the two Ti atoms). In 2_a, the two Ti atoms were evenly placed at opposite sides of the six-membered ring (6MR) of zeolite LTA [see Fig. 2(a)], separated by a $Ti–Ti$ distance of 6.2 Å. This ordering is similar to the one reported by Bare et al., who used EXAFS to study the substitution of Sn in $β$-zeolite and highlighted the non-random substitution of Sn in the zeolite framework.⁵⁹ We did not consider the presence of the SDA in this work; therefore, the stability of the 6MR substitution may be related not only to the employed SDA but also to the capacity of this configuration to minimize the framework deformations. Other important structures placed both Ti atoms in the same D4R, labeled as 0_a and 1_a in Fig. 2(a), with $Ti–Ti$ distances of 3.2 and 4.6 Å, respectively, and in two different D4Rs, labeled as 1_b and 2_b, with $Ti–Ti$ distances of 5.3 and 7.2 Å, respectively.

We generated two different Boltzmann distributions for the IP method. The first one, shown in Fig. 2(b), only considered the substitutions within the same unit cell, despite using a 2 $×$ 2 $×$ 2 supercell, to allow a direct comparison with the GGA results. The second distribution used the full spectrum of double substitutions (91 in total) and is shown in Fig. 2(c). A close agreement between the IP and GGA profiles can be noted from the comparison of Figs. 2(a) and 2(b), regardless of the difference in the level of theory and cell size. Again, 2_a was by far the most stable structure, although the intensity of 0_a decreased to almost zero and an additional configuration labeled as 2_c arose [see Fig. 2(b)]. The full spectrum of double substitution shown in Fig. 2(c) confirmed once more that 2_a is the most stable arrangement when two Ti are placed in the zeolite framework. In Fig. 2(c), 2_a was compared against structures with $Ti–Ti$ separations of more than 15 Å, which can be considered as two fully isolated species.

In the present work, we are not primarily concerned about the absolute stability of each type of substitution, but instead we have compared the relative stability of the double substitution against the single substitution by using the following equation:

where $ETinLTA$ is the energy of the n Ti-substituted zeolite, $ELTA$ is the energy of the pure-silica zeolite, and $ETi1LTA$ is the energy of the single Ti-substituted zeolite. The derivation of Eq. (11) from the absolute stabilities is shown in the supplementary material. This equation can be interpreted as the energetics of a process that moves n Ti atoms from n single substituted unit cells to a sole unit cell, leaving the other (n − 1) as pure-silica frameworks. By using Eq. (11) with the GGA energies and n = 2, the calculated values of $Δ2$ were −11, −9, −8, −18, and −7 kJ/mol for configurations 0_a, 1_a, 1_b, 2_a, and 2_b, respectively. The other seven substitutions yielded positive values. This means that the double Ti substitutions represented in Fig. 2(a) are more favorable than the single substitution per 1 $×$ 1 $×$ 1 unit cell. The IP method showed equivalent outcomes. These results are different from a similar analysis by Yang et al., who studied Ti substitution in $β$-zeolite by DFT and found that the single substitution is more stable than the double substitution of Ti in the 6MR.⁶⁰ Nevertheless, this apparent discrepancy may be related to the topology of the 6MR: the ring is planar in LTA whilst it is bent in $β$-zeolite, which translates into different capacities to support framework distortions. We considered pure-silica LTA, single Ti-substituted LTA, and double Ti-substituted LTA with configurations 0_a and 2_a, as fair representations of the structural diversity; hence, we have only used these structures to continue the analysis in this work.

The methodology developed from Eq. (11) can be used to determine the limit of Ti aggregation in zeolite LTA by adding increasing amounts of Ti atoms in the structure and evaluating the resulting stability of the n-substituted framework against the (n − 1) substitution. However, this analysis is beyond the scope of the present work and will be discussed in future publications.

The $Si–O$ bond length and the $O–Si–O$ and $Si–O–Si$ angles of Ti-substituted zeolite LTA did not manifest significant modifications when compared with the pure-silica structure (see Table S1 of supplementary material). The GGA- and IP-calculated lengths of the $Ti–O$ bonds were longer than the $Si–O$ bonds by approximately 0.2 Å. The GGA-calculated $Ti–O–Ti$ angle of the 0_a structure was 14° to 24° smaller than the average $Si–O–Si$ and $Si–O–Ti$ angles. This reduction was due to the simultaneous increase in length of the two $Ti–O$ bonds, which forced the bridging O atom away from the imaginary line that connects both Ti atoms, in consequence reducing the $Ti–O–Ti$ angle. In contrast, the IP method did not show marked differences among the $Si–O–Si$⁠, $Si–O–Ti$⁠, and $Ti–O–Ti$ angles (see Table S1 of the supplementary material). At the same time, upon Ti substitution, the difference between the lowest unoccupied band (LUB) and the Fermi energy decreased (see Table S2 of the supplementary material). In addition, the electronic charge density associated with the LUB was projected out of the D4Rs into the cavities following the Ti substitution, as shown in Fig. S1 of the supplementary material. These results point to the Lewis acid character that is induced in the structure by the Ti atoms.

Infrared (IR) spectroscopy is a powerful experimental technique that provides structural information specific to each material. The theoretical simulation of the IR spectrum enables the correlation of variations in the spectrum profile with explicit modifications in the material structure without undesired interferences, which then may complement and help us to understand the experimental measurements.

We have simulated the theoretical IR spectrum of pure-silica LTA and validated it against the equivalent experimental report by Huang and Jiang.⁶¹ Figure 3 shows the simulated IR spectra using the GGA and the IP method, with the vibrational frequencies of the main bands listed in Table II. The theoretical IR spectra of pure-silica LTA present four main bands for both the IP method and the GGA, resembling the experimental spectrum within the limits of the theoretical approximations.⁶¹ The most intense signal, located at 928 cm⁻¹ for IP and 1073 cm⁻¹ for GGA, corresponds to the asymmetric $Si–O$ stretching, which has an experimental value of 1096 cm⁻¹,⁶¹ i.e., very close to the GGA value. The other three bands are located within the range 650–400 cm⁻¹, reproducing the experimental signals with deviations of −31 to +6 cm⁻¹ for both the IP method and the GGA. Within this range, the experimental bands at 621 and 477 cm⁻¹ are related to the D4R, whilst the signal at 407 cm⁻¹ is associated with the 6MR.⁶¹ A close inspection of the calculated vibrations in Table II for pure-silica LTA indicates a regular underestimation of the frequencies compared with the experiment, a trend that has been noted previously by Huang and Jiang.⁶¹ 

The GGA was better than the IP method in reproducing the experimental spectrum of pure-silica LTA; therefore, we have only analyzed single and double Ti substitutions by GGA-simulated IR spectra. In addition, the IP method manifested a lack of sensibility towards the Ti substitution, showing practically no variations when one or two Ti atoms were included within the zeolite framework.

The single Ti substitution causes a red-shift in the band of the asymmetric $Si–O$ stretching by 6 cm⁻¹, which coincides with the softening of this band observed in the experiment when Si atoms are replaced by Al, shifting by almost 90 cm⁻¹ from pure-silica LTA to zeolite A.⁶¹ The stretching for $Ti–O$ appears at 978 cm⁻¹ with lesser shoulders at 951 and 920 cm⁻¹. This signal is shifted with respect to the asymmetric $Si–O$ stretch by 89 cm⁻¹, corresponding to the weaker and longer bond length of $Ti–O$ when compared to $Si–O$⁠. The rest of the spectrum does not show significant variations after the single substitution.

The double substitution represented by structure 0_a in Fig. 2 enhances the intensity of the asymmetric $Ti–O$ stretching as expected and produces an additional weak band at 755 cm⁻¹. The signal at 755 cm⁻¹ is related to the stretching of the $Ti–O$ bond associated with the O atom bridging the two Ti atoms. As described above, the angle $Ti–O–Ti$ in 0_a decreases when compared with the pure-silica value, moving the bridging O atom further away from the axis linking the two Ti atoms, which, combined with longer $Ti–O$ bonds, makes the $Ti–O(–Ti)$ stretching 217 cm⁻¹ weaker than $Ti–O(–Si)$⁠. The band at 755 cm⁻¹ disappears again in the spectrum of structure 2_a, where the two Ti atoms are separated by two Si atoms, in agreement with our interpretation.

The spectrum of 2_a presents weak peaks above the asymmetric Si–O band at 1168 and 1143 cm⁻¹. These signals are related to $Si–O(–Si)$ stretching modes associated with the D4R where the Ti atoms are substituted. Equivalent signals tend to appear in the spectrum of the single Ti-substituted LTA, but their intensities are very weak, almost negligible (see Fig. 3). However, the spectrum of structure 0_a lacks those bands, which may help us to differentiate structures $Ti–(Si)0–Ti$ from those with a configuration $Ti–(Si)n>0–Ti$⁠. Additionally, narrower $Ti–O–Si$ angles cause softening of the $Ti–O(–Si)$ stretching by 68 cm⁻¹, producing a very weak band at 905 cm⁻¹ in the spectrum of 2_a. For instance, while the $Ti–O–Si$ angle is approximately 165° for the band at 973 cm⁻¹, it is reduced to 138° for the signal at 905 cm⁻¹.

The Lewis acidity can be probed by studying the adsorption strength of small molecules at the active sites.^13,60,62 We have chosen water as a probe molecule owing to its presence during the synthesis of the zeolite,⁸ which may also influence the distribution of Ti atoms along the framework positions. Pure-silica LTA, single Ti-substituted LTA, and configurations 0_a and 2_a were used to analyze the adsorption of water. Table III lists the highest adsorption energies of a water molecule in the $α$-cage and the $β$-cage (see Fig. 4), including the Van der Waals (VdW) contribution to these energies, for pure-silica LTA and single and double Ti substitutions.

The adsorption of water in pure-silica LTA released energies of 21 and 32 kJ/mol for loadings in the $α$-cage and the $β$-cage, respectively, where the VdW contributions accounted for 66% and 68%, respectively. The water molecule did not have a preferred adsorption site, as expected from its interaction with a completely siliceous framework.⁶⁰ However, the reduced space of the $β$-cage, compared to the $α$-cage, produced stronger VdW and H-bond interactions between the water molecule and the surrounding framework atoms.

As expected, the Ti substitution increased the interaction strength of water with the zeolite. The exothermic energies for adsorption in the $α$-cage were estimated at 37–43 kJ/mol for single and double Ti substitutions, increasing to 48–52 kJ/mol for loadings in the $β$-cage (see Table III). The VdW contribution to the adsorption energy decreased from 66% to 68% for pure-silica LTA to 41%–48% for the Ti-substituted structures, which is anticipated for water-Ti interactions with more of a chemisorption character. The Ti–O_water distances ranged from 2.363 to 2.424 Å, with the $β$-cage producing smaller distances than the $α$-cage by 0.003–0.058 Å.

We have also analyzed the possibility of water dissociation on the Lewis site by breaking one of the $O–H$ bonds of water and transferring the proton to the $TiO4$ tetrahedron. However, in most cases, either the proton returned to the OH group binding the Ti atom, re-forming a water molecule, or the structure was more than 60 kJ/mol less stable than the non-dissociated system. Water dissociation was more likely when a two-membered ring (2MR) was formed, where both OH groups bound simultaneously the Ti atom and its nearest-neighbor Si atom from a different D4R, as shown in Fig. 5(b). We did not observe the formation of 2MRs when the Ti and the Si atoms belonged to the same D4R. We attribute this outcome to the inability of a single D4R to accommodate the distortions provoked by the formation of the 2MR. After observing the importance of the 2MR formation, we also optimized the pure-silica 2MR [Fig. 5(a)] and the double Ti-substituted 2MR [Fig. 5(c)] in order to analyze the dependence of the 2MR stability on the extent of the Ti substitution (for relative energies of the process, see Table IV).

The pure-silica 2MR (labeled as 2Si-2MR) has been experimentally detected by Hunger et al. using magic-angle spinning nuclear magnetic resonance (MAS-NMR).⁶³ The authors observed that when tetrapropylammonium (TPA), acting as the SDA, was used to synthesize zeolite ZSM-5, internal silanols were formed as defects in a concentration of up to 8%. Sokol et al., employing DFT under periodic boundary conditions, found a relatively low energy requirement for the hydrolysis of the Si-O-Si linkages of sodalite with the subsequent formation of 2Si-2MR.⁶⁴ They suggested that the 2MR may undergo dehydrogenation, creating triplet peroxide species, which transforms into a singlet peroxy bridge with an $O–O$ distance of approximately 1.58 Å. The authors further proposed that the migration of these peroxy bridges could be the reason for the increase in the Brønsted and Lewis acidity, which makes the redox chemistry of these materials even more complex.⁶⁴ 

In this study, the dissociation of adsorbed water that led to the formation of 2Si-2MR increased the energy of the system by +108 kJ/mol (see Table IV). The $Si–O$ bond lengths of 2Si-2MR remained between 1.795 and 1.866 Å, with an $O1–O2$ distance of 2.147 Å, in agreement with previous calculations.⁶⁴ 

The substitution of one Si by Ti within the 2MR (labeled as Si,Ti-2MR) increased the stability of the structure, which was apparent from the reduction of the energy difference relative to non-dissociated water. In this case, the Si,Ti-2MR was +50 kJ/mol less stable than the system with non-dissociated water, representing a marked decrease in energy difference when compared with the value of +108 kJ/mol found for 2Si-2MR (see Table IV). Nevertheless, Si,Ti-2MR was still hardly stable, with a mere relative energy of −2 kJ/mol against the standard of water in the gas phase plus the bare zeolite without adsorbate. In this regard, it has been reported that Si,Ti-2MR is the least stable configuration when compared with models that reflect different levels of hydration and inversion of the tetrahedra $SiO4$ and $TiO4$ such as tripodal, tetrapodal, and bipodal configurations.⁶⁵ 

When we considered the 2MR formed by two Ti atoms (labeled as 2Ti-2MR), the stability of 2Ti-2MR increased significantly. As shown in Table IV, 2Ti-2MR was only 15 kJ/mol less stable than the system with non-dissociated water, in contrast to the values of +108 and +50 kJ/mol for 2Si-2MR and Si,Ti-2MR, respectively. In addition, the adsorption energy of dissociated water forming the 2Ti-2MR structure was −56 kJ/mol, while the adsorption of non-dissociated water, which served as the precursor for the formation of 2Ti-2MR, released 71 kJ/mol. These adsorption energies were calculated as relative energies against water in the gas phase and the corresponding non-hydrated Ti configuration in 2Ti-2MR (labeled as 0_b in Fig. 2). Configuration 0_b was 19 kJ/mol less stable than 0_a; however, upon water adsorption in the $β$-cage, 0_b was 4 kJ/mol more stable than 0_a. This difference increases if the 2Ti-2MR formation is considered: 2Ti-2MR(0_b) was 32 kJ/mol more stable than 2Ti-2MR(0_a), which is due to the higher structural flexibility provided by locating the two Ti atoms in different D4Rs (0_b case) than in the same D4R (0_a case). Therefore, we reason that the presence of water during the synthesis process can stabilize 0_b, increasing its probability of occurrence against configurations that seem more stable in non-hydrated conditions, such as 0_a and 2_a. As such, the presence of configuration 0_b in zeolite LTA could have important consequences for the performance of this material as a catalyst for redox reactions: two Ti atoms, at close distance from each other and with the capacity to interact simultaneously with the same molecule, may adsorb and stabilize more effectively peroxide species than single Ti-substituted frameworks.¹⁷ 

We used the single Ti-substituted LTA and configuration 0_b, with water adsorbed in their respective $β$-cages [see Figs. 6(a) and 6(b)], as representative models to examine the electronic structure of hydrated Ti-LTA. The system water/0_b represented in Fig. 6(b) was the precursor for the formation of 2Ti-2MR(0_b), with an adsorption energy of −71 kJ/mol.

The electronic charge density difference (⁠$Δρ$⁠) is a simple but effective method to observe the movement of electronic charge within a system constituted by several individual components, such as water and Ti-LTA for the water/1Ti and water/0_b systems. Equation (12) is used to obtain (⁠$Δρ$⁠),

where $ρ(water∕ Ti−LTA)$ is the electronic charge density of the system water/Ti-LTA, and $ρ(water)$ and $ρ(Ti−LTA)$ are the densities of water and Ti-LTA, respectively, calculated separately but retaining the same box size and calculation conditions as for water/Ti-LTA. Figures 6(c) and 6(d) show $Δρ(water∕ 1Ti)$ and $Δρ(water∕ 0_b)$⁠, respectively, where we clearly observe the movement of charge density from the $O–H$ bonds of water towards the $H2O–Ti$ region. The determination of Bader atomic charges confirmed this charge transfer.^66–68 The water molecule adsorbed in 1Ti yielded a charge of +0.030 e⁻, whilst the value for water in 0_b was +0.042 e⁻. The increase in the positive electronic charge of water from 1Ti to 0_b is expected; water interacts with both Ti atoms simultaneously in 0_b, and the LUB of 0_b efficiently surrounds the water molecule [see Fig. 6(f)]. The shape of the positive section of $Δρ(water∕ 0_b)$ shows that the electronic charge moves evenly towards both Ti atoms.

Figure 7 shows the projected density of states (PDOS) onto the O atom of water in the gas phase and water in the systems water/1Ti and water/0_b. Also shown are the PDOS onto the Ti atoms of these systems. We observed that the main indications of hybridization between the electronic states of water and the Lewis acid site are provided by the analysis of the water molecule’s orbitals with symmetry 1b₁, which is the highest occupied molecular orbital (HOMO), 3a₁ (HOMO-1) and 1b₂ (HOMO-2) [see Fig. 7(a)].

The molecular orbital 1b₁ corresponds to the non-bonding lone pair of water [see charge density isosurfaces in Fig. 7(a)], and it shows the largest evidence of hybridization with the states of the Lewis acid. The peak that represents 1b₁ in the gas phase scatters over energy regions that span from −2.7 to 0.0 eV for water/1Ti and from −3.0 to 0.0 eV for water/0_b (with the Fermi energy shifted to 0.0 eV). The 3a₁ orbital also scatters but covering a smaller range of energy. These regions match the profiles of PDOS(Ti), which is confirmation of the states’ hybridization. For instance, the visualization of the electronic charge density corresponding to the most intense peaks within the 1b₁ region of systems water/1Ti and water/0_b (indicated with red arrows in Fig. 7) clearly shows continuous isosurfaces that go from the water molecule to the Ti atoms.

At the same time, the feature of 1b₂ as a single and intense peak is retained even after the adsorption of water on the Lewis site, in a section of energy where there are no peaks associated with PDOS(Ti). This means that 1b₂ remains practically intact, without mixing with the electronic states of the Lewis acid site.

We have simulated the IR spectra of the 2MRs with the GGA methodology and compared them with the profiles belonging to the equivalent structures without adsorbates, as shown in Fig. 8. The main feature related to the formation of the 2MRs, which may serve to identify this structure in any experiment, is the splitting of one of the bands associated with the D4R. For instance, the signal at 598 cm⁻¹ in the pure-silica spectrum evolves into two weak bands at 604 and 569 cm⁻¹ after the 2Si-2MR formation (indicated by small vertical arrows in the spectra). The same tendency is followed by the single Ti-substituted structure, where two weak bands at 604 and 571 cm⁻¹ are observed for Si,Ti-2MR. The spectrum of 2Ti-2MR shows a similar splitting but it tends to be slightly larger than in the previous two cases. The separation between the two weak bands is 35 and 33 cm⁻¹ for 2Si-2MR and Si,Ti-2MR, respectively, increasing to 54 cm⁻¹ for 2Ti-2MR.

We have performed density functional theory and interatomic potential-based calculations to study the pure-silica zeolite LTA and its Ti-substituted forms. The calculated mechanical properties of the pure-silica LTA predict that this material is more rigid than zeolite-A, which has the same framework type but with a Si/Al ratio of 1. The IP method tends to overestimate the value of the elastic moduli compared with the GGA, a trend that agrees with previous calculations on zeolite MFI.

The energetics of the Ti-substituted LTA suggests that the addition of a second Ti atom within the framework, at close distance from the first Ti, is a favorable process. In fact, both the GGA and the IP method predict that the global minimum corresponds to a double Ti substitution in a 6MR, with the two Ti atoms symmetrically opposing each other with a configuration $Ti–(Si)2–Ti$⁠.

The interaction of a water molecule with Ti has a calculated adsorption energy between −37 and −52 kJ/mol, where the adsorption in the $β$-cage is stronger than that in the $α$-cage, by at least 4 kJ/mol. The direct dissociation of water is only favored when it leads to the formation of a 2MR between T-sites of different D4Rs. The stability is especially high for the 2MR $Ti–(OH)2–Ti$⁠, with an energy that is comparable to non-dissociated water adsorbed at the Ti sites. The electronic analysis shows a net charge transfer from water to the Lewis acid site and a strong hybridization between their electronic states.

The GGA-simulated infrared spectrum of pure-silica LTA shows good agreement with experimental measurements, and it is capable of reproducing the four main bands within an error of −30 to −20 cm⁻¹. The Ti substitution causes the occurrence of new weak signals between 1200 and 700 cm⁻¹, depending on the configuration of the Ti atoms. The main feature related to the formation of the 2MRs is the splitting of the weak band at approximately 600 cm⁻¹, which is associated with the D4R building unit, into two weak bands that are separated by 33 to 54 cm⁻¹.

See supplementary material for the tabulated structural parameters of zeolite LTA after the Ti incorporation. In addition, the energies of the LUB for pure-silica LTA, the single Ti substitution and the structures 0_a and 2_a, and the isosurfaces of the electronic charge density associated with the LUB are shown. Furthermore, the derivation of Eq. (11) from the expression to evaluate the absolute stability of the Ti-substitution is presented.

This work was performed using the computational facilities of the Centre for High Performance Computing in Cape Town (CHPC). Via our membership of the UK’s HEC Materials Chemistry Consortium, which is funded by EPSRC [Grant No. EP/L000202], this work used the ARCHER UK National Supercomputing Service (http://www.archer.ac.uk). This work has also used the computational facilities of the Advanced Research Computing at Cardiff (ARCCA) Division, Cardiff University, and HPC Wales. The authors acknowledge the UK Economic and Social Research Council (ESRC) [Grant No. ES/N013867/1] and South Africa National Research Foundation for funding under the UK-SA Newton Fund. P.E.N. thanks the South African Research Chair Initiative of the Department of Science and Technology and the National Research Foundation. C.E.H.-T. acknowledges UCL and EPSRC [Grant No. EP/K009567/2] for a scholarship, and N.H.d.L. thanks the Royal Society for an Industry Fellowship. C.E.H.-T. also thanks Katlego Phoshoko and Kenneth Kgatwane for the invaluable assistance and Professor C. R. A. Catlow for a very fruitful discussion. All data created during this research are openly available from Cardiff University Research Portal at http://doi.org/10.17035/d.2016.0011512285.