8402
J. Phys. Chem. B 1997, 101, 8402-8409
ARTICLES Molecular Dynamics Simulation Studies of Zeolite-A. 5. Structure and Dynamics of Cations in H12-A and (CH3NH3)10Na2-A Using Rigid Dehydrated Zeolite-A Frameworks Song Hi Lee* Department of Chemistry, Kyungsung UniVersity, Pusan 608-736, Korea
Sang Gu Choi Department of Industrial Safety, Yangsan Junior College, Yangsan 626-800, Korea ReceiVed: March 3, 1997; In Final Form: July 4, 1997X
In the present paper a molecular dynamics simulation technique is applied to study the local structure and dynamics of H+, CH3NH3+, and Na+ ions in rigid dehydrated zeolite-A frameworks using a simple LennardJones potential plus Coulomb potential with Ewald summation. In the H12-A zeolite system, two structures appear, depending upon the choice of the Lennard-Jones parameter, σ, for the H+ ion. For the smaller values of σ, the 12th H+ ion is located on one of the 8-ring window sites which are already occupied by three H+ ions; for the larger values of σ, it is at one of the opposite-4-ring sites with the remaining 11 H+ ions almost fixed near their initial positions. In the (CH3NH3)10Na2-A zeolite system, the main structural differences from an X-ray crystallographic report are 4-fold: no facing CH3NH3+ ions through a 6-ring window, two ions in the β-cage, the appearance of a CH3NH3+ ion on one of the opposite-4-ring sites, and the lying of CH3NH3+ ions on the planes of the 8-ring window sites. Four kinds of time correlation functions for the CH3NH3+ ions show the dynamics of the ions, reflecting the different structural arrangements of the ions well. The analyses of hydrogen bond time correlation functions for the four nonequivalent CH3NH3+ ions indicate that about 0.8, 2.9, 0.6, 2.9, 1.6, 0.9, 1.7, and 2.0 hydrogen bonds formed between I and O(2), I and O(3), II and O(2), II and O(3), III and O(1), III and O(2), IV and O(1), and IV and O(3) are retained for 0.83, 1.89, 0.63, 1.61, 0.42, 0.31, 1.11, and 1.26 ps, respectively, before a breaking of the hydrogen bond occurs, leading to a significant exchange of O atoms hydrogen-bonded to the ion.
I. Introduction Computer simulation methods are playing an increasingly important role in the study of the behavior of small adsorbates in the channels and cavities of nanoporous materials like zeolites because many details of the dynamical properties of the adsorbates in porous media are inaccessible to experimental measurements. In particular, the molecular dynamics (MD) method has an advantage over the Monte Carlo (MC) method because not only static but also dynamic properties can be obtained. In our previous studies, we have performed MD simulations of four different zeolite-A systems to investigate the structure and dynamics of adsorbates: rigid dehydrated zeolite-A,1 two rigid Ca2+-exchanged zeolite-A systems,2 rigid hydrated zeoliteA,3 and rigid NH4+-exchanged zeolite-A.4 In the first paper,1 a best-fit set of electrostatic charges was chosen from the results of the MD simulation at 298.15 K, and the Ewald summation5 was used for the long-range character of Coulomb interaction. The calculated x-, y-, and z-coordinates of Na+ ions were in good agreement with the positions determined in an X-ray crystallographic study6 within statistical errors. Their random motions in different closed cages were well described by time correlation functions, and NaI ions were found to be less diffusive than NaII and NaIII. At 600.0 K, the * To whom correspondence should be addressed. X Abstract published in AdVance ACS Abstracts, September 15, 1997.
S1089-5647(97)00769-4 CCC: $14.00
unstable NaIII ion pushes one of the nearest NaI ions into the β-cage and sits on the stable site I instead, and the captured ion in the β-cage wanders over and attacks one of the eight NaI ions. In the second paper,2 the same method was applied to study the local structure and dynamics of cations in the rigid dehydrated Ca2+-exchanged zeolite-A systems Ca6-A and Ca4Na4-A. At 298.15 K, the calculated positions of Ca2+ and Na+ ions were in poor agreement with those determined by X-ray diffraction experiments,7,8 but this was reasonably explained by the large repulsive interactions between cations. A simple harmonic oscillation of Ca2+ ions in dehydrated Ca6-A zeolite and the simultaneous random motions for the same kinds of cations in the dehydrated Ca4Na4-A zeolite due to a positional symmetry were observed in their velocity autocorrelation functions and mean square displacements. MD simulation of the fictitious Ca8-A zeolite with Ca1.5+ confirmed our results on the behavior of Ca2+ ions in dehydrated Ca6-A zeolite. In the third paper,3 rigid hydrated zeolite-A at 298.15 K was studied by MD simulation to investigate the local structure and dynamics of Na+ ions and water molecules. The results of the MD simulation showed that the eight NaI ions remained at 6-ring windows, characterized by nondiffusive vibrational motion, and that four NaIV ions showed a diffusive mobility due to loose bindings to nearest framework atoms and full hydration with R-cage water molecules. β-cage water molecules formed a distorted tetrahedron with edge lengths in good agreement with © 1997 American Chemical Society
MD Simulation Studies of Zeolite-A
J. Phys. Chem. B, Vol. 101, No. 42, 1997 8403
TABLE 1: Initial Positions (Å) of Na+ Ions Used for the Positions of H+ Ions in the H12-A Zeolite ion
x
y
z
ra
ion
x
y
z
ra
H(1) H(3) H(5) H(7) H(9) H(11)
3.693 3.693 -3.693 -3.693 6.139 0.955
3.693 -3.693 -3.693 3.693 0.955 0.823
3.693 3.693 3.693 -3.693 0.823 6.139
6.396 6.396 6.396 6.396 6.259 6.259
H(2) H(4) H(6) H(8) H(10) H(12)
-0.693 3.693 3.693 -3.693 0.823 -3.536
3.693 3.693 -3.693 -3.693 6.139 -3.536
3.693 -3.693 -3.693 -3.693 0.955 0.000
6.396 6.396 6.396 6.396 6.259 5.001
a
r is the distance of the H+ ion from the center of the simulation box.
those predicted in an X-ray diffraction experiment.9 However, there was no clear evidence for the dodecahedral arrangement of 20 water molecules predicted in the X-ray diffraction experiment.9 The dynamical behavior of R- and β-cage water was more diffusive than that of bulk water due to the lack of hydrogen bonds in the closed cages. The rotational motions of water molecules in both R- and β-cages were characterized by slower librational motions and faster rotational relaxation than those of bulk water. An MD simulation in the true cell showed very little sensitivity for the structure and dynamics of Na+ ions and water molecules to the size of the unit cell. In the fourth paper,4 we reported an MD simulation study of NH4+ ions in rigid dehydrated zeolite-A at 298.15 K using the previously determined Lennard-Jones (LJ) potential parameters and electrostatic charges for the zeolite framework atoms1 and the newly determined ones for the NH4+ ion. During the preliminary equilibration, the unstable NH4(4) ion (a 12th ion) was pushed down to near a more stable 6-ring position in the R-cage that was already associated with an NH4(1) ion (a first ion) in the β-cage, which moved to another 6-ring position in the β-cage that was already associated with an NH4(2) ion (a seventh ion) in the R-cage. Calculated x-, y-, and z-coordinates of some NH4+ ions were in good agreement with those obtained from an X-ray diffraction experiment10 except that no NH4(4) ion was found and there were six NH4(2) ions instead of 0.5 and 5.5 occupancy. The analyses of calculated interatomic distances and time correlation functions of these ions indicated that the NH4(1-1) and NH4(3) ions were associated loosely with only one O(3) atom of the 6-ring and with only one O(1) atom of the 8-ring windows, respectively, while the NH4(1-2) and NH4(2) ions were associated strongly with two or three O(3) atoms of the 6-ring windows in the R- and β-cages, respectively. The analyses of hydrogen bond time correlation functions of these ions indicated about one, two or three, three, and one hydrogen bond of each NH4(1-1), NH4(1-2), NH4(2), and NH4(3) ion was retained for 1.4, 21, 75, and 1.4 ps, respectively, before a breaking of the hydrogen bond occurred, leading to a significant exchange of O atoms hydrogen bonded to the ion. The extent of ion exchange achievable in zeolite-A may be limited by ion size. Only small cations can penetrate the single 6-rings into the β-cage, especially if a rigid framework is used. Large organic cations cannot penetrate the 8-rings into the R-cage. Cation-cation repulsion may further limit the extent of exchange.11-12 Such considerations were used to explain why only Cs7Na5-A and Cs7K5-A were obtained when exhaustive Cs+ exchange was attempted with zeolites Na12-A and K12A.11-12 The order of decreasing selectivity for univalent ions in zeolite-A is Ag+ > Tl+ > Na+ > K+ > NH4+ > Rb+ > Li+ > Cs+ and for divalent ions, Zn2+ > Sr2+ > Ba2+ > Ca2+ > Co2+ > Ni2+ > Cd2+ > Hg2+ > Mg2+.13 CH3NH3+ ions are monopositive and relatively large, so only 10 CH3NH3+ ions would exchange with Na+ ions in zeolite-A.14 Formation of the hydroxyl group in place of the exchanged cation in the zeolite framework would necessarily result in the disruption of an Si, Al-O bond.15 However, there is evidence that the H+ form of zeolite-A framework can be reconstituted.16
A Ag12-A crystal lost its diffraction pattern upon exposure to H2 gas, apparently because Ag metal migrated out of the zeolite leaving H+ inside. Subsequent heating of this crystal in an O2 atmosphere caused most of the Ag atoms to be readsorbed and the single-crystal diffraction pattern to be fully restored. Another simple method to prepare H12-A involves exchange of the cation by aqueous NH4+ ions.15 Subsequent thermal treatment of the NH4+-exchanged zeolite results in the liberation of NH3 and the formation of H+ ions on the framework structure. To continue these MD simulation studies of zeolite-A, we have chosen to study next the rigid H+-exchanged and CH3NH3+-exchanged zeolite-A systems [H12-A and (CH3NH3)10Na2-A]. The primary purpose of this work is to investigate the local structure and dynamics of H+ ions and of CH3NH3+ and Na+ ions in the rigid dehydrated zeolite-A using the previously determined potential parameters of the zeolite framework atoms1 and the Ewald summation technique5 with the determination of the Lennard-Jones (LJ) parameters and electrostatic charges for the H+ and CH3NH3+ ions. In section II we present the molecular models and MD simulation method. We discuss our simulation results in section III and present the concluding remarks in section IV. II. Molecular Models and Molecular Dynamics Simulations H12-A and (CH3NH3)10Na2-A zeolite frameworks are modeled by the pseudo cell, (SiAlO4)12, using the space group Pm3m (a ) 12.2775 and 12.280 Å, respectively) which contain 12 H+ ions, or 10 CH3NH3+ and 2 Na+ ions, respectively. The zeolite frameworks are assumed to be rigid. The H12-A zeolite framework atoms are fixed in space at the positions determined by the X-ray diffraction experiment of Pluth and Smith6 for the dehydrated Na12-A zeolite, and the dehydrated (CH3NH3)10Na2-A zeolite framework atoms are fixed at those determined by the X-ray diffraction experiment of Jeong, Kim, and Seff,14 respectively. The LJ parameters and the electrostatic charges on Si and Al framework atoms are assumed to be equal because the Ewald summation for the pseudo cell is valid with this assumption. The cation positions obtained from these X-ray diffraction experiments are also used for the initial positions of the cations. Table 1 lists the initial Na+-ion positions used for the H+ ions in the H12-A zeolite, and Table 2 lists those of Na+ ions, and N and C atoms of the CH3NH3+ ions in the (CH3NH3)10Na2-A. In H12-A, eight and three H+ ions occupy each of the eight 6-ring and three 8-ring window sites, respectively, and the 12th H+ ion occupies one of the opposite-4-ring sites. In (CH3NH3)10Na2-A, two Na+ ions occupy 6-ring window sites diagonally, six CH3NH3+ ions occupy the remaining six 6-ring window sites with the C-N dipolar vectors pointing at the centers of the 6-ring windows, and a C-N vector (the first CH3NH3+ ion) in the β-cage faces the C-N vector of the second CH3NH3+ ion. The remaining three CH3NH3+ ions occupy 8-ring window sites also with the C-N dipolar vectors pointing at ideally the centers of the 8-ring windows. Figure 1 shows
8404 J. Phys. Chem. B, Vol. 101, No. 42, 1997
Lee and Choi
TABLE 2: Initial Positions (Å) of Na+ Ions and N and C Atoms of the CH3NH3+ Ions in the (CH3NH3)10Na2-A Zeolite ion/atom
x
y
z
ra
ion/atom
x
y
z
ra
Na(1) C(1) C(2) C(3) C(4) C(5) C(6) C(7) C(8) C(9) C(10)
3.770 -5.767 -2.185 2.185 -2.185 2.185 2.185 -2.185 4.766 0.000 0.000
-3.770 -5.767 -2.185 -2.185 -2.185 2.185 2.185 2.185 0.000 4.766 0.000
3.770 -5.767 -2.185 -2.185 2.185 2.185 -2.185 2.185 0.000 0.000 4.766
6.530 9.988 3.785 3.785 3.785 3.785 3.785 3.785 4.766 4.766 4.766
Na(2) N(1) N(2) N(3) N(4) N(5) N(6) N(7) N(8) N(9) N(10)
-3.770 -4.926 -3.026 3.026 -3.026 3.026 3.026 -3.026 6.222 0.000 0.000
3.770 -4.926 -3.026 -3.026 -3.026 3.026 3.026 3.026 0.000 6.222 0.000
-3.770 -4.926 -3.026 -3.026 3.026 3.026 -3.026 3.026 0.000 0.000 6.222
6.530 8.532 5.241 5.241 5.241 5.241 5.241 5.241 6.222 6.222 6.222
a
r is the distance from the center of the simulation box.
Figure 1. Stereoplot of the initial (CH3NH3)10Na2-A system viewed from the R-cage. H atoms are omitted.
TABLE 3: Lennard-Jones Parameters, σ (Å) and E (kJ/mol), and Electrostatic Charges (e) Used in This Work ion/atom
σ
charge
H+ ion Na+ ion N (CH3NH3+ ion) C (CH3NH3+ ion) HN (CH3NH3+ ion) HC (CH3NH3+ ion) Al (dSi) atom O(1) atom O(2) atom O(3) atom
a 1.776 1.824 1.908 0.6 1.1 4.009 2.890 2.890 2.890
0.0657 20.8466 0.7113 0.4577 0.0657 0.0657 0.5336 0.6487 0.6487 0.6487
0.55 0.55 -0.1062 -0.0206 0.1630 0.0626 0.6081 -0.4431 -0.4473 -0.4380
a
The value is not determined and HX is the hydrogen atom attached to the X atom.
the stereoplot, by ORTEP,17 of this initial (CH3NH3)10Na2-A structure viewed from the R-cage. The interaction potential for the cations and framework atoms is given as a sum of Lennard-Jones and Coulomb potentials. The LJ parameters, σ and , and electrostatic charges of the framework atoms and Na+ ions were determined in the previous work1 by comparing the interatomic distances of the dehydrated zeolite-A system. Those for the H+ and CH3NH3+ ions must be newly determined. For H+, the value for is chosen as 0.06569 kJ/mol from a previous study,18 but the value for σ is not determined until the preliminary MD simulations of the system are carried out for several values of σ and the results of energetics and structure are compared. 0.55 e is used for the electrostatic charge of H+. For CH3NH3+, a rigid Bernal-Fowler-like model for liquid ammonia19 is modified and used. The bond lengths C-N, C-H, and N-H and the bond angles ∠CNH () ∠NCH) and ∠HCH () ∠HNH) are fixed at 1.456, 1.059, and 1.033 Å and 109.5 and 107.9°, respectively. The total electrostatic charge on CH3NH3+ should be equal to 0.55 e, as for the Na+ ion in the
previous work,1 and is distributed as -0.1062, -0.0206, 0.1630, and 0.0626 e on N and C atoms and the H atoms attached to N and C atoms, respectively, according to ref 17. The LJ parameters and electrostatic charges used in this work are shown Table 3. The potential parameters used for the rigid zeolite-A framework atoms are already obtained in our previous study.1 The long-range Coulomb potential between the electrostatic charges on the cations and framework atoms is treated with the Ewald summation, which was discussed in detail in section II of refs 1 and 3. The interaction potential for cations and framework atoms is given as a sum of Lennard-Jones and Coulomb potentials between ions (or interaction sites of ions) and framework atoms as follows
1 Na
+
Φ)
2
Na+
CH3NH3 Na+
Na+ frame
φij + ∑ ∑ φij + ∑ ∑ ξij + ∑∑ i*j i j i j CH3NH3 frame
CH3NH3
∑i ∑j ξij + ∑i
CH3NH3
*
∑j
θij (1)
where the interaction energies φij, ξij, and θij are expressed as
φij ) 4ij
i
ξij )
∑k
[( ) ( ) ] σij rij
12
-
σij rij
6
+
qiqj rij
{ [( ) ( ) ] } σkj
4kj
σkj
12
-
rkj
qkqj
6
+
rkj
(2)
rkj
(3)
MD Simulation Studies of Zeolite-A i
θij )
j
∑k ∑l
{ [( ) ( ) ] } σkl
4kl
σkl
12
-
rkl
qkql
6
+
rkl
J. Phys. Chem. B, Vol. 101, No. 42, 1997 8405
rkl
(4)
with k and l representing the interaction sites of the CH3NH3+ ion. The LJ parameters and electrostatic charges of framework atoms were determined in a previous study1 by comparing the interatomic distances of the dehydrated zeolite-A system. A canonical ensemble of fixed N (number of particles), V (volume of fixed zeolite framework), and T (temperature) is chosen for the simulation ensemble. Gauss’s principle of least constraint20 is used to maintain the system at a constant temperature. The ordinary periodic boundary condition in the x-, y-, and z-direction and minimum image convention are applied for the Lennard-Jones potential with a spherical cutoff of radius equal to the half of each simulation box length. Gear’s fifth order predictor-corrector method21 is used to solve the equations of translational and rotational motions with a time step of 2.00 × 10-16 s. For the internal orientation of the rigid geometry of the CH3NH3+ ion, Evans’s quaternion method22 is used. The equilibrium properties are averaged over five blocks of 100 000 time steps, for a total of 500 000 time steps after 500 000 time steps to reach an equilibrium state. The configuration of each ion is stored every 5 time steps for further analysis. III. Results and Discussion A. H12-A. In order to determine the best-fit value for the Lennard-Jones parameter, σ, for the H+ ion in H12-A, 17 preliminary MD simulations were carried out for several values of σ, varying from 0.4 to 2.0 Å. As the value of σ increased, the total potential energy of the system increased monotonically. However, there were two distinct structural differences in the arrangement of the 12 H+ ions. For the larger values of σ (0.8 -2.0 Å), the 12th H+ ion is at one of the opposite-4-ring sites, which means that all ions are almost fixed near their initial positions, and for the smaller values of σ (0.4-0.7 A), the 12th H+ ion moves from one of the opposite-4-ring sites to one of the 8-ring window sites which are already occupied by 3 H+ ions (the 9th, 10th, and 11th) with the remaining 11 H+ ions almost fixed near their initial positions. For the second case, there are two H+ ions (the 9th and 12th) on the same 8-ring window; in the first case, the H+ ions are arranged nearly the same as the Na+ ions in dehydrated Na12-A.1 When the value of σ is reduced to 0.2 Å, a peculiar arrangement of 12 H+ ions is observed: two H+ ions occupy one of the 6-ring window sites, one inward into the R-cage and the other outward into the β-cage. We have selected two values of σ for the H+ ion, 0.6 and 1.1 Å, for our analysis of the energetic and structural properties of H+ ions in the H12-A zeolite system. These values of σ were suggested in a previous study,18 and the variance for the other LJ parameter, , is omitted for simplicity. Table 4 shows the average potential energy and average x-, y-, and z-coordinates of each H+ ion for σ ) 0.6 and 1.1 Å. The potential energy of each H+ ion for σ ) 0.6 Å is more negative than that for σ ) 1.1 Å. The relative positions of the first 11 H+ ions are almost the same; only the 12th ion is different. According to a direct determination of proton positions in D-Y and H-Y zeolite systems by a neutron powder diffraction experiment,23 the average O-H distance found crystallographically is 1.05 Å. From the analysis of our MD simulation results, the shortest O-H distance is 0.91 Å and the average O-H distance is about 1.17 Å in the case of σ ) 0.2 Å. However, as the value of σ increases, the O-H distance is monotonically
increased: for example, the average O-H distances for σ ) 0.3, 0.4, 0.5, 0.6, and 0.7 Å are 1.32, 1.43, 1.51, 1.55, and 1.58 Å, respectively. For the larger values of σ, the O-H distances are longer. B. (CH3NH3)10Na2-A. It is worth reporting the initial behavior of CH3NH3+ ions in the preliminary simulation run for equilibration. During that run of 500 000 time steps, the first CH3NH3+ ion in the β-cage, facing the second CH3NH3+ ion, moved to another 6-ring site in the β-cage to face the sixth CH3NH3+ ion. Then the sixth CH3NH3+ ion moved to one of the opposite-4-ring sites. Meanwhile, the fourth CH3NH3+ ion in the R-cage went through the 6-ring window and rotated its C-N vector so that it continued to point at the center of its 6-ring window. Figure 2 shows the stereoplots, by ORTEP,17 of (a) two Na+ and ten CH3NH3+ ions viewed from the R-cage and (b) two Na+ and six CH3NH3+ ions from the β-cage of the rigid (CH3NH3)10Na2-A zeolite system at 298.15 K after 500 000 time steps. The average positions of the Na+ ions, and of the C and N atoms of the CH3NH3+ ions, are given in Table 5. The final arrangement of cations in the (CH3NH3)10Na2-A zeolite system is different from the initial arrangement (discussed in section II): two Na+ ions occupy two 6-ring sites diagonally, two CH3NH3+ ions occupy two 6-ring sites diagonally in the β-cage, four CH3NH3+ ions occupy the remaining four 6-ring sites in the R-cage, three CH3NH3+ ions occupy three 8-ring sites, and the last one occupies one of the two opposite-4-ring sites in the R-cage between 6-ring sites which are occupied by a Na+ ion and a CH3NH3+ ion in the β-cage. All the C-N dipolar vectors at the 6-ring sites are pointing at the centers of the 6-ring windows, but the three C-N vectors at the 8-ring sites lie on the 8-ring planes because of the large area of the window. According to the positions of the CH3NH3+ ions, the ions can be classified into four typesstwo CH3NH3+ ions in the β-cage as type I, four on 6-ring window sites in the R-cage as type II, three on 8-ring window sites as type III, and the last ion opposite a 4-ring site as type IV. The ions of the same type show similar dynamic behaviors as discussed below. In NH4-A, the X-ray diffraction experiment10 found occupancies of 3.0, 5.5, 3.0, and 0.5 for NH4(1) (near a 6-ring in the β-cage), NH4(2) (near a 6-ring in the R-cage), NH4(3) (near an 8-ring), and NH4(4) (opposite a 4-ring), respectively. However, an MD simulation study4 reported three NH4(1), six NH4(2), and three NH4(3), but no NH4(4). It is interesting that in the case of the (CH3NH3)10Na2-A zeolite system the single-crystal X-ray diffraction experiment14 found occupancies of 1, 6, and 3 for types I, II, and III, respectively. Our results show occupancies of 2, 4, 3, and 1 for I, II, III, and IV and a different orientation for the type III ions. It is worth mentioning the disappearance of the opposite-4ring site in an MD simulation result4 as compared to the experimental result10 for NH4-A and the appearance of this site in this MD simulation result as compared to the experimental result14 for (CH3NH3)10Na2-A. The main structural difference between our MD simulation result and the experiment is 4-fold: no facing CH3NH3+ ions through a 6-ring window, two CH3NH3+ ions in the β-cage, the appearance of a CH3NH3+ ion on the opposite-4-ring site, and the lying of CH3NH3+ ions on the planes of the 8-ring window sites. Average potential energies for the Na+ and CH3NH3+ ions at 298.15 K for 500,000 time steps (100 ps) are listed in Table 6. The order of stability of these ions is Na+ followed by CH3NH3+ types II, III, I, and IV, with almost no difference between types I and III. In the NH4-A zeolite system, the order was
8406 J. Phys. Chem. B, Vol. 101, No. 42, 1997
Lee and Choi
TABLE 4: Average Potential Energy (kJ/mol) and Average Positions (Å) of Each H+ Ion for σ ) 0.6 and 1.1 Å at 298.15 Ka σ ) 0.6 Å ion H(1) H(2) H(3) H(4) H(5) H(6) H(7) H(8) H(9) H(10) H(11) H(12) a
σ ) 1.1 Å
E
x
y
z
r
ion
-149.0 (7.9 -158.1 (5.3 -152.0 -146.2 -153.1 -147.7 -159.3 -151.1 -135.5 -139.3 -138.2 -123.5
4.487 (0.280 -3.381 (0.326 4.053 2.815 -2.994 2.842 -3.067 -4.343 5.912 -0.728 1.995 -5.363
3.091 (0.366 2.891 (0.280 -4.298 4.121 -3.282 -4.443 4.720 -2.937 1.894 5.861 -0.674 -2.784
3.675 (0.371 4.718 (0.230 3.015 -3.973 4.763 -3.586 -3.294 -3.880 -1.270 1.945 5.961 0.465
6.572 (0.592 6.484 (0.487 6.633 6.379 6.513 6.378 6.522 6.522 6.337 6.217 6.322 6.060
H(1) H(2) H(3) H(4) H(5) H(6) H(7) H(8) H(9) H(10) H(11) H(12)
E
x
y
z
r
-139.3 (3.7 -143.8 (5.1 -141.5 -140.1 -133.7 -141.3 -145.4 -129.7 -122.3 -120.5 -117.7 -108.5
4.237 (0.582 -3.706 (0.509 4.080 3.259 -3.165 3.269 -3.254 -3.327 5.899 1.283 0.245 -4.391
3.182 (0.201 3.237 (0.389 -3.165 4.322 -3.431 -3.583 4.465 -3.284 0.848 5.922 1.064 -4.362
3.698 (0.526 4.092 (0.548 3.720 -3.621 4.348 -4.098 -3.255 -4.332 -1.645 1.344 5.951 0.040
6.461 (0.810 6.400 (0.843 6.364 6.513 6.379 6.350 6.413 6.373 6.183 6.207 6.050 6.189
The omitted standard deviations are of the same order of magnitude as the standard deviation reported.
Figure 2. Stereoplots of (a) 10 CH3NH3+ and 2 Na+ ions viewed from the R-cage and (b) 6 CH3NH3+(1, 2, 3, 4, 5, and 7) and 2 Na+ ions from the β-cage of the rigid dehydrated (CH3NH3)10Na2-A zeolite system at 298.15 K after 500 000 time steps.
NH4(3) (8-ring), NH4(2) (6-ring, R-cage), and NH4(1) (6-ring, β-cage). The higher energy of CH3NH3(1) and Na(2) as compared to the energies of CH3NH3(4) and Na(1) of the same type is attributed to the closeness of the unstable CH3NH3(6) ion. The slightly higher energy of CH3NH3(3) is also due to the CH3NH3(1) and CH3NH3(6) ions near to it. Average interatomic distances between N and C atoms of the CH3NH3+ ions, Na+ ions, and O atoms of the zeolite-A framework are listed in Table 7. The long N(I)-N(II) distance from our MD simulation result against the experimental result indicates the absence of two facing CH3NH3+ ions in the simulation. The Na-O distances tell us that the Na+ ions are located almost at the center of the 6-ring windows since the distance between the center of the 6-ring window and O(2) is 2.965 Å, and the distance between the center and O(3) is equal to 2.391 Å which is the experimental value. The short distances of N(I)-O(3) and N(II)-O(3) obtained from our MD simulation result against the experimental result indicate that the C-N dipolar vectors of the CH3NH3+ ions are perpendicular to and
come closer to the center of the 6-ring windows in the simulation than that predicted by the experiment.14 The short distances of N(III)-O also indicate the proximity of the ions to the O atoms, not perpendicular but lying on the planes of the 8-ring window sites which can be seen from almost the same distances of N(III)-O and C(III)-O. This kind of arrangement of the CH3NH3+ ion is also found in the distances of N(IV)-O and C(IV)-O at the opposite-4-ring site. Now we turn our attention to the dynamics of Na+ and CH3NH3+ ions. Average normalized velocity autocorrelation (VAC) functions and mean square displacements (MSD) of Na+ and CH3NH3+ ions are plotted in Figures 3 and 4. The velocities of Na+ ions vary very rapidly in the 6-ring window sites, but they diffuse not far away from their initial sites as seen in Figure 4. The CH3NH3+ ions of types I and II in the 6-ring window sites but in different cages have almost the same behavior in both the VAC function and the MSD, and when compared with those of Na+ in the 6-ring window sites, the VAC function is very different due to the mass of Na+, while the MSD has almost
MD Simulation Studies of Zeolite-A
J. Phys. Chem. B, Vol. 101, No. 42, 1997 8407
TABLE 5: Average Positions (Å) of Na+ Ions and N and C Atoms of the CH3NH3+ Ions in the (CH3NH3)10Na2-A Zeolitea type
ion/atom
x
y
z
rb
Na(1)
3.701 (0.182 4.829 (0.312 -2.927 (0.347 2.928 -4.604 2.963 0.107 -2.893 6.133 0.062 0.013
-3.671 (0.193 4.687 (0.315 -2.832 (0.357 -2.761 -4.561 2.680 2.943 2.719 0.033 -5.893 -0.014
3.668 (0.199 -4.648 (0.361 -2.684 (0.354 -2.716 4.574 2.725 -3.004 2.771 0.012 0.155 5.844
6.374 (0.332 8.179 (0.572 4.878 (0.611 4.855 7.932 4.836 4.207 4.842 6.133 5.895 5.844
I
C(1)
II
C(2)
II I II IV II III III III
C(3) C(4) C(5) C(6) C(7) C(8) C(9) C(10)
ion/atom Na(2) N(1) N(2) N(3) N(4) N(5) N(6) N(7) N(8) N(9) N(10)
x
y
z
rb
-3.930 (0.176 4.022 (0.179 -3.669 (0.181 3.643 -3.791 3.618 0.020 -3.641 6.135 -0.020 -0.019
3.836 (0.188 3.877 (0.194 -3.501 (0.222 -3.491 -3.749 3.503 3.655 3.501 0.042 -5.945 -0.093
-3.838 (0.167 -3.866 (0.243 -3.501 (0.206 -3.509 3.752 3.489 -3.655 3.484 0.048 0.092 5.953
6.700 (0.307 6.794 (0.359 6.162 (0.353 6.146 6.520 6.127 5.169 6.136 6.135 5.946 5.954
a
The omitted standard deviations are of the same order of magnitude as the standard deviations reported. b r is the distance from the center of the simulation box.
TABLE 6: Average Potential Energies (kJ/mol) of Na+ and CH3NH3+ Ions at 298.15 K for 500 000 Time Steps (100 ps) type
I II
III IV
ion
energy
〈energy〉
Na(1) Na(2) CH3NH3(1) CH3NH3(4) CH3NH3(2) CH3NH3(3) CH3NH3(5) CH3NH3(7) CH3NH3(8) CH3NH3(9) CH3NH3(10) CH3NH3(6)
-147.4 ( 5.5 -130.8 ( 5.0 -108.8 ( 5.9 -119.6 ( 5.9 -121.9 ( 2.4 -118.7 ( 4.1 -128.0 ( 3.4 -125.3 ( 3.4 -111.8 ( 3.3 -109.7 ( 4.6 -109.6 ( 2.8 -74.5 ( 5.6
-139.1 ( 8.3 -114.2 ( 0.4 -123.5 ( 3.5
-110.4 ( 1.0 -74.5 ( 5.6
TABLE 7: Average Interatomic Distances (Å) distance
expta
this work
distance
expta
this work
Na-C(I) 3.96 ( 0.23 Na-O(3) 2.39(1) 2.41 ( 0.04 Na-C(IV) 4.22 ( 0.06 Na-O(2) 2.98 ( 0.04 Na-N(I) 4.58 ( 0.21 N(I)-O(3) 3.15(8) 2.40 ( 0.06 Na-N(II) 5.06 ( 0.07 C(I)-O(3) 2.81 ( 0.09 Na-N(III) 5.24 ( 0.06 N(I)-O(2) 2.98 ( 0.07 Na-N(IV) 3.96 ( 0.02 C(I)-O(2) 3.15 ( 0.21 C(I)-C(I) 5.19 ( 0.01 N(II)-O(3) 2.72(1) 2.44 ( 0.08 C(II)-C(II) 4.29(3) C(II)-O(3) 2.97 ( 0.10 C(II)-C(III) 4.08(5) 5.08 ( 0.07 N(II)-O(2) 2.99 ( 0.09 C(III)-C(IV) 4.66 ( 0.03 C(II)-O(2) 3.36 ( 0.10 C(I)-N(II) 4.47 ( 0.02 N(III)-O(1) 3.30(2) 2.83 ( 0.35 C(II)-N(III) 5.05 ( 0.09 C(III)-O(1) 2.97 ( 0.15 C(IV)-N(I) 4.12 ( 0.01 N(III)-O(2) 3.55(1) 3.30 ( 0.11 C(IV)-N(III) 5.17 ( 0.04 C(III)-O(2) 3.41 ( 0.05 N(I)-N(II) 3.34(5) 5.01 ( 0.05 N(IV)-O(3) 2.45 ( 0.02 N(I)-N(III) 5.13 ( 0.07 C(IV)-O(3) 2.52 ( 0.04 N(I)-N(IV) 4.01 ( 0.02 N(IV)-O(1) 3.16 ( 0.06 N(II)-N(II) 5.01 ( 0.03 C(IV)-O(1) 3.19 ( 0.03 N(II)-N(III) 5.17 ( 0.02 N(III)-N(IV) 4.71 ( 0.02 a
Figure 3. Normalized velocity autocorrelation functions of the four types of CH3NH3+ ions at 298.15 K
Reference 14.
the same pattern. The CH3NH3+ ions of type III show the slowest VAC function and the largest MSD due to the wide 8-ring window sites. The unstable CH3NH3+ ion of type IV does not show any peculiar behavior in both the VAC and the MSD except rather rapidly changing VAC with relatively large MSD. The self-diffusion coefficients of Na+ and CH3NH3+ ions calculated from the slopes of MSDs in the range 0.05-0.25 ps by the Einstein relation are 1.49, 2.52, 2.15, 4.51, and 2.46 × 10-4 cm2/s, respectively. In Figures 5 and 6, angular velocity autocorrelation(AVAC) and cosine autocorrelation (CAC) functions of four types of CH3NH3+ ions are plotted. The AVAC function of type III shows a relatively slow change of angular velocity due to the
Figure 4. Mean square displacements of the four types of CH3NH3+ ions at 298.15 K.
wide 8-ring window sites which is very different from those of the remaining three types. The CAC function, 〈cos φ(t)〉, of a rigid molecule is defined as 〈cos φ(t)〉 ) 〈ux(0)ux(t) + uy(0)uy(t) + uz(0)uz(t)〉 where φ(t) is the angle between the molecular axis and the z-axis at time t and ui(t) is the i-th component of the direction vector of the molecular axis with unit magnitude at time t. 〈cos φ(t)〉 is also related to the rotational relaxation time: 〈cos φ(t)〉 ) exp(-t/τ). Calculated τ’s for four types of CH3NH3+ ions in the range 0-1 ps for types I, II, III, and IV are 2.43, 1.89, 0.87, and 1.47 ps, respectively.
8408 J. Phys. Chem. B, Vol. 101, No. 42, 1997
Lee and Choi
Figure 8. Normalized hydrogen bond time correlation functions of types III and IV of CH3NH3+ ions at 298.15 K. Figure 5. Normalized angular velocity autocorrelation functions of the four types of CH3NH3+ ions at 298.15 K.
Figure 6. Normalized cosine autocorrelation functions of the four types of CH3NH3+ ions at 298.15 K.
Figure 7. Normalized hydrogen bond time correlation functions of types I and II of CH3NH3+ ions at 298.15 K.
To estimate the lifetime of hydrogen bonds between H atoms (attached to the N atom) of CH3NH3+ ions and O atoms of the zeolite-A framework, we consider the hydrogen bond time correlation function defined by
〈Ri(t)〉 )
〈
θi(0) θi(t)
〉
θi(0) θi(0)
(5)
where θi(t) is the Heaviside unit step function which is 1 if the distance between the H atoms of the NH4+ ions and O atoms of zeolite-A is less than or equal to 2.1 Å at time t and 0 otherwise, and the subscript i indicates the i-th hydrogen bond. The zero value of θi(0) is excluded in the average of 〈Ri(t)〉. The 〈Ri(t)〉 functions are determined by the strength of the hydrogen bonds and the stability of their sites. This time correlation function is somewhat similar to the residence time correlation function of water near a silica surface24 and that of water in the first shell around an ion.25
Average number, distance, and interaction energies of hydrogen bonds are collected in Table 8. This table gives a good summary for the behavior of hydrogen bonds for each type of CH3NH3+ ion: Type I ions in the β-cage and type II ions in the R-cage show almost the same behaviors2.9 short hydrogen bonds to O(3) atoms in the 6-ring window sites and 0.6-0.8 relatively longer H-bond to O(2) atom. Type III ions are bound more strongly to O(1) atoms in the 8-ring window sites than to O(2) atoms. In the case of type IV ions, the degree of H-bonding to O(1) and O(3) atoms is nearly the same due to the small size of the 4-ring window site. The ratio of the Lennard-Jones potential energy to the Coulomb energy is about 0.2-0.3%, which tells us that the hydrogen bonds are characterized mostly by the ionic aspect rather than the covalent one. Figures 7 and 8 show the time dependence of 〈Ri(t)〉 for the hydrogen bonds of four different types of CH3NH3+ ions at 298.15 K. A characteristic decay time, τ, of the hydrogen bond is obtained by fitting the time correlation function to an exponential decay 〈Ri(t)〉 ) exp(-t/τ). The characteristic decay times of H-bonds for each type of CH3NH3+ ion calculated from the time correlation function of the H-bond (the fitting region: 0-1 ps for types I, II, and IV and 0-0.5 ps for type III) are listed in Table 8. The similar behavior of H-bonds between type I-O(2) (or O(3)) and type II-O(2) (or O(3)), and that between type IV-O(1) and type IV-O(3) are also observed. In the case of type III, the loosely bound H-bond of type IIIO(2) decays more rapidly than that of type III-O(1). The results of our simulation indicate that about 0.8, 2.9, 0.6, 2.9, 1.6, 0.9, 1.7, and 2.0 hydrogen bonds formed between I and O(2), I and O(3), II and O(2), II and O(3), III and O(1), III and O(2), IV and O(1), and IV and O(3), are retained for 0.83, 1.89, 0.63, 1.61, 0.42, 0.31, 1.11, and 1.26 ps, respectively, before a breaking of the hydrogen bond occurs leading to a significant exchange of O atoms hydrogen bonded to the ion. IV. Concluding Remarks An analysis of the structural features for the rigid dehydrated H+-exchanged zeolite-A has shown two distinct structural differences of the arrangement of 12 H+ ions according to the choice of an LJ parameter, σ. For the larger values of σ, the 12th H+ ion occupies an opposite-4-ring site. For the smaller values of σ, the 12th H+ ion moves from the opposite-4-ring site to one of the 8-ring window sites which are already occupied by 3 H+ ions (9th, 10th, and 11th). The first case has almost the same arrangement of H+ ions as that of the Na+ ions in dehydrated Na12-A.1 When the value of σ is further reduced to 0.2 Å, a third arrangement of 12 H+ ions is observedstwo H+ ions in the same 6-ring window site, one inward into the R-cage and the other outward into the β-cage. On the basis of the results obtained from our MD simulation of rigid dehydrated CH3NH3+-exchanged zeolite-A, the follow-
MD Simulation Studies of Zeolite-A
J. Phys. Chem. B, Vol. 101, No. 42, 1997 8409
TABLE 8: Average Number, Distance (Å), Interaction Energies (kJ/mol), and Characteristic Decay Time (ps) of hydrogen bonds CH3NH3+ ion-O atom
number of H-bond
distance of H-bond
Coulomb energy
LJ energy
decay time (τ)
I-O(2) I-O(3) II-O(2) II-O(3) III-O(1) III-O(2) IV-O(1) IV-O(3)
0.80 ( 0.09 2.90 ( 0.01 0.64 ( 0.04 2.89 ( 0.01 1.60 ( 0.11 0.92 ( 0.17 1.71 ( 0.06 1.99 ( 0.11
1.94 ( 0.01 1.72 ( 0.01 1.93 ( 0.01 1.70 ( 0.01 1.80 ( 0.00 1.84 ( 0.02 1.85 ( 0.03 1.81 ( 0.01
-26.51 ( 0.33 -32.42 ( 0.38 -26.91 ( 0.19 -32.17 ( 0.18 -30.32 ( 0.01 -29.48 ( 0.49 -30.30 ( 0.26 -29.86 ( 0.18
-0.045 ( 0.002 -0.091 ( 0.003 -0.048 ( 0.001 -0.097 ( 0.001 -0.072 ( 0.000 -0.065 ( 0.003 -0.062 ( 0.005 -0.071 ( 0.002
0.83 1.89 0.63 1.61 0.42 0.31 1.11 1.26
ing structural and dynamic features are seen. During the preliminary equilibration, some CH3NH3+ ions make large moves. The final positions of the CH3NH3+ ions are much different from those of the X-ray crystallographic experiment: 14 no facing CH NH + ions through a 6-ring window, two ions 3 3 in the β-cage, the appearance of a CH3NH3+ ion at an opposite4-ring site, and the lying of CH3NH3+ ions on the planes of the 8-ring window sites. The most interesting feature is a CH3NH3+ ion at an opposite-4-ring site with no facing CH3NH3+ ions through a 6-ring window. On the other hand, in the MD simulation of NH4+-exchanged zeolite-A,4 the results have shown that there are two facing NH4+ ions through a 6-ring window with the disappearance of NH4+ ion at an opposite-4ring sites against an X-ray diffraction experimental result.10 According to the results of our MD simulations, as the size of the cation increases from Na+ to NH4+ and CH3NH3+, an opposite-4-ring site is occupied, then is not, and then is again occupied. However, this relationship depends not only on the size of the adsorbate but also on the number of ions in the β-cage and the heterogeneous ionic species. The four kinds of time correlation functions for CH3NH3+ show the dynamics of the ions which reflect the different structural arrangements of the ions well. The hydrogen bond time correlation functions of four types of CH3NH3+ ions decay more rapidly than those of NH4+ ions.4 From the analysis of the characteristic decay times for the H-bonds, about 0.8, 2.9, 0.6, 2.9, 1.6, 0.9, 1.7, and 2.0 hydrogen bonds formed between I and O(2), I and O(3), II and O(2), II and O(3), III and O(1), III and O(2), IV and O(1), and IV and O(3), are retained for 0.83, 1.89, 0.63, 1.61, 0.42, 0.31, 1.11, and 1.26 ps, respectively, before a breaking of the hydrogen bond occurs, leading to a significant exchange of O atoms hydrogen bonded to the ion. The rigid framework approximation adopted for simplicity in the present work is physically less adequate in the study of cation-exchanged zeolites. The inclusion of the intraframework interaction of zeolite-A is essential to take account of energy exchange between the adsorbed molecules and the framework atoms and with dynamical couplings of the adsorbate with framework vibrations, as well as the flexibility of the host lattice, and must be considered in the next study. Acknowledgment. This work was supported in part by a research grant to S.G.C. from Yangsan Junior College (1996)
and in part by a research grant (KOSEF 901-0303-032-1) to S.H.L. from the Korea Science and Engineering Foundation. The authors thank the Computer Center at Korea Institute of Science and Technology for the access to the Cray-C90 system. We gratefully acknowledge fruitful discussions with Prof. Karl Seff at the University of Hawaii. References and Notes (1) Moon, G. K.; Choi, S. G.; Kim, H. S.; Lee, S. H. Bull. Korean Chem. Soc. 1992, 13, 317. (2) Moon, G. K.; Choi, S. G.; Kim, H. S.; Lee, S. H. Bull. Korean Chem. Soc. 1993, 14, 356. (3) Lee, S. H.; Moon, G. K.; Choi, S. G.; Kim, H. S. J. Phys. Chem. 1994, 98, 1561. (4) Choi, S. G.; Lee, S. H. Mol. Simul. 1996, 17, 113. In Figures 2, 3, 4, and 5, NH4(2) and NH4(3) should be replaced by each other. (5) (a) de Leeuw, S. W.; Perram, J. W.; Smith, E. R. Proc. R. Soc. London 1980, A373, 27. (b) Anastasiou, N.; Fincham, D. Comput. Phys. Commun. 1982, 25, 159. (6) Pluth, J. J.; Smith, J. V. J. Am. Chem. Soc. 1980, 102, 4704. (7) Pluth, J. J.; Smith, J. V. J. Am. Chem. Soc. 1983, 105, 1192. (8) Jang, S. B.; Han, Y. W.; Kim, D. S.; Kim, Y. Korean J. Crystallogr. 1990, 1, 76. (9) Gramlich, V.; Meier, W. M. Z. Kristallogr. 1971, 133, 134. (10) McCusker, L B.; Seff, K. J. Am. Chem. Soc. 1981, 103, 3441. (11) Vance, T. B.; Seff, K. J. Phys. Chem. 1975, 79, 2163. (12) Firor, R. L.; Seff, K. J. Am. Chem. Soc. 1977, 99, 6249. (13) Breck, D. W. Zeolite Molecular SieVes; Wiley: New York; 1974; p 538. (14) Jeong, M. S.; Kim, Y.; Seff, K. J. Phys. Chem. 1994, 98, 1878. (15) Breck, D. W. Zeolite Molecular SieVes; Wiley: New York; 1974; p 570. (16) Kim, Y.; Seff, K. J. Phys. Chem. 1978, 82, 921. (17) Johnson, C. K. ORTEP. Report ORNL-3794; Oak Ridge National Laboratory: Oak Ridge, TN, 1972. (18) Meng, E. C.; Cieplak, P.; Caldwell, J. W.; Kollman, P. A. J. Am. Chem. Soc. 1994, 116, 1201. (19) Sprik, M.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1985, 83, 5802. (20) Gauss, K. F. J. Reine. Angew. Math. 1829, IV, 232. (21) Gear, C. W. Numerical initial Value problems in ordinary differential equation; Prentice-Hall: Englewood Cliffs, NJ, 1971. (22) (a) Evans, D. J. Mol. Phys. 1977, 34, 317. (b) Evans, D. J.; Murad, S. Mol. Phys. 1977, 34, 327. (23) Czjzek, M.; Jobic, H.; Fitch, A. N.; Vogt, T. J. Phys. Chem. 1992, 96, 1535. (24) Lee, S. H.; Rossky, P. J. J. Chem. Phys. 1994, 100, 3334. (25) (a) Lee, S. H.; Rasaiah, J. C. J. Chem. Phys. 1994, 101, 6964. (b) Lee, S. H.; Rasaiah, J. C. J. Phys. Chem. 1996, 100, 1420.
