Potential Functions for Silica and Zeolite Catalysts Based on ab

Leanna K. Toy, Scott M. Auerbach, Thom Vreven, and Michael J. Frisch .... Alexandra Simperler, Robert G. Bell, Martin D. Foster, Aileen E. Gray, D...
19 downloads 0 Views 299KB Size
Download PDF
J. Phys. Chem. 1996, 100, 11043-11049

11043

Potential Functions for Silica and Zeolite Catalysts Based on ab Initio Calculations. 3. A Shell Model Ion Pair Potential for Silica and Aluminosilicates† Klaus-Peter Schro1 der and Joachim Sauer* Max Planck Society, Research Unit “Quantum Chemistry” at the Humboldt UniVersity, Ja¨ gerstrasse 10/11, D-10117 Berlin, Germany ReceiVed: NoVember 20, 1995X

A shell model potential for silica and protonated zeolite catalysts is presented. The potential parameters are fitted exclusively to ab initio data generated by nonempirical quantum chemical calculations on small molecular models made of SiO4 and AlO4 tetrahedra. The Hartree-Fock method has been used with a basis set of double zeta + polarization quality on Si, Al, and H atoms and of valence triple zeta + polarization quality on O atoms. Comparison is made with an ab initio molecular mechanics force field previously derived from the same data and with an empirical parametrization of the shell model potential. The power of the new potential for predicting the crystal structures and the properties is demonstrated for a variety of silica and zeolite structure types. Cell parameters of dense and microporous silica are predicted with a mean error of 1.8%. Compared to earlier potentials, substantial progress is made in predicting dynamical properties. Examples are the phonon dispersion curves of R-quartz and the infrared spectrum of zeolite ZSM-5. As a first application to acidic zeolite catalysts the local structures and vibrational frequencies of the different bridging hydroxyl groups in faujasite are calculated. The results do not show the artifacts exhibited by previous potentials.

1. Introduction Zeolites and related microporous materials often have such complex unit cells that their structures cannot be solved without assistance from modeling. Academic1,2 and commercial codes3,4 are available which predict periodic structures by lattice energy minimization and vibrational spectra by lattice dynamics. Finite temperature effects can be included by free energy minimization2,5 or molecular dynamics. Moreover, once the structures are obtained, diffraction patterns can be calculated and NMR spectra simulated.6 The success of all these predictions critically depends on the reliability of the potential functions used to describe the interactions between the atoms. Traditionally, empirical potential functions are used, but their applicability to zeolite catalysts is hampered by the problem of the availability of relevant and reliable observed data for deriving potential parameters. For example, little is known from experiments about the local structure of the active sites of acidic zeolite catalysts, the bridging hydroxyl groups, SiO(H)Al, and their location in the framework. This is connected with their low concentration (typically 1 Al among 30-50 Si atoms) and with the fact that they are not ordered into unit cells. This is a strong incentive for using quantum chemical ab initio techniques for structure predictions. However, a direct ab initio description of the full periodic structure of a zeolite catalyst is an enormous computational task because the unit cells are very large and show little point symmetry. Hence periodic ab initio calculations on zeolites are limited to energy calculations for given structures with not too many atoms in the unit cell.7-10 Only for some silica polymorphs with only a few atoms in their asymmetric units have the structures been determined by energy minimization.11,12 A way out of this dilemmasthere are not enough observed data to generate reliable empirical potential functions, but * Corresponding author. Telephone: +49-30-20192300. FAX: +4930-20192302. E-mail: [email protected]. † See refs 13 and 14 for parts 1 and 2 of this series. X Abstract published in AdVance ACS Abstracts, June 1, 1996.

S0022-3654(95)03405-8 CCC: $12.00

periodic ab initio calculations are not feasiblesis the following strategy:13,14 Quantum chemical ab initio calculations are performed for finite molecular models of the periodic zeolite structures. The results obtained are used to fit the parameters of the potential functions. These parameters are assumed to be transferable and used for calculating the energies of the periodic lattices and the forces on its atoms. This approach has several advantages. The set of data that may be included is in principle unlimited, and there is a one-to-one correspondence between the “predicted” (by the potential function) and the “observed” (by the ab initio calculation) data. The energy and first and second derivatives with respect to the displacement of all nuclei serve as “observables”. They can be generated not only for the energy minimum structure of the models but also for an, in principle, unlimited number of distorted structures. Two previous papers in this series described the data base generated for silica13 and H-zeolites.14 As the functional form of the potential functions a molecular mechanics type force field was selected, which is typically used for organic molecules with more or less covalent bonds. Specifically, the “consistent force field” (cff)15 was adopted which has the advantage that parameters are already available for large classes of organic molecules, e.g., hydrocarbons15,16 and polycarbonates.17 This allows the simulation of adsorbed molecules inside the zeolite pores without much additional parametrization work. The bonding in silica and zeolites has a significant ionic component, which makes an ion pair model the natural choice for a potential function. Both rigid ion pair potentials18-21 and shell model potentials21-25 which take the polarizability of the oxygen anions into account have been used in the past. This paper reports a fit of the shell model potential to a subset of the same data as previously used for fitting the consistent force field.14 The structures and spectra predicted with the ab initio shell model potential obtained are compared with the corresponding results of the ab initio consistent force field. This allows one to judge the influence that the analytical form of a potential has on the results. The predictions made with the ab © 1996 American Chemical Society

11044 J. Phys. Chem., Vol. 100, No. 26, 1996

Schro¨der and Sauer is an adjustable parameter the ion charge is not. It is most convenient to assign full formal charges to the ions. This choice allows one to treat defects or materials with different chemical compositions without problems. Moreover, it has been recently demonstrated that a shell model cannot properly describe the dynamics of quartz and sapphire unless a formal charge of -2 is chosen for the oxygen ions.30 The full potential expression is

E ) Eelec + Ecore-shell + Eshort-range + E3-body

(1)

The first three terms are pair interactions. The Coulomb energy, Eelec, is a sum that runs over all point charges (cores and shells),

Eelec ) ∑qiqjrij-1

(2)

i,j

where rij is the distance between the point charges i and j. The core-shell interaction runs over all ions for which the shell model is used and acts between core and shell:

Ecore-shell ) ∑kiri2

(3)

i

Figure 1. Molecular models used in the ab initio calculations.

initio parametrized shell model potential are also compared with those made with a basically empirically parametrized shell model potential.21-23 This allows one to study the influence of the data base and the parameters on the results.

ri is the actual distance between the ith core and its shell. We use the shell model only for anions; this means here for oxygen ions only. We use a purely repulsive term between shells and cations to describe the short-range interaction of ions:

Eshort-range ) ∑Akl exp(-rklFkl-1)

(4)

k,l

This function is frequently extended by a term which describes the dispersion energy,

-∑Cklrkl-6

2. Calculations

(5)

k,l

2.1. The Quantum Chemical Data Base. Previously, a large quantum chemical data base has been established for molecules that are built of SiO4 and AlO4 tetrahedra.13,14 The SCF method has been used with a basis set of double zeta plus polarization quality on Si, Al, and H and of valence triple zeta plus polarization quality on O.26 For the Si,Al/O/H atoms the (11s, 7p)/(9s, 5p)/(4s) Gaussian basis sets were taken from Huzinaga27,28 and contracted into [6s, 4p]/[5s, 3p]/[2s] according to the following pattern {5 2 1 1 1 1, 4 1 1 1}/{5 1 1 1 1, 3 1 1}/{3 1}. Polarization functions with the following exponents were added: 0.4 (Si), 0.3 (Al), 1.2 (O), 0.8 (H). For the purposes of the present work only a part of these data has been adopted since the shell model potential includes substantially fewer adjustable parameters than the consistent force field. We selected exclusively ring structures (Figure 1), since microporous aluminosilicates are built mainly of four-, five-, and sixmembered rings. To include data that correspond to aluminumrich zeolites a four-membered ring with two Si and two Al atoms was considered. For all models the optimized structures with their zero gradients of the energy with respect to the positions of the nuclei were used in the fit, and for the Al-free four- and five-membered rings and the 4-membered ring with 2 Al atoms the force constants were included. 2.2. The Shell Model Potential. The shell model potential, as introduced by Dick and Overhauser,29 is an ion pair potential which represents the ions (often the anions only) by a pair of point charges, the positive core and the negative and massless shell, which are connected by a harmonic spring. All shortrange interactions are defined between shells rather than cores which is consistent with the assumption that the short-range repulsion is caused by electron repulsion. While the shell charge

We neglect this term since the Hartree-Fock method was used to generate the data base which does not provide any dispersion contributions. Between ions of equal charge no short-range repulsion term is defined since the electrostatic repulsion dominates. For tetrahedrally coordinated atoms often a so-called threebody interaction term is added to the potential expression:

E3-body )

1

kib(θijk - θ0)2 ∑i ∑ j,k

2

(6)

The indices j and k run over the shells of all oxygen atoms that are bonded to the T atom i. This term is typical for molecular mechanics potential functions rather than for ionic ones. It gives the TO4 tetrahedron additional stiffness. The reference angle θ0 usually is 109.47°. It is assumed that such a term is needed to take into account that the T-O bond is partly covalent and that the hybridization of the T atom is sp3.20,22 2.3. Parameter Fitting. The shell model is designed for ionic crystals. However, the data base was generated for molecule-type models of partially ionic materials. Hence, it was necessary to use the shell model to describe molecules. In particular the large formal charges could cause problems. The question was how to treat the bridging hydroxyl groups and the “boundaries” of the moleculessthe terminal hydroxyl groups. It turned out to be impossible to obtain a good fit with only one type of oxygen ion. It was necessary to distinguish between a 2-fold coordinated oxygen ion in terminal OH groups or between two T ions, O, and a 3-fold coordinated one in bridging OH groups, Ob (Figure 2). The shell model was

Potential Functions for Silica and Zeolite Catalysts

J. Phys. Chem., Vol. 100, No. 26, 1996 11045 TABLE 1: Parameters of the ab Initio Shell Model Potential charge (e) Si Al O Ob Hb Ht

Figure 2. Definition of ion types used in the ab initio shell model potential demonstrated for the four-membered ring model containing one Al ion.

applied to both of them. All other ions were treated as rigid (cores only). Different hydrogen ion types were assumed for bridging and terminal hydroxyl groups, Hb and Ht, and different short-range parameters were fitted for the Ob-Hb and O-Ht bonds. Since the nonbonded interaction between O and Hb is described by parameters different from those of the “bonding” Ob-Hb interaction the potential is not suited to describe the proton transition between different framework oxygen sites. The previous empirical shell model potential for H-zeolites21 deviates from the original shell model approach in three respects: (i) The oxygen and the hydrogen ion of the bridging hydroxyl obtained partial charges (O, -1.426 e, H, 0.426 e). (ii) The bridging oxygen ion was assumed to be nonpolarizable, i.e., it was treated as a core only. (iii) The interaction between the O and the H ion of the bridging hydroxyl groups was described by a Morse potential. The empirical shell model potential is, therefore, less consistent than the present ab initio potential in both the functional form as well as the source of the data. Weighting different types of data by different factors during a fit is a standard procedure to obtain a well-balanced parametrization. Most of the oxygen ions in the molecular models considered are located in terminal hydroxyls. These groups are needed to terminate the models used in quantum chemical calculations, but they are not present in the periodical crystal to which the derived potential functions are eventually applied. To limit the influence of the terminal hydroxyls on the results of the fit, all data connected with them were weighted with a factor of 0.1. Since we were interested in good predictions of both structures and vibrational spectra we decided to weight the ab initio force constants with a factor of 0.01 since their number in the data base was much larger than the number of gradient components. Energies were not used in the fit. Other authors24 report on multiple step procedures in least squares fits. This can be useful if there is a group of parameters that is independent of all other parameters or that depend on a subset of data only. In our case, however, it seemed unnecessary to adopt such an approach, since the number of adjustable parameters was low (20 in the final run) and since we could not exclude interdependencies. All parameters have been adjusted in one step. The final values are given in Table 1. 2.4. Lattice Energy Minimization. To obtain static and dynamical properties of silica and zeolites we used standard lattice modeling techniques.31 On the basis of the shell model potential, the energy of a unit cell in an infinite periodical lattice has been calculated by accurate summation of the individual pair and three-body interactions. The Ewald summation technique has been used for the long-range Coulomb interaction. A cut-off radius of 10 Å has been chosen for the summation of the short-range interactions. The fractional coordinates of all cores and shells as well as the cell parameters have been varied to minimize the lattice energy. Properties like elastic and dielectric constants as well as phonon dispersion curves and

core

shell

4.0a 3.0a 1.06237 1.23944 1.0a 1.0a

-3.06237 -3.23944

short-range repulsion Si-O Si-Ob Al-O Al-Ob Ob-Hb O-Hb O-Ht

A (eV)

F (Å)

1550.950 2078.349 1068.711 887.969 452.466 6758.796 1034.329

0.30017 0.28130 0.32260 0.32377 0.21143 0.20511 0.16886

core-shell interaction k (eV Å-2) O Ob

112.7629 137.0634 three-body interaction

Si Al a

kb (eV rad-2)

θ0 (deg)

0.18397 0.64984

109.47a 109.47a

Not adjusted in the fitting procedure.

infrared spectra have been determined for this equilibrium structure. The calculations have been performed with the METAPOCS1 and GULP2 codes. 3. Results and Discussion 3.1. Structure of Silica Polymorphs. Silica shows a rich polymorphism. A reliable zeolite potential is expected to yield a good description of both dense and microporous silica. Table 2 shows the unit cell parameters of selected silica modifications. The deviation of the calculated from the experimental values is generally small, in most cases less than 2%. The average deviation of cell parameters for all silica modifications in Table 2 is 1.8%. Our shell model results agree somewhat better with experimental findings than predictions obtained with the molecular mechanics potential fitted to the same ab initio data (3.1% average deviation). For the dense modifications, quartz and cristobalite, the low-temperature forms are obtained as expected. By contrast, the cff potential yields the hightemperature β forms which possess higher symmetry. Both ab initio derived potentials predict larger cell parameters than the empirical shell model potential, which is in even better agreement with observations (0.7% average deviation of cell parameters). The source of the data base, namely the HartreeFock method, seems to be the reason for this systematic deviation. This is confirmed by calculations on the structure of pure-silica chabazite. Results12 of a periodical Hartree-Fock calculation using a split-valence basis set with polarization functions are compared to those obtained by a lattice energy minimization with the ab initio shell model potential (Table 3). Although slightly different basis sets have been used, cell parameters and fractional coordinates are almost identical. The empirical shell model potential yields a somewhat smaller cell parameter of 919.9 pm but very similar fractional coordinates. Experimental values are not available since naturally occurring chabazite has an Si:Al ratio of approximately 2.

11046 J. Phys. Chem., Vol. 100, No. 26, 1996

Schro¨der and Sauer

TABLE 2: Cell Parameters of Microporous and Dense SiO2 Modifications (in Picometers and Degrees) empirical shell ab initio ab initio model21 cff14 shell model obsda

modification A)B)C A)B)C A B C mordenite A B C ZSM-5 A B C R zeolite rho A)B)C R-quartz A)B C β-quartz A)B C R-cristobalite A ) B C β-cristobalite A ) B ) C coesite A B C β faujasite sodalite theta-1

2423 882 1382 1739 500 1802 2004 743 1998 1974 1332 90.8 1477 484 535 500 550 497 701 731 703 1229 712 122.5

2484 906 1438 1806 529 1834 2071 769 2066 2060 1384 90 1520 b 517 566 c 745e 729e 1273e 746e 120.6e

2463 895 1414 1777 516 1830 2049 758 2043 2021 1363 90 1501 499 551 510 562 513 727 739 719 1254 726 122.8

2426 883 1386 1742 504 1810 2038 749 2011 1988 1337 90.7 1485 492 541 500d 546d 498 695 717 714 1237 717 120.3

a Faujasite, ref 32; sodalite, ref 33; theta-1, ref 34; mordenite, ref 35; ZSM-5, ref 36; zeolite rho, ref 37; R-quartz, ref 38; R- and β-cristobalite, ref 39; coesite, ref 40. b Relaxes to β-quartz. c Relaxes to β-cristobalite. d Calculated from data in ref 41. e Reference 42.

TABLE 3: Cell Parameters (in Picometers and Degrees), Fractional Coordinates, Bond Lengths (in Picometers) and Bending Angles (in Degrees) of Pure-Silica Chabazite Optimized by a Crystal Orbital Method and by Lattice Energy Minimization Using the ab Initio Shell Model Potential Apra´ et al.12

this work

Apra´ et al.12

this work

a

932.6

933.5

γ

94.7

94.5

xSi ySi zSi xO(1) xO(2)

0.106 0.333 0.878 0.255 0.146

0.106 0.333 0.877 0.257 0.147

xO(3) zO(3) xO(4) zO(4)

0.250 0.897 0.023 0.326

0.249 0.894 0.022 0.325

rSi-O(1) rSi-O(2) rSi-O(3) rSi-O(4) 〈rSi-O〉

161.0 160.5 160.7 161.4 160.9

161.5 161.2 161.0 161.7 161.3

∠SiO(1)Si ∠SiO(2)Si ∠SiO(3)Si ∠SiO(4)Si 〈∠SiOSi〉

153 148 152 151 151

151 147 152 151 150

Table 4 compares observed and calculated bond lengths and angles for R-quartz and silica faujasite, the end member of the zeolite X and Y family. The averaged calculated Si-O bond length is about 1 pm longer than the observed one. For faujasite the three potentials agree in predicting rSi-O(4) being an Si-O bond of average length. This finding is at variance with the result of structure refinement from a neutron diffraction experiment32 that Si-O(4) was clearly the longest bond. The overestimation of Si-O-Si angles (by about 4°) is mainly responsible for the too large calculated unit cell volume. The biggest deviation from experiment has been obtained for the Si-O-Si angles within the six-membered rings of the sodalite unit of faujasite (at O(2) and O(3)). This is valid for both ab initio potentials, although the deviation is much larger for the molecular mechanics force field than for the shell model potential. The empirical parametrization of the shell model potential predicts Si-O bond lengths very close to the experimental ones. For faujasite the Si-O-Si bending angles are exactly the same as reported from the neutron diffraction

TABLE 4: Bond Lengths and Si-O-Si Bending Angles in r-Quartz and Pure-Silica Faujasite (in Picometers and Degrees) empirical shell model21 R-quartz

faujasite

a

rSi-O ∠SiOSi rSi-O(1) rSi-O(2) rSi-O(3) rSi-O(4) 〈rSi-O〉 ∠SiO(1)Si ∠SiO(2)Si ∠SiO(3)Si ∠SiO(4)Si

160.8 161.5 139 161.4 159.9 160.9 160.8 160.8 138 149 146 141

ab initio cff14

ab initio shell model

obsda

163.3 162.3 162.8 162.6 162.8 135 165 162 145

161.4 161.7 148 162.2 160.9 161.7 161.6 161.6 141 154 152 145

160.5 161.4 144 160.7 159.7 160.4 161.4 160.6 138 149 146 141

For refs, see Table 2.

TABLE 5: Elastic (in 1010 N/m2) and Dielectric Constants of r-Quartz empirical shell model21

ab initio shell model

obsd38

C11 C33 C44 C66 C14 C13 C12

9.47 11.61 5.01 3.82 -1.45 1.97 1.84

8.45 9.63 4.11 3.65 -1.37 2.39 1.15

8.69 10.60 5.83 3.99 -1.81 1.19 0.70

11 33 ∞

4.74 5.01 2.12

4.07 4.42 1.76

4.52 4.64 2.40

experiment.32 For R-quartz the Si-O-Si bending angle is underestimated by 5°. This indicates that the empirical shell model potential provides better structure predictions for microporous than for dense silica modifications. This is somewhat surprising, since its parameters were adjusted to R-quartz data. In contrast, the quality of the ab initio shell model potential appears to be the same for the whole class of silica polymorphs. The predictions of all considered potentials are reasonable. A further test of the potential is the comparison of calculated and measured elastic constants of R-quartz (Table 5). The calculated constants are of the same quality as those obtained with the empirical potential, although the parameters of the latter have been fitted to the experimental values. As a consequence of the large core-shell spring constant the polarizability of the oxygen anions is underestimated and the dielectric constants are too small. A full neglect of polarizability in rigid ion or molecular mechanics potentials leads, however, to a highfrequency dielectric constant of exactly 1. The advantage of the shell model approach is evident. 3.2. Infrared Spectra of Silica Polymorphs. One main goal of the new potential is to model the dynamical properties of zeolites to support spectroscopists in assigning infrared and Raman bands to vibrational modes. None of the previous potentials is able to satisfactorily predict vibrational spectra of silica or zeolites. Figure 3 shows calculated phonon dispersion curves for R-quartz in the [0 0 ξ] direction. In addition to the shell model potentials considered so far, we consider for comparison Kramer’s rigid ion potential20 and a shell model potential for R-quartz which is based on density functional calculations.25 The present ab initio shell model potential provides a remarkable improvement over the previously published rigid ion and shell model potentials. The two previous parametrizations of the shell model potential functions are not able to adequately describe the Si-O bond stretching vibrations in the range between 1000 and 1300 cm-1. The rigid ion

Potential Functions for Silica and Zeolite Catalysts

J. Phys. Chem., Vol. 100, No. 26, 1996 11047

Figure 5. The four different types of bridging hydroxyl groups in silicon-rich faujasite. For the purpose of clarity, oxygen atoms are left out unless they are part of the hydroxyl group. Figure 3. Experimental43 ([) and calculated (solid lines) phonon dispersion curves of R-quartz along the [0 0 ξ] direction. Calculated values were obtained with the present ab initio shell model potential (a), Kramer’s rigid ion potential20 (b), the empirical shell model potential21 (c), and Purton’s DFT-derived SiO2 shell model potential25 (d). An extended zone scheme is used to plot the dispersion curves for the three different irreducible representations separately.

TABLE 6: Relative Energies and OH Bond Stretching Frequencies of the Bridging OH Groups in Silicon-Rich H-Faujasite (Results Obtained with the Empirical Shell Model Potential23 in Parentheses) νOH (cm-1) O(1)H O(2)H O(3)H O(4)H a

Figure 4. Experimental44 infrared spectrum of ZSM-5 (a) and infrared spectra of silicalite calculated with the ab initio shell model potential (b and c), with the empirical shell model potential (d), and with the ab initio consistent force field (e). The curves c, d, and e have been generated from the discrete spectra assuming Gaussian bandshape.

potential, which describes the spectrum in this region better, has flaws in the lower frequency region. Only the ab initio shell model potential provides phonon dispersion curves, which are in very good agreement with experiment over the whole range. Small deviations from experiments are found in the middle region between 700 and 900 cm-1 only. It is important to know whether this promising predictive power for vibrational spectra also extends to microporous polymorphs with a broader range of Si-O-Si angles and of Si-O bond lengths. Such a case is the complex lattice of silicalite, the aluminum-free form of zeolite ZSM-5, which has 96 Si and 192 O atoms in the unit cell. Figure 4 shows infrared spectra calculated at the zone center (zero wave vector). The ab initio shell model potential reproduces the main features of the measured spectrum well. The maxima of three of the five main bands are found at the same positions (at about 550, 1100, and 1225 cm-1), the other two are only slightly shifted (from observed 441 to 402 cm-1 and from 798 to 861 cm-1). The largest deviation is found in the region of the symmetric Si-O stretches between 700 and 900 cm-1. This seems to be a systematic feature of the parameter set since it was also found for R-quartz. The empirical shell model potential also yields

∆Erel (kJ mol-1) calcd

calcda

obsd46

8.8 (5.3) 18.0 (19.8) 0.0 (0.0) 27.9 (23.7)

3626 (3622) 3506 (3552) 3569 (3586) 3572 (3601)

3623 3550

See text.

all five major bands of the observed spectrum. Below 600 cm-1, in the region of bending modes, both shell model potentials yield almost identical shapes of the spectra. However, the asymmetric Si-O stretching bands at 1100 and 1225 cm-1 are predicted at much too low frequencies (912 and 1081 cm-1). The ab initio molecular mechanics force field is obviously not suited for predicting infrared spectra, although its parameters have also been fitted to force constants. It yields only four bands in the range from 700 to 1100 cm-1. 3.3. Structure of H-Zeolites. For applications as catalysts the active H-forms of silicon-rich zeolites are of major interest. As already mentioned, many details of their structures are not known from experiments. Among the different H-zeolites, silicon-rich H-faujasite has become the favorite subject of both experimental and theoretical studies. The reason is that macroscopic observations of H-faujasite are well understood at the atomistic level. Faujasite has a large unit cell of high symmetry (space group Fd3hm). All tetrahedral atoms (Si and Al) are crystallographically equivalent as long as no difference is made between T ) Si and T ) Al, but because none of the symmetry elements passes through the tetrahedral sites there is no local symmetry and, as a consequence, all four oxygen sites in the TO4 tetrahedron are distinct. For the high-silica form of H-faujasite, where the Al atoms are isolated from each other, only four types of bridging hydroxyl groups can emerge (Figure 5). Only two of them, the O(1)H and the O(3)H groups, are actually detected in the infrared and 1H NMR spectra. The Al-H distances in these two bridging hydroxyl groups have been obtained by a numerical analysis of 1H MAS NMR sideband patterns.45 From a reliable potential one should expect that this still incomplete experimental information available is reproduced and that additional details of the structure are predicted which cannot be derived from experiments. Indeed, calculations with both shell model parametrizations for single bridging hydroxyl groups in the primitive cell of faujasite (one Al per 48 T sites) give the result that the O(1)H and O(3)H bridging OH groups are energetically preferred compared with the remaining two. Furthermore, it is fully confirmed by both potentials that the O(3)H hydroxyl group with the stronger bonded hydrogen atom corresponds to the low-frequency infrared band, whereas the

11048 J. Phys. Chem., Vol. 100, No. 26, 1996

Schro¨der and Sauer

TABLE 7: Geometries of the Four Bridging Hydroxyl Groups in Silicon-Rich Faujasite: Distances in Picometers, and Angles in Degrees rAl-O(H)

rAl-H

rSi-O(H)

rO-H

∠SiO(H)Al

∠SiOH

∠AlOH

〈rAl-O〉

〈rSi-O〉

O(1)H

ab initio shell model ab initio cff14 empirical shell model23 ab initio (HF)47

191.4 185.8 191.0 193.5

247.6 233.3 238.6 243.2

169.6 167.1 169.4 172.0

95.4 95.2 100.0 95.6

130.2 133.7 131.1 129.3

114.5 117.9 123.0 119.1

115.2 107.9 105.9 109.9

175.8 174.4 173.8 176.7

162.3 162.4 160.3

O(2)H

ab initio shell model ab initio cff14 empirical shell model23

189.6 185.7 190.3

234.7 222.6 229.6

168.9 168.1 168.7

96.9 95.9 100.4

141.8 141.3 142.4

111.9 113.6 117.6

105.4 99.5 99.7

175.8 175.7 174.3

162.6 163.3 160.7

O(3)H

ab initio shell model ab initio cff14 empirical shell model23

194.7 189.7 193.0

239.7 231.2 233.2

169.1 167.3 169.7

96.1 95.5 100.2

140.2 138.7 138.7

113.8 117.9 120.7

106.0 103.3 100.6

176.4 175.2 174.2

162.6 163.2 160.7

O(4)H

ab initio shell model ab initio cff14 empirical shell model23

189.4 184.0 190.6

241.6 225.1 236.3

169.6 167.7 168.8

96.1 95.6 100.1

132.9 136.5 134.0

112.7 115.1 121.4

111.4 102.6 104.4

175.9 175.6 174.3

162.4 162.8 160.3

O(1)H group is the source of the high-frequency band (Table 6). This is contrary to the sometimes assumed correlation between OH stretching frequency and acidity. To achieve more than qualitative agreement with experiment, the calculated harmonic OH bond stretching frequencies have to be multiplied by a scale factor, which refers to the anharmonicity of the vibration and to the known overestimation of OH bond stretching force constants in Hartree-Fock calculations. A factor of 0.86 has been chosen for scaling the frequencies obtained with the ab initio shell model potential. The frequencies obtained with the empirical shell model potential have been corrected for anharmonicity only.23 The final relative energies and vibrational frequencies are similar for both potentials. Table 7 shows the bond distances and angles predicted by different potentials for the bridging hydroxyl groups in siliconrich faujasite. The ab initio shell model potential yields structures that are consistent with results of a large-scale Hartree-Fock calculation on a molecular model which is large enough to reflect the typical connectivity of the faujasite framework. Its central Al atom is part of three four-membered rings and a six-membered ring.47 The Al-O(H) bond lengths predicted by the potential vary by about 5 pm, which is close to the range of 188-199 pm found for this distance in the molecular models. From 1H MAS NMR sideband patterns Al-H distances of 248 ( 4 and 240 ( 4 pm for O(1)H and O(3)H groups, respectively, have been deduced.45 This compares well with the predictions of the ab initio shell model potential (247.6 and 239.7 pm, respectively). Previous potentials for modeling H-forms of zeolites exhibit typical drawbacks in predicting the detailed structure of bridging hydroxyl groups (Table 7). The ab initio molecular mechanics potential shows not enough flexibility to model the variation of the Al-O(H) bond. As a consequence the Al-O(H) bond lengths in zeolitic bridging hydroxyl groups are systematically underestimated. An ionic shell model appears better suited to describe such critical cases. This is demonstrated for both parametrizations. However, the empirical shell model potential clearly overestimates the O-H distance. Lattice energy calculations were also performed for Hfaujasite with a Si:Al ratio of 2.43, typical for zeolite Y. A sample of the same module but with a small amount of sodium cations was recently studied by Czjzek et al.48 in a neutron powder diffraction experiment. The composition was Na3H53Al56Si136O384. We distributed the Al ions over the lattice sites according to a model suggested by Klinowski et al.,49 and the protons were distributed among the different oxygen sites as O(1):O(2):O(3):O(4) ) 8:2:4:0, which is close to the population deduced from the neutron diffraction data. The relaxed structure with space group P1 cannot be compared directly with the experimental one. Since the Al ions are not really ordered into

TABLE 8: Calculated and Observed Mean Bond Lengths and Angles as Well as the Unit Cell Parameter of Zeolite H-Y with a Si:Al Ratio of 2.43 (in Picometers and Degrees). In Parentheses: Deviations from the Values for Pure-Silica Faujasite empirical shell model23 〈rT-O(1)〉 〈rT-O(2)〉 〈rT-O(3)〉 〈rT-O(4)〉 〈∠TO(1)T〉 〈∠TO(2)T〉 〈∠TO(3)T〉 〈∠TO(4)T〉 a

168.9 (7.5) 163.0 (3.1) 165.6 (4.7) 162.3 (1.5) 136 (-2) 152 (3) 146 (0) 143 (2) 2484 (61)

ab initio cff14

ab initio shell model

obsd48

168.3 (5.0) 169.4 (7.2) 167.7 (7.0) 165.7 (3.4) 165.0 (4.1) 163.2 (3.5) 167.5 (4.7) 166.7 (5.0) 165.4 (5.0) 165.0 (2.4) 164.5 (2.9) 163.6 (2.2) 139 (4) 138 (-3) 136 (-2) 153 (-12) 156 (2) 145 (-4) 147 (-15) 150 (-2) 140 (-6) 147 (2) 144 (-1) 144 (3) 2531 (47) 2522 (59) 2477 (51)

the unit cell of zeolite Y, the refinement of the diffraction data resulted in averaged fractional coordinates only, without distinguishing between Si and Al sites and with identical fractional coordinates for protonated and unprotonated oxygen sites. The space group Fd3hm was assumed. Mean T-O distances are shown in Table 8. For comparison, we calculated for each of the four oxygen types the average of its T-O bond lengths, regardless of whether the particular O ion was protonated or not. This is not the same type of averaging as was involved in the refinement of the neutron diffraction data. Nevertheless, a comparison with the predictions made by the different potentials is useful. The results of the different calculations are similar to each other. For the mean T-O bond lengths the sequence O(1) > O(3) > O(2) > O(4) is predicted, which reflects the different degrees of protonation of the oxygen sites. This is reasonable, since the protonation of an oxygen lengthens its T-O bonds and decreases its T-O-T angle (see also Table 7). Thus, the oxygen site with the highest occupancy should show the largest mean T-O bond length. As a deviation from this rule the observed mean T-O bond is larger for O(4) than for O(2). This means that the occupation of 1/12 of all O(2) sites is not enough to make the mean T-O(2) bond longer than the mean T-O(4) bond, which is longer for pure-silica faujasite. The observed increase of the mean T-O bond lengths by between 2.2 pm for O(4) and 7.0 pm for O(1) is adequately described by the two shell model potentials. There is a linear dependence of the lengthening on the number of protonated sites. This is not valid for the molecular mechanics potential which underestimates the lengthening of the averaged T-O(1) bond. The situation is much less clear for T-O-T bond angles. By the neutron diffraction experiment it was found that the angles for all protonated sites are smaller than in pure-silica faujasite, with the largest changes for the less occupied sites O(2) and O(3) which are located in six-membered rings. None

Potential Functions for Silica and Zeolite Catalysts of the tested potentials confirms this trend. The shell model potentials predict only minor changes of bond angles compared with the pure-silica faujasite. In contrast, the molecular mechanics potential predicts large differences between the mean bond angles at O(2) and O(3) in the H-faujasite and the puresilica faujasite. 4. Conclusions An ab initio shell model potential for silica and H-zeolites has been derived by fitting its parameters to data obtained by quantum chemical calculations on molecular models which represent typical subunits of zeolites. The potential successfully models the crystal structures of silica and H-zeolites and appears to be a valuable tool for simulating infrared spectra. Comparison with two previous zeolite potentials provides information about the influence of the form of the potential expression as well as of the nature of the data base used for the fit. Comparison of the two ab initio potentials fitted to the same quantum chemical data base leads to the following conclusions: (i) Both ab initio potentials reproduce known crystal structures well. The shell model potential yields improved unit cell parameters and individual bond lengths and angles. In particular, it is better suited to describe strongly varying parameters like the T-O-T bond angle or the Al-O(H) bond length. We see the reason in the Taylor expansion type expression of the molecular mechanics potential which leads to large terms when an internal coordinate is far from its reference values. (ii) The shell model potential predicts infrared spectra of pure-silica zeolites in very close agreement with experimental findings. In contrast, the spectra predicted by the consistent force field are disappointing in view of the basically correct structure predictions and the large number of adjusted parameters (336). The reason is not yet clear. The influence of the data base used in the fit can be estimated comparing the results of both shell model potentials: (i) For silica polymorphs the empirical potential yields structures that are very close to observed ones, whereas the predictions of the ab initio potential show small systematic deviations from observations which are connected with known systematic errors of the Hartree-Fock method used in generating the data. The detailed structure of bridging hydroxyl groups in H-zeolites is more realistically described by the ab initio potential. Both parametrizations agree in describing the details of energetic and dynamical differences between bridging hydroxyl groups at different framework sites in H-faujasite. Generally, the ab initio potential shows more consistent results due to the uniform data base that was used for parametrization. (ii) The ab initio potential has significant advantages over the empirical potential in predicting lattice vibrations. Simulated and experimental spectra agree over large frequency ranges. Between 700 and 900 cm-1 the frequencies are shifted in a systematical way. The main improvement compared to the empirical potential is observed in the asymmetric Si-O stretching region above 1000 cm-1. We believe that the good prediction of vibrational spectra is due to inclusion of force constants into the quantum chemical data base used to fit the parameters of our potential. We conclude that the success of the present potential is based on the independent way of data generation by ab initio calculations and on the adequate analytical form of the shell model potential, which guarantees a proper description of the long-range Coulomb interaction and of the polarizability of the oxygen anions. Work is in progress which extends this approach to elements used for isomorphous substitutions of tetrahedral ions, e.g., Ti. Acknowledgment. This work has been supported by the Deutsche Forschungsgemeinschaft and the Fonds der Chemischen Industrie.

J. Phys. Chem., Vol. 100, No. 26, 1996 11049 References and Notes (1) Program METAPOCS: Catlow, C. R. A.; Cormack, A. N.; Theobald, F. Acta Crystallogr. B 1984, 40, 195. (2) GULP (the General Utility Lattice Program) written and developed by J. D. Gale, Royal Institution/Imperial College, UK, 1992-1994. (3) Program DISCOVER, Biosym Technologies, Inc., San Diego, CA, 1992. (4) Program CERIUS2, Molecular Simulations, Inc., Cambridge, 1993. (5) Program PARAPOCS: Parker, S. C.; Price, G. D. AdV. Solid State Chem. 1989, 1, 295. (6) Engelhardt, G.; Radeglia, R. Chem. Phys. Lett. 1984, 271, 108. (7) White, J. C.; Hess, A. C. J. Phys. Chem. 1993, 97, 6398. (8) White, J. C.; Hess, A. C. J. Phys. Chem. 1993, 97, 8703. (9) Anchell, J. L.; White, J. C.; Thompson, M. R.; Hess, A. C. J. Phys. Chem. 1994, 98, 4463. (10) Nicholas, J. B.; Hess, A. C. J. Am. Chem. Soc. 1994, 116, 5428. (11) Dovesi, R.; Roetti, C.; Freyria-Fava, C.; Apra´, E.; Saunders, V. R.; Harrison, N. M. Philos. Trans. R. Soc. London A 1992, 341, 203. (12) Apra´, E.; Dovesi, R.; Freyria-Fava, C.; Pisani, C.; Roetti, C.; Saunders, V. R. Modell. Simul. Mater. Sci. Eng. 1993, 1, 297. (13) Hill, J.-R.; Sauer, J. J. Phys. Chem. 1994, 98, 1238. (14) Hill, J.-R.; Sauer, J. J. Phys. Chem. 1995, 99, 9536. (15) Maple, J. R.; Hwang, M.-J.; Stockfish, T. P.; Dinur, U.; Waldman, M.; Ewig, C. S.; Hagler, A. T. J. Comput. Chem. 1994, 15, 162. (16) Hwang, M.-J.; Stockfish, T. P.; Hagler, A. T. J. Am. Chem. Soc. 1994, 116, 2515. (17) Sun, H.; Mumby, S.; Maple, J. R.; Hagler, A. T. J. Am. Chem. Soc. 1994, 116, 2978. (18) Tsuneyuki, S.; Tsukada, M.; Aoki, H.; Matsui, Y. Phys. ReV. Lett. 1988, 61, 869. (19) van Beest, B. W. H.; Kramer, G. J.; van Santen, R. A. Phys. ReV. Lett. 1990, 64, 1955. (20) Kramer, G. J.; de Man, A. J. M.; van Santen, R. A. J. Am. Chem. Soc. 1991, 113, 6435. (21) Jackson, R. A.; Catlow, C. R. A. Mol. Simul. 1988, 1, 207. (22) Sanders, M. J.; Leslie, M.; Catlow, C. R. A. J. Chem. Soc., Chem. Commun. 1984, 1271. (23) Schro¨der, K.-P.; Sauer, J.; Leslie, M.; Catlow, C. R. A.; Thomas, J. M. Chem. Phys. Lett. 1992, 188, 320. (24) de Boer, K.; Jansen, A. P. J.; van Santen, R. A. Chem. Phys. Lett. 1994, 223, 46. (25) Purton, J.; Jones, R.; Catlow, C. R. A.; Leslie, M. Phys. Chem. Miner. 1993, 19, 392. (26) Sauer, J.; Horn, H.; Ha¨ser, M.; Ahlrichs, R. Chem. Phys. Lett. 1990, 173, 26. (27) Huzinaga, S. Approximate atomic waVefunctions I, II; University of Alberta: Edmonton, 1971. (28) Huzinaga, S. J. Chem. Phys. 1965, 42, 1293. (29) Dick, B. G.; Overhauser, A. W. Phys. ReV. 1958, 112, 90. (30) Schober, H.; Strauch, D. J. Phys.: Condens. Matter 1993, 5, 6165. (31) Catlow, C. R. A.; Cormack, A. N.; Theobald, F. Acta Crystallogr. B 1984, 40, 195. (32) Hriljac, J. A.; Eddy, M. M.; Cheetham, A. K.; Donohue, J. A.; Ray, G. J. J. Solid State Chem. 1993, 106, 66. (33) Richardson, J. W., Jr., Pluth, J. J.; Smith, J. V.; Dytrych, W. J.; Bibby, D. M. J. Phys. Chem. 1988, 92, 243. (34) Marler, B. Zeolites 1987, 7, 393. (35) Rudolf, P. R.; Garce´s, J. M. Zeolites 1994, 14, 137. (36) van Koningsveld, H.; Jansen, J. C.; van Bekkum, H. Zeolites 1990, 10, 235. (37) Fischer, R. X.; Baur, W. H.; Shannon, R. D.; Staley, R. H.; Vega, A. J.; Abrams, L.; Prince, E. J. Phys. Chem. 1986, 90, 4414. (38) Levien, L.; Previtt, C. T.; Weidner, D. J. Am. Mineral. 1980, 65, 920. (39) Peacor, D. R. Z. Kristallogr. 1973, 138, 274. (40) Geisinger, K. L.; Spackman, M. A.; Gibbs, G. V. J. Phys. Chem. 1987, 91, 3237. (41) Grimm, H.; Dorner, B. J. Phys. Chem. Solids 1975, 36, 407. (42) Hill, J.-R. Personal communication, 1995. (43) Strauch, D.; Dorner, B. J. Phys.: Condens. Matter 1993, 5, 6149. (44) Serrano, D. P.; Li, H. X.; Davis, M. E. J. Chem. Soc., Chem. Commun. 1992, 745. (45) Fenzke, D.; Hunger, M.; Pfeifer, H. J. Magn. Reson. 1991, 95, 477. (46) Anderson, M. W.; Klinowski, J. Zeolites 1986, 6, 455. (47) Haase, F.; Sauer, J. J. Am. Chem. Soc. 1995, 117, 3780. (48) Czjzek, M.; Jobic, H.; Fitch, A. N.; Vogt, T. J. Phys. Chem. 1992, 96, 1535. (49) Klinowski, J.; Ramdas, S.; Thomas, J. M.; Fyfe, C. A.; Hartman, J. S. J. Chem. Soc., Faraday Trans. 2 1982, 78, 1025.

JP953405S