CATIVIC: Parametric quantum chemistry package for catalytic reactions: I

CATIVIC: Parametric Quantum Chemistry Package for Catalytic Reactions: I ´ NCHEZ,2 G. MARTORELL,1 C. GONZA ´ LEZ,3 F. RUETTE,1 M. SA ˜ EZ,1 A. SIERRAALTA,1 L. RINCO ´ N,4 C. MENDOZA5 R. AN 1

Laboratorio de Quı´mica Computacional, Centro de Quı´mica, Instituto Venezolano de Investigaciones Cientı´ficas, Apartado 21827, Caracas, Venezuela 2 Departamento de Quı´mica, Instituto Universitario de Tecnologı´a Federico Rivero-Palacio, Apartado 40347, Caracas, Venezuela 3 Physical and Chemical Properties Division, National Institute of Standards and Technology, Gaithersburg, MD 20899, USA 4 Facultad de Ciencias, Departamento de Quı´mica, ULA, Me´rida, Venezuela 5 Laboratorio de Fı´sica Computacional, Centro de Fı´sica, Instituto Venezolano de Investigaciones Cientı´ficas, Apartado 21827, Caracas, Venezuela Received 1 November 2002; accepted 16 May 2003 DOI 10.1002/qua.10719

ABSTRACT: A quantum chemistry package for catalytic reactions, referred to as CATIVIC and based on simulation techniques of parametric Hamiltonians, is presented. We describe in detail the computational procedures for modeling adsorption on a catalytic substrate, especially the parameterization scheme using examples of atomic Al and AlOX bonds (X ⫽ H, N, O, Si, Fe) for atomic and molecular parameters. The code features are illustrated with the adsorption of NO on models of ZSM-5 zeolite doped with Fe. © 2003 Wiley Periodicals, Inc. Int J Quantum Chem 96: 321–332, 2004 Key words: parametric method; semiempirical; catalysis; zeolite; CATIVIC

1. Introduction

T

he modeling of catalytic processes entails the calculation of several reaction steps, namely, physisorption, chemisorption, diffusion, surface reCorrespondence to: F. Ruette, e-mail: [email protected]

International Journal of Quantum Chemistry, Vol 96, 321–332 (2004) © 2003 Wiley Periodicals, Inc.

action, and desorption. In addition, surface coverage plays an important role in each of them. A realistic and accurate representation of the catalytic system often leads to calculations that are intractable for the standard ab initio quantum chemistry methods. For instance, in a catalyst formed by the dispersion of relatively small metallic clusters on a support the adsorbed clusters usually have a differ-

RUETTE ET AL. ent structure from their free or bulk states. Hence, cluster–support interactions are relevant to delineate an adsorption site, and several of the latter plus the support are required to mimic the adsorbate– surface and catalyst–support interactions. Further, model size can be of importance to evaluate activation-barrier heights and the stability of intermediates because the electric field in the adsorption site environment may be of significance in the kinetics and reactivity of the catalytic system, especially for supports of ionic nature. Given the size and complexity of such systems, it is necessary to develop quantum chemical methodologies based on parametric methods that treat with some degree of reliability catalysis in complex and amorphous systems, where it is essential to understand their bases, e.g., parametric functionals, to introduce improvements in a systematic fashion. Several parametric methods have been recently improved and successfully extended [1–7] to transition metals. The main differences between these methods and the one presented here are related to parameterization technique adopted, basic parametric functionals employed, and the conception of parametric methods. In CATIVIC, we follow the MINDO/3 [8] philosophy; i.e., a separation between atomic and diatomic parameters. A specific set of diatomic parameters for each pair of atoms based, as a first approach, on data form as diatomic molecules and atomic parameters adjusted from spectroscopic data of atoms. The form of the functionals, in this first version, corresponds to those employed in MINDO/SR [9 –12]. In addition, we introduce special tools to facilitate the analysis of adsorbate–substrate interactions, as shown below. Given that the foundations of parametric methods have been previously established [13–17], we developed a versatile computer code to treat catalytic processes based on parametric functionals. In earlier work [13, 14], comparisons between semiempirical MINDO/SR functionals [9 –12] with the analytic ab initio HF/STO-3G method led to the possibility of improving parametric methods by considering a suitable selection of electron– electron repulsion functionals. A simple parameterization scheme has also been carried out based on simulation techniques of the binding energy functionals [15]. They use as reference an exact electronic energy functional, that is, a Born–Oppenheimer and nonrelativistic configuration interaction energy functional (CI-EF) instead of the Hartree–Fock functional. Reasonable parameter transferability was found for small molecules containing light and

322

heavy atoms, e.g., hydrogen, carbon, and first-row transition metals. More recently, applications of functional analysis concepts have been employed to understand the nature of energy parametric functional in terms of the Riesz and Weiestrass’s theorems [16]. Functionals and the parameters have been constrained by considering the virial, hypervirial, and Hellmann–Feynmann theorems. Basis sets have been examined by defining optimal transformed minimum basis sets (OTMBSs, extended minimum basis set implicitly included in parametric functionals) and considering different environments involved in the electronic interactions. Finally, the general conditions of the energy functional component space and their properties have been studied by considering sesquilinear forms [17] and disjoint sets of basic functionals. The minimax and variational principles were defined in terms of a functional set that integrates the total electronic energy functional represented by singleelectron bases. Previous applications of parametric Hamiltonians in complex chemical systems have been in general successful [1–7, 18 –21]. Further, these methods have the ability to reproduce at least qualitatively the location of the critical points of the electronic and spin densities at the atomic valence shell [22]. It is well known that these points are closely related with the active sites in materials and molecules [23], in particular those related to the spin density [22, 24]. The aim of this article is to present the first version of a parametric code, referred to as CATIVIC, formally based on simulation techniques and customized for catalytic applications. A brief description of the theoretical background is presented in Section 2. Examples of atomic and molecular parameterization schemes are given in Section 3. The organic structure of the package is detailed in the context of chemisorption on a model catalytic surface (Section 4). Applications to the adsorption of NO on zeolite–Fe systems are exemplified in Section 5, including the use of tools for data analysis. Finally, comments, conclusions, and suggestions are discussed in Section 6.

2. Theoretical Foundations 2.1. PARAMETRIC HAMILTONIANS According to operator simulation techniques [14, 15], the optimal parametric Hamiltonian (Hpa) can be obtained by minimization of the distance be-

VOL. 96, NO. 4

PARAMETRIC QUANTUM CHEMISTRY PACKAGE tween the exact (Hexa) and a family of parametric operators (F) using the following expression: min Hpa僆F储Hexa ⫺ Hpa储 ⫽ minHpa僆F ⫻

冉 冘 兩共⌿ , H I

冊

1/ 2

⌿J兲 ⫺ 共⌿I, Hpa⌿J兲兩2

exa

I,J僆S

,

(1)

where Hpa properly describes the physical problem of interest, and ⌿I belongs to the subspace S of the Hilbert space (᐀) spanned by Hexa. The I and J subscripts indicate states of an atom or molecule. For a system of n-electrons, the n-electron wave function is expanded in terms of a set of linearly independent one-electron wave functions {␾␮}僆 ᐀. Note that (⌿I,Ht ⌿J) (t ⫽ exa, pa) corresponds to the energy functionals. Because our main interest is determining the state energies, Eq. (1) can be expressed as min Epa僆具F典

冉 冘 兩E

冊

1/ 2

I exa

I 2 ⫺ Epa 兩

I僆S

,

(2)

[27, 28]. The diatomic energy ␧XY can be expressed in terms of resonance energy integrals {␤X␮Y␯(RXY)}, electron–nucleus attractions {VX␮Y␯ (RXY)}, interatomic electron– electron repulsions {␥X␮Y␯ (RXY)}, and the nucleus–nucleus repulsion energy funcc tions {fXY(␣XY, ␥XY (RXY), CXY,. . .)}. Note that ⑀XY depends on the coordinates, {RXY}: ␧ XY ⫽ f共兵 ␤ X ␮Y ␯其, 兵 ␥ X ␮Y ␯其, 兵V X ␮Y ␯其, 兵 f XY 其兲.

(5)

Thus, the total energy is a functional that depends on {RXY} and parameters: J ⫽ f共兵 jX␮␯其, 兵kX␮␯其, 兵UX␮␯其, 兵 fXY其, 兵␤X␮Y␯其, E pa

兵VX␮Y␯其, 兵␥X␮Y␯其, RXY兲.

(6)

The parametric functional components of the total energy used in CATIVIC are similar to those employed in MINDO-SR [9 –12] and MINDO-3 [8]. Nevertheless, new parametric functionals, as suggested in Refs. [14, 15], are included in the code.

where 具F典 is a family of parametric energy funtionals.

3. Parameterization 3.1. ATOMIC PARAMETERS

2.2. ENERGY FUNCTIONAL In parametric methods, the total energy functional of a molecule can be partitioned in monoatomic and diatomic energy terms [25, 26]: E pa ⫽

冘␧ ⫹ 冘␧ X

X

.

XY

(3)

XY

The monoatomic energy ␧X depends on three kinds of terms: (1) core integrals {UX␮␯} that include the kinetic and electron–nucleus interactions; (2) intraatomic Coulomb electron– electron repulsions {jX␮␯}, and (3) the exchange {kX␮␯} integrals. The two latter terms are usually expressed as a function of Slater–Condon integrals (Fk and Gk). Here, the core integral and Coulomb integral are defined as UX␮␯ ⫽ (␺␮ ,Hcore ␺␯) and jX␮␯ ⫽ (␺␮␺␮,1/r12 ␺␯␺␯), respectively. Thus, the monoatomic energy is a functional that depends on constant functions for each pair of ␺␮␺␯ one-electron orbitals: ␧ X ⫽ f共兵 j X ␮␯其, 兵k X ␮␯其, 兵U X ␮␯其兲.

(4)

Interaction between atoms is given by diatomic parametric functionals that obey the ZDO approach

Because the estimate of binding energies requires the calculation of atomic energies, parameterization of functionals for single atoms is performed considering the corresponding integrals ({ jX␮␯}, {kX␮␯}, {UX␮␯}) shown in expression (4). Substraction and addition of the atomic ground-state 0 energy (Eexa ) in the first and second terms of Eq. (2) lead to an expression that depends on energy differences: min Epa僆具F典

冉 冘 兩⌬E

I exa

I僆S

冊

1/ 2

I 2 ⫺ ⌬Epa 兩

,

(7)

I I 0 I where ⌬Epa ⫽ Epa ⫺ Eexa . Values of ⌬Eexa are obtained from atomic spectroscopic data [29] because they are related to the atomic excitation energies relative to the corresponding ground state. The simulating annealing (SA) program developed by Gonza´lez [30] has been employed to locate the global minima in Eq. (7). The energy difference I between the different configurations (⌬Epa ) is comI pared with the experimental values (⌬Eexa ), taken from Moore’s tables [29]. Calculations for states associated to the electronic configurations were performed maintaining the corresponding electronic

INTERNATIONAL JOURNAL OF QUANTUM CHEMISTRY

323

RUETTE ET AL. TABLE I ______________________________________________________________________________________________ Calculated values of total energy (ENERGY), calculated (EEpa) and experimental (EEexp) excitation energies, and their differences (EEexp ⴚ EEpa) for the Al atom. State 2

P: 3s23p1 P: 3s13p2 2 S: 3s13p2 2 P: 3s13p2 1 b S : 3s2 3 b P : 3s13p1 1 b P : 3s13p1 4

a b

ENERGY (a.u.)

EEexpa (kcal/mol)

EEpa (kcal/mol)

EEexp ⫺ EEpa (kcal/mol)

⫺2.07760 ⫺1.94716 ⫺1.83911 ⫺1.82047 ⫺1.85771 ⫺1.68558 ⫺1.58685

0.0 81.88 149.70 161.40 138.02 246.07 308.04

0.0 82.97 147.97 161.95 138.03 244.94 309.15

0.0 ⫺1.09 1.73 ⫺0.55 ⫺0.01 1.13 ⫺1.11

From Ref. [29]. Positively charged atom.

configuration fixed during the self-consistent field (SCF) process. The calculated valence-electrons energy (ENERGY), the experimental excitation energies (EEexp), the theoretical energy differences (EEpa) relative to the ground state, and the (EEexp ⫺ EEpa) difference are presented in Table I for the Al atom. Results show that a small discrepancy is reached by the adjustment of seven parameters 0 [core integrals (Uss, Upp) and Slater–Condon (Fss , 0 0 2 1 Fpp, Fsp, Fpp, and Gsp)] using seven different electronic states. Selected states include the 2P(3s23p1) ground state and the 4P(3s13p2), 2S(3s13p2), and 2 P(3s13p2) excited states of the neutral atom. In addition, cation states were also considered, namely, 1 S(3s2), 1P(3s13p1), and 3P(3s13p1), because of the scarce number of states in the neutral Al atom spectrum [29]. The selection of Al cation spectrum may be relevant for modeling positive charged molecules. The core integrals (UAl␮␮ ) and one-center Coulomb and exchange integrals ( jAl␮␯, kAl␮␯) (␮ ⫽ s, p), evaluated from the Slater—Condon parameters [31], are displayed in Table II. Details of the atomic parameterization for several types of atoms (Be, C, K, V, Ga, Mo) are analyzed elsewhere [32].

3.2. MOLECULAR PARAMETERS As shown previously [15], simulation techniques of Eq. (2) for the bond formation of a diatomic molecule can be used to obtain molecular parameters. Consider the formation of the XY diatomic molecule: X ⫹ Y 3 XY; the binding energy can then be obtained by the subtraction and addition of X and Y atomic energies in both terms of Eq. (2), min BEXY XY pa 僆兵 f 其

冉冘 兩BE

XY exaJ

J

冊

1/ 2

XY 2 ⫺ BEpa 兩 J

,

where {fXY} is a family of binding energy parametric functionals for the set of J molecular states considering states K and L for atoms X and Y, respectively. Thus, BEXY j is given by the expression XY BE J3K⫹L ⫽ 具⌿ JXY 兩H XY 兩⌿ JXY 典 ⫺ 具⌿ XK 兩H X 兩⌿ XK 典

⫺ 具⌿ LY 兩H Y 兩⌿ LY 典. TABLE II ______________________________________ Aluminium atomic parameters. Parameter U␮␮ (eV) J (a.u.) K1 (a.u.)

324

ss

pp

30.5927 22.2421 0.37375 0.33623 — —

pp⬘

sp

— 0.31760 0.00932

— 0.34997 0.09873

(8)

(9)

The simulation process was performed considering several constraints: (1) Expression (8) is minimized at the equilibrium bond distance (Req); (2) at this distance, parameters and functionals are adjusted to reproduce dissociation energy; (3) the energy gradient at the Req is required to be zero; and (4) the electronic configuration and multiplicity at Req are found to be the most stable.

VOL. 96, NO. 4

PARAMETRIC QUANTUM CHEMISTRY PACKAGE III. Note that only parameters involving Al have been calculated with the method described above. Other parameters have been obtained from the MINDO-SR [9 –12] and MINDO-3 [8] methods.

4. Software Structure

FIGURE 1. Scheme for searching molecular parameters.

The algorithm employed in molecular parameter optimization is presented in Figure 1. It starts with the reading of the diatomic molecule (XY) properties (equilibrium bond distance, dissociation energy, and multiplicity) obtained from a database [33]. The latter has been built up from both experimental and theoretical data available in the literature. For example, data for AlOH [34, 35], AlON [34, 35], AlOO [35], AlOAl [34, 36], and AlOSi [37] have been used for the parameterization of the corresponding bonds. In the cases such as FeOAl, where data is not available, the parameters are computed by ab initio methods. Calculations have been carried out with the Gaussian 94 suite of programs [38] at the B3LYP-6-311 ⫹ G(d,p) level. The resulting FeOAl equilibrium bond distance, dissociation energy, and ground-state symmetry are 2.445 Å, ⫺42.9 kcal/mol, and 4⌺, respectively. Then, the molecular parameters ␣FeAl and ␤FeAl are determined using an iterative procedure until the gradient is close to zero and the bond energy near to ⫺42.9 kcal/mol, at the bond distance of 2.445 Å. With these parameters, different electronic configurations are evaluated considering single and double excitations as well as in ␣ and ␤ occupations. If a more stable configuration is found, a new search of parameters is carried out, using in some cases a fixed configuration, until the cycle is consistent. The molecular parameters for a system that contains H, N, O, Al, Si, and Fe are displayed in Table

As mentioned above, although CATIVIC is a general-purpose code it is specially targeted for catalytic reactions. It is thus outfitted with several tools to facilitate the modeling of such systems. In Figure 2 the implemented flow diagram for the adsorption step is shown. Because most catalysts have unpaired electrons, due to the presence of transition metals in their composition, optimization of the substrate spin configuration is required. The code therefore spans a range of multiplicities after geometry optimization, selecting the most stable one. Calculations starting with the aufbau configuration do not always reach the most stable state.

TABLE III _____________________________________ Set of molecular parameters employed in zeolite-Fe-NO systems. Atomic pair (XOY) FeOSia FeOAlb FeOOb FeONa FeOHa SiOHc SiONa SiOOa SiOAlb SiOSic AlOHb AlONb AlOOb AlOAlb OOOc OONc OOHc NONc NOHc HOHc

INTERNATIONAL JOURNAL OF QUANTUM CHEMISTRY

␣XY

␤XY

0.677390 0.830510 1.071625 0.815285 0.318123 0.289647 0.404333 0.536648 0.319236 0.291703 0.306647 0.448987 0.517709 0.478300 0.659407 0.458110 0.417759 0.377342 0.360776 0.244770

3.001750 4.619700 1.743737 1.595467 1.168525 0.497834 0.601946 0.692882 2.385445 0.486004 1.685439 1.129520 0.643442 6.099146 0.813430 0.991584 0.478901 1.074007 0.589380 0.788168

a

Ref. [9 –12]. This work. c Ref. [8]. b

325

RUETTE ET AL. Once a possible adsorption site orientation has been found, a potential energy curve (PEC) is calculated by means of a relaxed scan along the line that joins critical points of adsorbate and substrate. For each point on the PEC the whole system is reoptimized, except the adsorbate–surface distance. Finally, properties such as diatomic binding energies (DBEs) [42], diatomic energies (DEs), bond orders (BOs), charges, and molecular orbital correlation diagrams are computed. As mentioned elsewhere [42, 43], the binding energy between all pair of atoms (DBE) can be evaluated approximately by a single calculation of the system. Routines for energy partitioning have been implemented in the code with the following expressions: DBE共X ⫺ Y兲 ⫽ ␧XY ⫹ fX共X ⫺ Y兲␧X ⫹ fY共X ⫺ Y兲␧Y (10) and FIGURE 2. Strategy employed in calculations of interf X 共X ⫺ Y兲 ⫽ ␧ XY

action adsorbate–substrate (adsorption process).

冒冘

␧ XW ,

(11)

W⫽X

Therefore, several techniques have been implemented to find the most stable states associated with the electronic configurations of the substrate and substrate–adsorbate. After calculation of the adsorption of an arbitrary adsorbate (usually a hydrogen atom) at different sites, partitions of the substrate Fock and density matrices are stored as an initial guess for subsequent calculations of the free substrate. This usually results in substrate density relaxation, leading to an adequate initial guess for the calculation of the most stable substrate state. Further, small variations of the substrate coordinates will also contribute in this search. In the case of an adsorbate–substrate system, there is also the possibility of selecting different initial guess configurations by considering only the shifting of the molecular orbitals. Bader’s theory on reactivity can be a powerful approach to determine the location and properties of adsorption sites, as well as the optimal orientation of an adsorbate on a surface [39, 40]. Previous applications of this theory to parametric methods have proved that it is possible to reproduce qualitatively ab initio theoretical results [22]. Several modifications have been carried out with the EXTREM program [41] to improve the performance of calculations of critical points of PCs for big systems.

326

where ⑀XY and ⑀X are the diatomic and monoatomic energies, respectively, as shown in Eq. (3). The molecular orbital diagrams can be evaluated by expressing molecular orbitals of the adsorbate– substrate system in terms of molecular orbitals of the isolated adsorbate and substrate [44], i.e., one can obtain the adsorbate and substrate orbitals that participate in the formation of an orbital of the system. Thus, the eigenvector matrix (␾) of adsorbate–substrate can be spanned as well as from the basis set {␺␮} or from the set of eigenvectors of the isolated substrate and adsorbate {␾MOF}.

␾ ⫽ ␸ C ⫽ ␾ MOF C MOF ,

(12)

MOF ␾ MOF ⫽ ␸ C block .

(13)

where

C is the coefficient matrix of MO from the whole system and the CMOF matrix is defined from Eqs. (12) and (13) as MOF ⫺1 C MOF ⫽ 共C block 兲 C,

(14)

where

VOL. 96, NO. 4

PARAMETRIC QUANTUM CHEMISTRY PACKAGE MOF C block ⫽

冉C0

ads

0

冊.

Csubs

(15)

Cads and Csubs are coefficient matrices of MOs for isolated adsorbate and substrate, respectively.

5. NO Adsorption on Iron–Zeolite The modeling of catalytic processes has been of importance in recent years, in particular on zeolites doped with transition metals. A good example is the NO catalytic elimination due to atmospheric pollution from automobile exhaust and chemical industries. The modeling of reactions on zeolites has been usually carried out in systems with a small number of atoms. In this work, we carried out calculations for different zeolite cluster sizes (models A, B, and C) as shown in Figure 3. Comparisons with ab initio methods were performed beforehand for model A (Si2AlO4H8 and Si2AlO4H8Fe) using the B3LYP density functional theory (DFT) method from the Gaussian 94 code [38] with a pseudopotential for the Fe atom: a 6-311G(d, p) for N and O and 6-31G(d, p) basis sets for Al, Si, and H. Results are presented in Table IV. The comparison is reasonable with respect to the optimized geometry although in the ab initio method symmetry restrictions have been imposed. An improved adjustment can be reached, but it requires the fitting of more parameters and modifications of the core– core functionals [14, 15]. Computations have been carried out for models A, B, and C for the following multiplicities: 1, 3, 5, and 7, in all cases the most stable being 5. The optimization was performed first using the H atoms as bond saturators, and then for all the rest of the atoms keeping fixed the already optimized hydrogen atoms. A search for the most stable electronic configuration was then performed, followed by estimates of the critical points (CPs) of the Fe valence density in the A zeolite model as well as for N and O on the NO molecule. The topology of the Laplacian of the electronic density gives details of active sites on surfaces [39], i.e., the region where electronic charge is locally depleted [(3,⫹1) CPs] or concentrated [(3,⫺3) CPs]. The local (3,⫺3) and (3,⫹1) CPs of ⌬␳(r) can be employed to settle the physical basis of the Lewis acid– base reaction model, that is, regions of local charge concentration in the Lewis base react with regions of local charge depletion

(Lewis acid). Results for values of (3,⫺3) and (3,⫹1) CPs for Fe, N, and O are shown in Figure 4. The N and O atoms present (3,⫺3) CPs above these atoms [spheres labeled L in Fig. 4(a)], while Fe shows (3,⫹1) ones (spheres labeled X in Fig. 4) that are above and below the plane perpendicular to the O-Fe-O plane of the surface model. The topological analysis of the Laplacian density shows that adsorption of NOO for N and O atoms must be favored above the Fe atom in the end-bond mode rather than the side-bond mode. In the case of the substrate–adsorbate system a similar search of mutiplicities was performed (2, 4, 6, and 8 multiplicities), 4 being the most stable. Potential energy curves were calculated according to the location of the CPs of Fe on the substrate, O, and N. Results show an attractive potential well for the NO end-bond mode adsorption interaction through N or O atom above Fe, as shown in Figure 5. Adsorption by N is more stable than for the O atom. On the other hand, the side-bond interaction is not stable. These results confirm the NO adsorption mode expected by the location of the valence– shell critical points. Once the characterization of the adsorption was performed, analysis of the adsorbate–substrate interactions was carried out for the NO adsorption by the N atom, considering correlation diagrams according with Eqs. (12)–(15). The tedious selection of interacting orbitals is performed by the program and results for ␣ and ␤ electrons are shown in Figures 6(a) and 6(b). To simplify the analysis of interacting orbitals, the selection was performed considering values of the CMOF matrix coefficient [see Eq. (12)] higher or equal to 0.15. Several features are brought forth from the energy orbital location and correlation lines. Energy levels with a strong contribution from the Fe cation are more stable than NO ones. This suggests a charge transfer from NO to Fe⫹. The strongest interaction occurs between the 2s-2pz orbitals of NOO with the 4s-4pz of the substrate. In addition, the interaction of the dz2 orbitals of Fe indicates a small participation. The interaction of the frontier orbitals is mainly of ␲-type symmetry (px orbitals) of the NO and substrate. Finally, values of DBE, DE, BO, charge, and bond distances for A, B, and C models are displayed in Tables V–VII. The program is able to select the interactions according to a minimal adsorbate–substrate distance criteria. Thus, values of bond distances (FeOO, FeON, NOO) and angles (Fe-N-O) are presented in Table V for the different models

INTERNATIONAL JOURNAL OF QUANTUM CHEMISTRY

327

RUETTE ET AL.

FIGURE 3. Model employed for NO chemisorption on a doped zeolite doped with Fe: (A) (Si2AlO4H8Fe, 16 atoms); (B) (Si9AlO12H20Fe, 43 atoms); and (C) (Si31AlO43H44Fe, 120 atoms).

328

VOL. 96, NO. 4

PARAMETRIC QUANTUM CHEMISTRY PACKAGE TABLE IV _____________________________________ Comparison between CATIVIC and DFT calculations for model A with and without Fe. Distance (Å) SiOO AlOO AlOO(OH) SiOH OOH FeOO FeOAl

DFT

CATIVIC

1.61 (1.73) 1.78 (1.91) 1.78 (1.70) 1.50 (1.48) 0.96 (0.96) — (1.90) — (2.80)

1.67 (1.79) 1.80 (2.13) 1.79 (1.71) 1.50 (1.48) 0.94 (0.93) — (1.96) — (2.96)

Values in parentheses correspond to the model with Fe.

shown in Figure 3. Results depict an increase of the FeOO equilibrium bond distances due to the N–O interaction. These distances change with the model size and in model C is mainly through a single FeOO bond. An inverse effect is observed with the FeON distance, i.e., the FeON is strengthened because the FeOO is weakened. The NOO distance does not change with the size of the model and is slightly shorter than in free NOO. The Fe-N-O

FIGURE 4. Critical points of valence shell Laplacian density of Fe in Si2AlO4H8Fe and of N and O in NO.

FIGURE 5. Potential energy curves for NOO adsorption on Fe adsorption site.

angle shows some deviation from 180° due to asymmetrical bonding of Fe with O atoms. Values of BO and DE are presented in Table VI. The same trends observed in the equilibrium bond distances are reproduced for the bond strength: the weakening of a FeOO bond due to NOO chemisorption and the strengthening of the FeON bond with the increase of model size. The NOO bond is stronger than in the free molecule and it does not change with the model size. This stronger bond in the adsorbed NOO can be explained by analyzing the frontier orbitals of the NO molecule. The NO highest occupied molecular orbital (HOMO) is an antibonding orbital and a transfer of electron density to the adsorption site (Fe⫹) will strength the NOO bond. Qualitative information can be obtained from the charge changes in different atoms of the substrate and adsorbate, as presented in Table VII. The Fe⫹ atom gains electrons because they show positive values for the charge change. This value increases with model size. The same trend is observed with the zeolite oxygen atoms directly bonded to Fe. An inverse charge change (loss of electrons) is observed in the O and N atoms of NO. Adsorption of NOO produces an electronic charge transfer from NOO to Fe⫹, as was previously anticipated. On the other hand, charges on the O and N atoms of NO increase, giving a significant NOO charge separation with respect to the free molecule. The positive

INTERNATIONAL JOURNAL OF QUANTUM CHEMISTRY

329

RUETTE ET AL. TABLE V ______________________________________ Bond distances (Å) and angles for models with and without NO. Bond distances Model A B C NO (free)

Angle

FeOOa

FeON

NOO

FeONOOb

2.04 (1.96) 2.04 (1.96) 2.07 (2.06) 1.99 (1.98) 2.06 (2.03) 2.47 (2.00) —

2.02

1.15

173.7

1.98

1.14

173.4

1.95

1.16

167.0

—

1.16

—

NO is adsorbed through the N atom. a Values in parentheses correspond to systems without NO. b Oxygen of Zeolite.

charge on the N atom increases as the size of the substrate increases. These results seem to indicate that the cluster size model may play an important role in the characteristics of the adsorption site.

6. Conclusions and Comments We proceed to itemize the main conclusions and comments:

FIGURE 6. Correlation diagrams of NO interaction with Si2AlO4H8Fe. (a) ␣ orbitals. (b) ␤ orbitals.

330

1. Results indicate that parametric methods, such as CATIVIC, can be powerful in the analysis of catalytic processes, in particular the implementation of subroutines that automate calculations in different steps of surface reactions to facilitate the interpretation of quantum chemical data. 2. A natural application of this type of method is the study of complex systems that require a large number of atoms to simulate the catalytic site. 3. It is possible to systematically improve the parametric methods due to the generalization of simulation techniques. The methodology showed here may be improved including more complex basic functionals that better describe the physics of the problem. 4. The use of qualitative information of topological properties of the density Laplacian can be of importance for the location of adsorption site and adsorption mode. 5. Qualitative results of the application of CAT-

VOL. 96, NO. 4

PARAMETRIC QUANTUM CHEMISTRY PACKAGE TABLE VI _____________________________________________________________________________________________ Diatomic energies (DEs) and bond orders for models with and without NO. Bond orders a

Model A B C NO (free)

DE (a.u.) a

FeOO

FeON

NOO

FeOO

FeON

NOO

0.94 (1.03) 0.89 (1.03) 0.86 (0.94) 0.96 (1.14) 0.83 (1.01) 0.32 (1.12) —

0.75

2.07

⫺0.380

⫺1.036

0.81

2.08

⫺0.433

⫺1.043

0.80

2.02

⫺0.529

⫺1.039

—

2.05

⫺0.442 (⫺0.505) ⫺0.431 (⫺0.505) ⫺0.389 (⫺0.447) ⫺0.452 (⫺0.531) ⫺0.344 (⫺0.463) ⫺0.090 (⫺0.494) —

—

⫺0.984

Values in parentheses correspond to systems without NO. a Oxygen of Zeolite.

IVIC to chemisorption of NO on a Fe-ZSM-5 show that it is important to use an aggregate of a sufficient size to model the metal-ZSM-5 catalyst and the adsorption site. 6. Other techniques to the code have been incorporated, such as embedding procedure [45], Monte Carlo simulated annealing for global geometric optimization [46], and an orbital energy partition [26]. 7. A generalization of a similar process for diffusion, surface reaction, coverage effects, surface reconstruction, and desorption is in progress. 8. A molecular builder is already incorporated to the program to facilitate the input and analysis of results.

TABLE VII ____________________________________ Charges changes on Fe, O, O,a N, and NO for models A, B, and C due to adsorption of NO. ⌬Chargeb (a.u.) Model

Fe

O

Oa

N

NOO

A

0.35

⫺0.17

⫺0.31

⫺0.48

B

0.38

⫺0.18

⫺0.37

⫺0.55

C

0.50

⫺0.02 ⫺0.02 ⫺0.01 0.02 0.04 0.11

⫺0.05

⫺0.45

⫺0.50

a

b

Oxygen of NO. Values of ⌬Charge on atoms Fe and O are referenced with respect to isolated systems (zeolite–Fe and NO).

b

ACKNOWLEDGMENTS This research has been supported by CONICIT under Contracts G-9700667 and CONIPET 97003734. The authors also thank Lic. Luis Rodrı´guez for help with the figures.

References 1. Bredow, T.; Geudtner, G.; Jug, K. J Comput Chem 2001, 22, 861, and references therein. 2. Pearl, G. M.; Zerner, M. C. J Am Chem Soc 1999, 121, 399, and references therein. 3. Cundari, T. R.; Deng, J. J Chem Info Comput Sci 1999, 39, 376. 4. Sizova, O. V.; Baranovski, V. I. J Comput Chem 1995, 16, 586. 5. Jomoto, T.; Lin, J. J Mol Struct Theochem 2002, 577, 143, and references therein. 6. Ma, X.; Schobert, H. H. J Mol Catal A 2000, 160, 409, and references therein. 7. Voityuk, A. A. ; Ro¨sch, N. J Phys Chem 2000, 104, 4089. 8. Dewar, M. J. S.; Lo, D. H. J Am Chem Soc 1971, 93, 7201. 9. Blyholder, G.; Head, J.; Ruette, F. Theor Chim Acta 1982, 60, 429. 10. Rodrı´guez, L. J.; Ruette, F.; Castro, G. R.; Ludena, P. V.; Herna´ndez, A. J. Theor Chim Acta 1990, 77, 39. 11. Rodrı´guez, L. J.; Ruette, F.; Rossa-Brussin, M. J Mol Catal 1990, 62, 199. 12. Rinco´n, L.; Ruette, F.; Herna´ndez, A. J Mol Struct 1992, 254, 395. 13. Romero, M.; Primera, J.; Sa´nchez, M.; Sierraalta, A.; Brassesco, S.; Bravo, J.; Ruette, F. Folia Chim Theor Latina 1995, XXIII, 45. 14. Primera, J. R.; Romero, M.; Sa´nchez, M.; Sierralta, A.; Ruette, F. J Mol Struct Theochem 1999, 469, 177. 15. Romero, M.; Sa´nchez, M.; Sierraalta, A.; Rinco´n, L.; Ruette, F. J Chem Info Comput Sci 1999, 39, 543.

INTERNATIONAL JOURNAL OF QUANTUM CHEMISTRY

331

RUETTE ET AL. 16. Ruette, F.; Gonza´lez, C.; Octavio, A. J Mol Struct Theochem 2001, 537, 17. 17. Ruette, F.; Marcantognini, A. M.; Karasiev, V. J Mol Struct Theochem (in press). 18. Ruette, F.; Sierraalta, A.; Herna´ndez, A. In: Ruette, F., ed. Quantum Chemistry Approaches to Chemisorption and Heterogeneous Catalysis; Kluwer: Dordrecht, The Netherlands, 1992; p 320 and references therein. 19. Ruette, F.; Poveda, F. M.; Sierraalta, A.; Sa´nchez, M.; Rodrı´guez-Arias, A. N. J Mol Catal 1997, 119, 335. 20. Lathiotakis, N. N.; Andriotis, A. N.; Menon, M.; Connolly, J. J Chem Phys 1996, 104, 992. 21. Semiempirical Methods: Is There a Future? J Mol Struct Theochem 1997, 401 (special issue), 193. 22. Sierraalta, A.; Ruette, F.; Machado, E. Int J Quantum Chem 1998, 70, 113. 23. Bader, R. F. Atoms in Molecules: A Quantum Theory; Clarendon Press: Oxford, UK, 1990. 24. Sierraalta A.; Ruette, F. J Mol Catal A 1996, 109, 227. 25. Sa´nchez, M.; Ruette, F. J Mol Struct Theochem 1992, 254, 335. 26. Ruette, F.; Gonza´lez, C.; Sa´nchez, M. J Mol Struct Theochem 2001, 537, 17. 27. Pople, J. A.; Beveridge, D. L. Approximate Molecular Orbital Theory; McGraw-Hill: New York, 1970. 28. Neymeyr, K; Seelig, F. F. Int J Quantum Chem 1995, 53, 515, and references therein. 29. Moore, C. E. Atomic Energy Levels, vols. 1–3, National Standard Reference Data Series; National Bureau of Standards: Washington, DC, 1971. 30. Gonza´lez, C. In: Alvarez, L. J., ed. Me´todos Computacionales en la Cata´lisis; UNAM, CYTED: Mexico (in press). 31. Slater, J. C. Quantum Theory of Atomic Structure, vols. I and II; McGraw-Hill: New York, 1960. 32. Ruette, F.; Sa´nchez, M.; .Mendoza, C.; Sierraalta, A.; Martorell, G.; Gonza´lez, C. Int J Quantum Chem (in press).

332

33. Bermudez, A. J.; An˜ez, R.; Sierraalta, A.; Sa´nchez, M.; Ruette, F. (to be submitted). 34. Handbook of Chemistry and Physics, 53rd ed.; CRC Press: Cleveland, OH, 1973; p F-183. 35. Rosen, B. Selected Constants Ppectroscopic Data Relative to Diatomic Molecules; Pergamon Press: New York, 1970. 36. Koseki, S.; Gordon, M. S.; Schmidt, M. W.; Matsunaga, N. J Phys Chem 1995, 99, 12764. 37. Ornellas, F. R.; Iwata, S. Chem Phys 1998, 232, 95. 38. Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Gill, P. M. W.; Johnson, B. G.; Robb, M. A.; Cheeseman, J. R.; Keith, T.; Peterson, G. A.; Montgomery, J. A.; Raghavachari, K.; AlLaham, M. A.; Zakrzewski, V. G.; Ortiz, J. V.; Foresman, J. B.; Cioslowski, J.; Stefanov, B. B.; Nanayakkara, A.; Challacombe, M.; Peng, C. Y.; Ayala, P. Y.; Chen, W.; Wong, M. W.; Andres, J. L.; Replogle, E. S.; Gomperts, R.; Martin, R. L.; Fox, D. J.; Binkley, J. S.; Defrees, D. J.; Baker, J.; Stewart, J. P.; Head-Gordon, M.; Gonza´lez, C.; Pople, J. A. Gaussian 94, revision D.4; Gaussian, Inc.: Pittsburgh, PA, 1995. 39. Bader R. F. Atoms in Molecules: A Quantum Theory; Clarendon Press: Oxford, UK, 1990. 40. Aray, Y.; Rodrı´guez, J.; Vega, D.; Rodrı´guez-Arias, E. N. Angew Chem 2000, 112, 3968, and references therein. 41. Biegler-Konig, F. W.; Bader, R. F. W.; Tang, T. H. J Comput Chem 1982, 3, 317. 42. Ruette, F.; Poveda, F. M.; Sierraalta A.; Rivero, J. Surf Sci 1996, 349, 241. 43. Ruette, F.; Poveda, F. M. In: Ilisca, E.; Makoshi, K. Electron Processes at Solid Surfaces; World Scientific: Singapore, 1996; p 303–327. 44. Rinco´n, L.; Almeida, R. Chem Phys Lett 2000, 323, 416. 45. Poveda, F. M.; Ferna´ndez-Sanz, J.; Ruette, F. J Mol Catal A 2003, 191, 101. 46. Ruette, F.; Gonza´lez, C. Chem Phys Lett 2002, 359, 428.

VOL. 96, NO. 4