Molecular dynamics simulation study: The decryption of bi and tri aromatics behavior with NaX zeolite

Journal of Molecular Graphics and Modelling 69 (2016) 61–71

Contents lists available at ScienceDirect

Journal of Molecular Graphics and Modelling journal homepage: www.elsevier.com/locate/JMGM

Molecular dynamics simulation study: The decryption of bi and tri aromatics behavior with NaX zeolite M. Haghighi Asl a,c,∗ , F. Moosavi b , J. Sargolzaei a , Kh. Sharifi d a

Department of Chemical Engineering, Ferdowsi University of Mashhad, P.O. Box 9177948944, Mashhad, Iran Department of Chemistry, Ferdowsi University of Mashhad, P.O. Box 9177948974, Mashhad, Iran Department of Engineering, Senior Process Engineer, South Pars Gas Complex, P. O. Box 311/75391, Assalouyeh, Iran d Research Institute of Petroleum Industry (RIPI), Azadi Sport Complex, P.O. Box 14665-1998, Tehran, Iran b c

a r t i c l e

i n f o

Article history: Received 11 May 2016 Received in revised form 21 July 2016 Accepted 8 August 2016 Available online 9 August 2016 Keywords: Molecular dynamics simulation Zeolites NaX Naphthalene Anthracene Diffusion Isosteric heats of adsorption

a b s t r a c t Molecular dynamics simulations have been carried out to provide an atomic description of the behavior of naphthalene and anthracene in NaX zeolite. The force field parameters, which were selected in this process, were chosen carefully to examine dependence of the self-diffusion coefficient of sorbates over a wide range of loading, temperature, and pressure. The simulated adsorption isotherm and calculated adsorption energies at low concentration were found to be in a reasonable qualitative and quantitative agreement with the corresponding scarce experimental data which can evaluate the effectiveness of proposed calculation method and force field parameters. The simulations provided new insights into the simulated concentration dependence of the atomic behavior of bi and tri aromatics inside NaXmicropores. Collective effects of the mutual interactions of sorbates molecules competing for the most preferable sites at the supercages (SCs) were found to be key factors responsible for the observed behavior of the adsorption isotherms, heat of adsorption, activation energy, and self-diffusivity. All involved calculations were performed in time period 6 ns and repeated calculations have been done at least two times to confirm the results and use the average values, which made the results being reliable. © 2016 Elsevier Inc. All rights reserved.

1. Introduction Zeolites are crystalline alumino silicates with micropores of uniform size extensively used in adsorption applications [1–7]. These materials may possess variety of framework structures originated from distinct microporous structures. The size and shape of the zeolite micropores are tremendously important for adsorption considering the organic compounds as adsorbates have molecular sizes similar to the pore dimension [8]. Literature survey demonstrates that slight variations in the size of the pores as well as the surface chemical composition are able to impact drastically on the performance of zeolitic materials as adsorbents and Molecular Sieves (MS) [9–20]. Faujasite (FAU) zeolites have extensive applications in adsorption, separation, and catalysis especially in the oil refining and petrochemical industries, mainly because of their high surface area, high hydrothermal stability, and appropriate acidity. In-depth

∗ Corresponding author at: Department of Chemical Engineering, Ferdowsi University of Mashhad, P.O. Box 9177948944, Mashhad, Iran. E-mail address: [email protected] (M.H. Asl). http://dx.doi.org/10.1016/j.jmgm.2016.08.002 1093-3263/© 2016 Elsevier Inc. All rights reserved.

knowledge of adsorption of sorbates inside FAU zeolites, in all series of applications, is prerequisite to clarify the underlying mechanism of these processes as well as exploring further applications of the zeolites [21]. MS 13X (faujasite framework-structure code FAU) is sodium form of the type X crystal with larger opening pores than a type A crystal. Molecules with a kinetic diameter less than 0.74 nm are adsorbed by MS 13X and those, which are larger, are excluded. In addition, MS 13X has the highest theoretical capacity amongst common adsorbents and possesses high standard mass transfer rates. Molecular sieve type 13X can be regenerated by either heating in the case of thermal swing processes or by lowering the pressure in the case of pressure swing processes [22]. Understanding the diffusion mechanism of guest molecule in these porous materials shed lights on enhancing the efficiency of industrial processes. In addition, it can also help designing transport optimized catalysts and membranes. Self-diffusivity (DS) is measured under equilibrium conditions by microscopic techniques, such as pulsed field gradient (PFG), NMR, and quasi-elastic neutron scattering (QENS) [23,24]. It has been shown that the diffusivities depend on the relevant space and time scales of observation as well as on the physical state under which the measurements are carried

62

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

out [25]. From the other side of view, thermal effect, pore size distribution, and surface roughness are other important factors that must be taken into account in diffusion measurements. However, real zeolite samples always contain imperfections, whereas “perfect” crystals without defects are used in molecular simulations [17,26]. Given the difficulties in experimentally obtaining accurate and reproducible diffusion coefficients, molecular simulation methods have been shown to be particularly powerful to study diffusion in zeolites [10,15,27–29]. Current available atomic models, force fields, and algorithms allow accurate calculations for diffusion coefficients in different systems. Molecular dynamics (MD) is the most common method to carry out simulations for diffusional behavior in porous adsorbents [30]. It is, therefore, interesting to investigate the adsorption and diffusion of aromatic hydrocarbons within the pores of such aluminosilicate zeolites. As described earlier, the diffusion of molecules can be studied experimentally; however, details of the diffusion process are very difficult to be achieved by experiments. Atomistic MD simulations based on analytical potential functions provide an appropriate way to study the microscopic details of such processes as well as being effective in modeling the diffusion of organic molecules in zeolites [10,11,31,32]. Considering these various important points on adsorption, the present study focuses on the diffusion of naphthalene and anthracene in zeolite 13X. Though there exists numerous simulations and experimental studies on diffusion of single aromatics in zeolites including FAU [15,21,29,33], the current research investigates the adsorption of polyaromatics on NaX zeolites. As the authors are aware, this study is performed for the first time. In contrast to the growing body of data related to adsorption and diffusion of aromatics in zeolites at dilute loading, such as benzene, toluene, xylene, and 1,3,5-trimethylbenzene (TMB), the research on the adsorption mechanism of larger aromatics at different loadings in zeolites requires more attention. Limited information is available on the entire loading range in the case of the same materials. The complexity of the system may significantly increase if pore filling approaches saturation values [33,34]. For example, anomalous adsorption behavior of cyclic hydrocarbons (benzene, p-xylene, and so on) has been reported in silicalite-1 at high sorbate loadings [15,29,33]. This work further explores the adsorption behavior of bi and tri aromatics (naphthalene and anthracene) on aluminasilicate FAU zeolite by force-field-based MD simulations, with a focus on evaluating the size-related effects of sorbates, if any, on the loading dependence of the adsorption behavior. The naphthalene/anthracene trajectories are less restricted at higher temperatures, which suggest displacements are limited to the length scale of one SC, the diameter of which is approximately 12 Å. We also investigate for the first time the effect of temperature and loading on self-diffusivities determined using MD for bi and tri aromatic hydrocarbon molecules adsorbed on X type zeolite. All results are for diffusion over the temperature range of 298–573 K and the loading range of 2–12 molecules per unit cell (molecules/UC). Self-diffusivities decrease markedly with loading for both sorbates at all temperatures examined. The calculated self-diffusivity of naphthalene molecules is in the order of 10−7 and anthracene 10−8 -10−7 cm2 /s. To the best of our knowledge, a detailed microscopic interpretation of the concentration dependence of the self-diffusivity and adsorption of naphthalene and anthracene in NaX zeolite has never been unraveled and remains to be established. In this work, we have also applied MD simulation qualitatively to obtain microscopic interpretation of the scarce experimentally observed concentration dependence of the bi and tri aromatics diffusion in NaX zeolite.

Fig. 1. Snapshots of FAU zeolite and extra framework cation distribution.

The activation energies of diffusion over a specified range of loadings and temperature are also reported. Clarifying whether an adsorption mechanism exists that is universally applicable for FAU zeolites can shed light on the separation processes, which are of interest to industries that generally use different types of aromatics. We continue this paper with adsorption isotherms compared quantitatively at low sorbate loadings with experimental data and also the loading dependence of the energetic property of adsorption. The isosteric heats of adsorption increases with increasing kinetic diameters of sorbates over the loading range 4–8 naphthalene molecules/UC and 4–16 anthracene molecules/UC and decreases at higher loadings. 1.1. Model and simulation details The structure of FAU was taken from X-ray and neutron powder diffraction data [35]. Each of two Supercages (SCs) was connected by a 12-membered window(12-T ring) with a diameter of 7.4 Å.The framework of FAU zeolite was shown in a neutron diffractionstudy which has a cubic unit cell with celllength a = 24.8536 Å. The accessible space for aromatic molecules wasthe SC with a diameter of ∼13 Å formed by sodalite cages, as shown in Fig. 1. Each unit cell contains 96 Si, 96 Al, 96 Na, and 384O atoms with Si/Al ratio equalt to one. As a result, NaX structure in a unit cell obeys the composition of Na96 Si96 Al96 O384 . As it is well known, the kinetic diameter is the smallest diameter thatallows the adsorbate molecule to enter the inner cavity (SC in this case).Thus, the largest sorbate employed in the present study wasnaphthalene, which had a kinetic diameter comparable with the size ofthe 12-T ring. All MD simulations were performed in NVT ensemble usingNosé-Hoover thermostat implemented in DL POLY programpackage version 2.17 [36]. At the beginning of the simulations, sorbate molecules were randomly distributed around NaX surface. The simulations were performed for the loadings inthe range between 1–8 molecules per unit cell and the highest loading of 36 molecules per unit cell was followed by the production run of 6 ns. Atomic trajectories were stored every 1000 fs for further analyses. The equations of motion were integrated by using the Verlet leapfrog algorithm with a time step of 1.0 fs. Periodic boundary conditions were applied to all three dimensions. The long-range electrostatic interactions were computed using Ewald summation method and short-range van der Waals interactions were evaluated up to a cutoff radius of 12.5 Å, the half of simulation box edge. The framework structure was built with strict alternation of the SiO4 and AlO4 tetrahedral and the extra framework Na+ ions were

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

63

Table 1 Types of interaction potentials and mathematical expression used in the analytical potential Vtotal . Type of terms

Expression

Used in

Harmonic stretching

Vtwo−body = 1/2kb (r − r0 )2

Vsorbate

Harmonic bending

Vthree−body = 1/2k  − 0

Quadratic bending

Vfour−body = A 1 + cos m − ı

Vsorbate

Lennard-Jones

VLennard−Jones = 4ε

Vsorbate, VZeolite Vsorbate- sorbate, VZeolite- sorbate

Coulomb

VLennard−Jones =







2

Vsorbate



 12  6  ␴ ␴ rij

−

rij

qi qj

Vsorbate-sorbate, VZeolite- sorbate

rij

Fig. 2. Numbering of carbon and hydrogen atoms of naphthalene (left) and anthracene (right).

distributed to different location sites as established by the previous powder neutron diffraction study shown in Fig. 1 which was fully occupied by Na+ ions according to model described by [27] and [37]. The extra framework cations were allowed to move while all other framework atoms were kept fixed at the crystallographic positions during the simulations. The total potential energy of the zeolite–sorbate complex consists of three terms in the potential function: Vtotal = Vsorbate + Vsorbate sorbate + VZeolite sorbate where Vsorbate is the intramolecular potential energy of naphthalene or anthracene, Vsorbate-sorbate is the potential energy between sorbate molecules and VZeolite-sorbate is the potential energy between the naphthalene/anthracene molecules and the zeolitic atoms. The mathematical expressions used in the potential functions are shown in Table 1. As Table 1 illustrates, the intramolecular potential energy of the sorbatemolecule is taken from the work of Rungsirisakun et al. [10]. Itconsists of five terms for bond stretching, bending, and torsionas well as van der Waals (vdw) and Coulomb terms: Vsorbate = Vbond + Vbend + Vtorsion + Vvdw + VCoul Bond stretching and bending terms are modeled byharmonic potentials. The torsion potential is represented bya harmonic cosine function. The van der Waals interactions aredescribed by a Lennard–Jones potential. The sorbate–sorbate and zeolite–sorbate interactions aredescribed by Lennard–Jones plus Coulomb terms: VNaphthalene−Naphthalene = VLennard-Jones + VCoul VZeolite−Naphthalene = VLennard-Jones + VCoul The parameters of the sorbate–sorbate and zeolite–sorbate interactions as well as the partial charges are reported in Table 2 described by [10,13,38]. The partial charges of sorbates were calculated according to DFT method at B3LYP/6–311 + +G(d,p) level of theory; after optimizing each electronic structure, by CHELPG method. The indication numbers of carbon and hydrogen atoms of naphthalene and anthracene molecules in order to use for partial charge calculation has been shown in Fig. 2. The partial charges of sorbates obtained by electrostatic surface potential (ESP) method with

Table 2 Potential parameters used for sorbate molecules with NaX (see equations from Table 1) and partial charges used for model zeolite. Two-body parameters C C C H Three-body parameters C C C C C H Four-body parameters C C C C C C C H Lennard–Jones parameters O Na Na-Na C C C H H H O C O H C Na H Na charge q [e]

A (kcal/mol) 3.6250 3.6250

O −1.2

kb (kcal/mol)

r0 (Å)

938.00 680.00 k␪ (kcal/mol)

1.37249 1.08425 0 (◦ )

126.00 70.00 m

120.289 120.289 ı

2.000 2.000 ␧ (kcal/mol)

180.0 180.0 ˚ ␴(A)

0.0725 0.3920 0.0860 0.0359 0.0150 0.0762 0.0318 0.1836 0.0767 Al Si +1.4+2.4

3.400 2.116 3.400 3.367 2.918 3.170 2.770 2.758 2.517 Na +1.0

Table 3 Calculated partial charges according to B3LYP/6–311 + +G (d, p) for naphthalene atoms. Type-number of atom C-1 C-2 C-3 C-4 C-5 C-6 H-7 H-8 H-9

Partial charge -0.04712 -0.26905 0.21651 0.21651 -0.26905 -0.04712 0.12931 0.07861 0.12931

Type-number of atom C-10 C-11 H-12 H-13 C-14 C-15 H-16 H-17 H-18

Partial charge -0.26905 -0.26905 0.12932 0.07860 -0.04712 -0.04712 0.12931 0.07860 0.07860

6–311 + +G (d,p) basis set werelisted in Table 3 for naphthalene and in Table 4 for anthracene.

64

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

Table 4 Calculated partial charges according to B3LYP/6–311 + +G (d, p) level of theory for anthracene atoms. Type-number of atom C-1 C-2 C-3 C-4 C-5 C-6 C-7 C-8 C-9 C-10 C-11 H-12

Partial charge -0.2180 -0.1590 -0.0690 -0.0630 -0.1820 -0.1850 -0.1380 -0.1380 -0.0690 -0.0630 -0.1820 0.2020

Type-number of atom C-13 C-14 C-15 H-16 H-17 H-18 H-19 H-20 H-21 H-22 H-23 H-24

Partial charge -0.1850 -0.2180 -0.1590 0.2000 0.2050 0.2020 0.2020 0.2040 0.2000 0.2040 0.2050 0.2020

Fig. 4. Radial distribution functions for Cn Na at temperature and pressure of 298 K and 1 bar for different naphthalene (a) and anthracene (b) loadings.

Fig. 3. Representative snapshots of (top) naphthalene, (bottom) anthracene in FAU model.

2. Results and discussion Fig. 3 shows simulated snapshotof a typical configuration of naphthalene and anthracene molecules in oneSC of FAU, respectively. The sorbate molecules prefer to accommodate in parallel to the surface of SC. The planes of the aromatic rings are perpendicular with respect to FAU pore structure, frequently found in this zeolite at various loadings molecules/SC. For simplicity, some abbreviations are applied to indicate sorbate and zeolite atoms. Si, Al, O and Na refer to the zeolite atoms and Cn and Hn are abbreviations applied for naphthalene and anthracene atoms. Radial distribution functions (RDFs) analysis give us information about positions and structural characteristics of sorbates in the pore of zeolites. Fig. 4(a,b) shows RDFs of Cn Na at the temperature and pressure of 298 K and 1 bar for different naphthalene and anthracene loadings, respectively. Though the height of the peak maximum in RDFs shows no special order with loading, the

maximum peak for anthracene is higher and sharper than naphthalene. It is obviously observed that the probability of vicinity of Cn Na for greater molecule (anthracene) is further than smaller one (naphthalene). At the condition of 2–16 molecules/UC, as shown in Fig. 4, the RDF of Cn Na for both sorbates is similar with the peak maximum at 2.7 Å, which confirms naphthalene and anthracene are adsorbed ideally at the overall range of loading. One can see that the position of the first maximum of RDFs does not change whether there are 2–16 molecules of sorbate present. For all loadings of FAU, the first peak position does not experience any change and has its first maximum at about 2.7 Å. A graphical analysis of these configurations reveals a range of flat-shaped structures. Fig. 3 is a snapshot of typical configuration of sorbate molecules in one SC of FAU. The sorbate molecules prefer parallel position with respect to the surface of the sodalite cage above the six-membered-ring. A general point of view is that the height of the peak maximum in RDF decreases with increasing naphthalene loading and vice versa is observed for anthracene. Actually, there is an increasing disordering between the heights of RDF’s intensities with different sorbate loadings. It may be related to the rearrangement of adsorbed molecules in a parallel configuration as a function of loading and steric hindrance. However, as the loading increases, the sorbate molecules diffuse into straight channel to be stacked parallel to each other (Fig. 3). This configuration, including some roughly “parallel stacked and displaced” structure is also observed at high loadings of 16 molecules/UC where the movement of the molecules is also highly correlated. This behavior may also be responsible for the distorted height of the peaks. Furthermore, itcan be seen that the intensity of the first peak of g(r) for anthracene is generally greater than naphthalene. This may be consistent with the more extended structure of anthracene rather than naphthalene and more probability of Na-Cn vicinity. Fig. 5(a,b) shows the radial distribution functions of Cn with other zeolite atoms at 298 K and 1 bar for 6 molecules/UC. The position of the first maximum of Cn Na pair is located at the distance

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

Fig. 5. Radial distribution functions of Cn with zeolite atoms at 298 K and 1 bar for the 6 naphthalene/UC (a) and 6 anthracene/UC (b).

of 2.725 Å and in contrast, the first peak of Cn with Al, Si, and O is located almost over 4 Å. An analysis of the first maximum of the cations surrounded by adsorbed molecules in comparison with other zeolite atoms reveals that the cations are the most preferable sites for naphthalene and anthracene adsorption which is completely in qualitative agreement with the adsorption position of the same sorbates in NaY zeolite [39]. Considering the influence of temperature, inthe case of pressure 1 bar, on the adsorption process as shown in Fig. 6(a,b), the distances amongCn Na is about 2.725 Å at different temperatures (298, 373, 473, and 573 K) though the RDF intensity decreased gradually while temperature increased. This result is in agreement with the kinetic energy of the molecules that increases with increasing temperature and mean free path of Cn and Na atoms that enhances. In fact, some decrease of Cn congestion around Na ions is observed. It might be interesting that RDF peaks of naphthalene have been completely broadened at higher temperatures but the corresponded anthracene peaks keep almost constant in sharpness with a slight decrease in height. This behavior may be corresponded to the steric hindrance against freely motion of anthracene molecules according to their length scale and higher adsorption tendency. Fig. 7(a,b) shows the RDF of the O Hn for different naphthalene and anthracene loadings plotted at temperature and pressure of 298 K and 1 bar. It can be seen that the intensity of the first peak of g(r) gradually decreases as loading increases, which indicates the vicinity of oxygen atoms becomes increasingly less occupied by Hn atoms. As Fig. 7 demonstrates, the main peaks stay almost at the same position over the entire loading range. As discussed earlier at high loading, the more aggregation of Cn Na pairs occurs; consequently, the segregation of O Hn may be observed by considering the steric hindrance. Mean square displacements (MSD) and the velocity auto correlation function (VACF) obtained from the simulations can be used

Fig. 6. The effect of temperature on Cn lene/UC (a) and 6 anthracene/UC (b).

65

Na distribution at 1 bar and 6 naphtha-

Fig. 7. RDF of the O Hn for different naphthalene (a) and anthracene (b) loadings plotted at the temperature and pressure of 298 K and 1 bar.

in analyzing the self-diffusion coefficients following the Einstein and Green-Kubo equations [40].

2 ∂MSD 1 1   rj (t) − rj (0)  lim = 6 t→∞ ∂t 6Nm Nm

D=

j=1

(1)

66

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

Fig. 8. MSD of the Cn atoms of the naphthalene (a and c) and anthracene (b and d) molecules for different loadings (2–8 molecules per UC) plotted as a function of time and various temperatures.

1 D= 3Nm

∞  Nm 

vj (t) .vj (0) dt

(2)

j=1

0

wherer(t) and  (t) are three-dimensional sorbate position and velocity at time t, respectively, and Nm is the component loading. Einstein and Green–Kubo relationships are equal in theory. The term

Nm  j=1

rj (t) − rj (0)

Nm

2

in Einstein relation is regarded as MSD,

v .v in Green–Kubo relation is usually while the term j=1 j (t) j (0) called VACF. Fig. 8(a,d) reports the MSDs for naphthalene and anthracene in NaX at the various loadings investigated (2, 4, 6 and 8 molecules/UC) and for different temperatures at the range 298–573 K. The superior statistics obtained for self-diffusion exhibits non linearity of long-time MSDs. Furthermore, especially for the lowest loading (2 molecules/UC), the MSDs appear to be reasonably nonlinear over a broad time domain indicating restricted three dimensional diffusion within the zeolite pore system. The values of MSD are greater than 20 Å2 over 500 ps for naphthalene and up to 16 Å2 for anthracene at the same time of simulation, which denote hindered diffusion across just one SC. In somecases, it can also be seen that diffusivity does not always increase smoothly with temperature as, for instance, is evident from the noticeable “gap” between the plots for 298 and 373 K in Fig. 8d. Similar behavior can be discerned in the MSD plots for 2 and 4 naphthalene molecules/UC (Fig. 8b). For higher loadings (6 and 8 molecules/UC), the MSD plots shown in Fig. 8(a,b) appear more linear. In addition, the naphthalene/anthracene trajectories are less restricted (MSDs reach more than 100 Å2 for naphthalene and over 20 Å2 for anthracene) at high temperatures, which suggest displacements are again limited to the length scale of one SC, the diameter of which is approximately 12 Å. Similarly, at lower loadings and temperatures, the average displacement of the molecules is less than the dimensions of a SC. In target FAU zeolite, in spite that small molecules have freedom to move in all directions from one SC to another cage (intercage), this is not performed for larger molecules such as naphthalene and anthracene. Due to smaller size, naphthalene molecules (8

molecules/UC) have a better ability to move inside one SC than anthracene as shown in Fig. 9. Movement tendency of naphthalene molecules towards the center of phase space shows a better diffusion of naphthalene molecules within the zeolite pores rather than anthracene molecules. In fact, as revealed qualitatively in Fig. 9, almost all anthracene molecules were adsorbed on the external surface of the zeolite. This diffusion tendency which is different from self-diffusivity leads a very limited diffusion of larger molecules into the NaXmicropores obtained by [41] in the range of 10−8 and 10−15 cm2 /s for naphthalene and anthracene, respectively. As interpreted qualitatively in Fig. 9, anthracene molecules could not diffuse toward the center NaX zeolite atoms even at higher temperatures. At high loadings such as 8 molecules/UC and lower temperatures, reduced translational mobility due to a reduced effective volume in observed for the individual molecule. The sorbate–sorbate interaction further contributes to these movements. Similar results were also found for alkane diffusion in FAU zeolite [11]. Self-diffusion of naphthalene and anthracene molecules in zeolite 13X was calculated at 298, 373, 473, and 573 K, and at loadings of 2, 4, 6, 8, and 12 molecules/UC. A direct comparison between results for those two sorbates is shown in Fig. 10. As expected, the diffusivity increases with temperature in all cases. The loading dependency of naphthalene and anthracene on self-diffusivities in zeolite 13X from which it is clear that self-diffusivities decrease markedly with loading for both sorbates at all temperatures examined. Several studies on self-diffusion of light alkanes and single aromatics such as benzene in zeolites by Monte Carlo or MD exist in the literature [10,11,15–17,26,27,31,42]. To the best of our knowledge, no experimental results for naphthalene/anthracene self-diffusion in FAU have been reported. As a result, it can be found that the present results open new windows for considering polyaromatics adsorption and diffusion on NaX. We compared the estimated temperature dependence of our work with benzene’s self-diffusivity in NaX. As expected for a guest-zeolite system where guest motion is a size related phenomenon, the self-diffusivity trends obeys the size dependence

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

67

Fig. 9. Trajectory of the pathway of naphthalene (a & c) and anthracene (b & d) molecules during the simulation at the condition of 1 bar and 298 K (top) and 1 bar and 573 K (bottom) for 8 molecules/UC.

behavior. The self-diffusivity of benzene molecules are in the range of 10−5 -10−4 cm2 /s [27]; however, the self-diffusivity of naphthalene molecules is in the order of 10−7 and anthracene 10−8 10−7 cm2 /s. The data for 2, 4, 8, and 12 naphthalene and anthracene molecules/UC are shown as Arrhenius plots in Figs. 11 and 12. The activation energies, which are characteristic of long-range motion, were derived from linear least squares fits to the Arrhenius plots. It can be observed that these activation energies follow an overall increasing trend as the loading increases, even though with an unexpected decrease between loadings of 2 and 4 naphthalene molecules/UC. This trend can be interpreted as follows: during the initial stage of naphthalene loading, the value of the activation energy is mainly governed by the strong interaction between the sorbate molecules and the Na+ extra framework cations. As the loading increases, up to 4 naphthalene molecules/UC, this dominant interaction, which tends to maintain the naphthalene molecules in the vicinity of the Na+ cations, is more and more attenuated by the increasing number of surrounding naphthalene molecules, thereby leading to a decrease of the activation energy to about 6.53 kJ/mol at 4 naphthalene molecules/UC. Of course, it should be recalled that the diffusivities generally rise up to 12 molecules/UC and that cooperative interactions between naphthalene molecules also play an important role. The properties of adsorption isotherms are important for understanding the adsorption behavior of sorbates in micropores. Fig. 13

shows the adsorption isotherms developed for the first time for naphthalene on the NaX Zeolite obtained theoretically at the loading of 8 molecules/UC and different temperatures. Adsorption isotherms are simulated for pure naphthalene on NaX at 298, 373, 473, and 573 K and pressure range up to 300 bar. The results presented in Fig. 13 show the effect of the temperature on the adsorption capacity of naphthalene on NaX. The loadings of naphthalene decrease with increasing temperature because of thermal motion enhancing with temperature. The adsorption isotherms obtained theoretically are of type I, as shown in Fig. 13, with a maximum loading of 0.63, 0.56, 0.50 and 0.36 mmol/g at 298, 373, 473, and 573 K, respectively. However, absolutely scarce experimental data available for polyaromatic adsorption isotherms on NaX limits comparison. Interestingly, one can observe simulated isotherms by MD are in good agreement with those obtained by experiment in terms of the shape and the amount of adsorbed sorbates at 298 K and low concentrations. This allows validation of microscopic model for both host and guest as well as force field parameters used for describing the system throughout this work.The sorbates exhibit no inflection point, in which the curvature changes from positive (increasing slope) to negative (decreasing slope). This observation is similar to experimental results reported by [43] for polyaromatics in FAU system and might be related to the fact that there is no change of adsorption mechanism with the pressure (loading) of naphthalene at low concentrations. In order to further explore the effect

68

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

Fig. 10. Self-diffusion coefficients of naphthalene (top) and anthracene (bottom) in NaX from MD simulations.

of size-related packing on adsorption mechanisms, the spatial distributions of aromatics on the NaX model need to be carefully examined. For comparison of adsorption mechanisms of naphthalene and anthracene in more real atmosphere, a simulation with higher concentrations (36 molecules/UC) has been carried out and the results of adsorption isotherms are illustrated in Fig. 14. This observation for naphthalene and anthracene is in reasonable order agreement with other molecules adsorption in FAU-type zeolites X [37,44–49]. Based on the adsorption mechanisms, the loading dependence of the energetic property of adsorption for the aromatics in the NaX model can be improved significantly. The isosteric heats of adsorption of aromatics (Q) are calculated as follows:



Q = RT −

∂ (Uad − Uintra ) ∂Nad



(3)

by the zeolite cations attains its maximum value. The sharp drop corresponds to a saturation of the energy minima and generation of configurations at less energetically favorable positions. Since for lower loadings, electrostatic interactions are not in minimal conditions, the zeolite cavity serves as an almost inhomogeneous surface and the isosteric heat shows stepwise increase. The average occupancy of the zeolite cavity varies from 8 to 36 molecules for naphthalene and from 16 to 36 molecules for anthracene and the adsorbed phase exhibits moderate and sharp negative deviations for naphthalene and anthracene, respectively. This adsorption heats were then found from the temperature dependence of the Henry constant described by [50] calculated 25.5 kcal/mol and 34.0 kcal/mol for naphthalene and anthracene, respectively which illustrate good agreement at low concentration with our simulation results.

V,T

where R is the gas constant, Nad is the loading of sorbate, Uad is determined by summing all interactions between pairs of sorbate molecules (Uads-ads ) as well as all sorbate interactions with the zeolite framework (Uads-zeo ), and Uintra is the intramolecular energy of the sorbate molecules. Thus, a positive correlation is observed between Q and the absolute value of Uad . Fig. 15 shows Q of naphthalene and anthracene on NaX at 298 K. Two dominant features of Q are revealed. First, Q increases with increasing kinetic diameters of sorbates over the loading range 4–8 naphthalene molecules/UC and 4–16 anthracene molecules/UC and decreases at higher loading. Second, one changing point separates the data of Q into two ranges (low and high loading) which is distinguished by blue and red dashed lines at about 8 naphthalene molecules and 16 anthracene molecules/UC. It is interesting that the isosteric heats of naphthalene and anthracene show a step upward at coverage of 8 and 16 molecules/UC, respectively, and then downward at coverage of higher loadings. For adsorption of molecules, the model used in simulations has some major sites corresponding to the positions where the electric field generated

3. Conclusions This study reports the adsorption of di- and tri-aromatics in FAU zeolite, which is the first case concerning the basic law of adsorption mechanism at various loadings, temperatures, and pressures. The change of adsorption mechanism due to rearrangement is observed for naphthalene and anthracene, which is believed to be a general characteristic of adsorption in FAU zeolites for polyaromatics. Besides, self-diffusion of naphthalene and anthracene in the NaX zeolite was investigated by using MD technique. The calculated self-diffusion coefficients at various temperatures and loadings are almost in the range of 10−7 and 10−8 -10−7 cm2 s−1 for naphthalene and anthracene, respectively, which are in the range of experimental values. Tracking the trajectories reveals that naphthalene molecules (8 molecules/UC) have a better ability to move inside one SC than anthracene with the same loading. Strongest adsorption sites in NaX (extra framework cations) are well separated, so stabilizing interactions between neighboring sorbate molecules are relatively strong. This leads to the observed

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

69

Fig. 11. Arrhenius plot of estimated loading-dependent activation energies for long-range motion with 2 (a), 4 (b), 8 (c) and 12 (d) naphthalene molecules/UC.

Fig. 12. Arrhenius plot of estimated loading-dependent activation energies for long-range motion with 2 (a), 4 (b), 8 (c) and 12 (d) anthracene molecules/UC.

increase of the isosteric heat of adsorption on increasing sorbate loading. It can be observed that the activation energies follow an increasing trend from 5 to 28 kJ/mol as loading increases from 2 to 8 for anthracene, even though the mentioned property experiences an unexpected decrease for naphthalene molecules/UC by increasing loadings of 2 and 4 which can be interpreted because of mobility hindrance of naphthalene molecules.

Moreover, adsorption isotherms of the naphthalene and anthracene in NaX zeolite were considered. From the obtained data, finite adsorption capacities and active site occupation probabilities are demonstrated. The simulated adsorption isotherm and calculated adsorption energies at low concentration were found to be in a qualitative and quantitative agreement with the corresponding scarce experimental data evaluating the effectiveness of proposed computation method and force field parameters.

70

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71

[2]

[3]

[4]

[5] [6]

[7]

[8] Fig. 13. Adsorption isotherms developed for naphthalene on the NaX Zeolite obtained theoretically at the loading of 8 molecules/UC and different temperatures.

[9]

[10]

[11]

[12]

[13]

[14]

[15]

[16] Fig. 14. Adsorption isotherms developed for high concentration 36 (naphthaleneanthracene)/UC on the NaX Zeolite obtained at 298 K and 1 bar according to amount adsorbed.

[17]

[18]

[19]

[20]

[21]

[22] [23]

[24]

Fig. 15. Calculated Q of sorbates at room temperature in the NaX zeolite at different loadings.

[25] [26]

Acknowledgment We gratefully acknowledge support from the High Performance Computing Centre of Ferdowsi University of Mashhad.

[27]

[28]

References [29] [1] P.A. Russo, M.M.L.R. Carrott, P.J.M. Carrott, P. Matias, J.M. Lopes, M. Guisnet, F.R. Ribeiro, Characterisation by adsorption of various organic vapours of the

porosity of fresh and coked H-MCM-22 zeolites, Micropor. Mesopor. Mat. 118 (2009) 473–479. Y. Nakasaka, T. Tago, H. Konno, A. Okabe, T. Masud, Kinetic study for burning regeneration of coked MFI-type zeolite and numerical modeling for regeneration process in a fixed-bed reactor, Chem. Eng. J. 207–208 (2012) 368–376. H. Yi, H. Deng, X. Tang, Q. Yu, X. Zhou, H. Liu, Adsorption equilibrium and kinetics for SO2 NO, CO2 on zeolites FAU and LTA, J. Hazard. Mater. 203-204 (2012) 111–117. H. Deng, H. Yi, X. Tang, Q. Yu, P. Ning, L. Yang, Adsorption equilibrium for sulfur dioxide nitric oxide, carbon dioxide, nitrogen on 13X and 5A zeolites, Chem. Eng. J. 188 (2012) 77–85. G. Weber, J.P. Bellat, F. Benoit, C. Paulin, Adsorption equilibrium of light mercaptans on faujasites, Adsorption 11 (2005) 183–188. M. Granato, M. Jorge, T.J.H. Vlugt, A.E. Rodrigues, Diffusion of propane, propylene and isobutane in 13X zeolite by molecular dynamics, Chem. Eng. Sci. 65 (2010) 2656–2663. G. Weber, F. Benoit, J.-P. Bellat, C. Paulin, P. Mougin, M. Thomas, Selective adsorption of ethyl mercaptan on NaX zeolite, Microporous Mesoporous Mater. 109 (2008) 184–192. J. Cejka, H.v. Bekkum, A. Corma, F. Schüth, Introduction to Zeolite Science and Practice, Elsevier, Amsterdam, 2007. M. Trejda, M. Ziolek, P. Decyk, D. Duczmal, The radical species and impurities present in mesoporous silicas as oxidation active centres, Micropor. Mesopor. Mat. 120 (2009) 214–220. R. Rungsirisakun, T. Nanok, M. Probst, J. Limtrakul, Adsorption and diffusion of benzene in the nanoporous catalysts FAU, ZSM-5 and MCM-22: A molecular dynamics study, J. Mol. Graphics Modell. 2 (2006) 373–382. P. Wongthong, R. Rungsirisakun, M. Probst, J. Limtrakul, Adsorption and diffusion of light alkanes on nanoporous faujasite catalysts investigated by molecular dynamics simulations, Microporous Mesoporous Mater. 100 (2007) 160–166. H. Yi, H. Deng, X. Tang, Q. Yu, X. Zhou, H. Liu, Adsorption equilibrium and kinetics for SO2 NO, CO2 on zeolites FAU and LTA, J. Hazard. Mater. 203-204 (2012) 111–117. A.M. Granato, N.L. Thijs, J.H. Vlugt, A.E. Rodrigues, Adsorption equilibrium of isobutane and 1-Butene in zeolite 13X by molecular simulation, Ind. Eng. Chem. Res. 47 (2008) 6166–6174. Z. Yongping, J. Shengui, X. Weihong, C. Changlin, Adsorption of mercaptan from model gasoline on 13X loaded with Zn2+, Can. J. Chem. Eng. 86 (2008) 186–191. T. Mosell, G. Schrimpf, J. Brickmann, Diffusion of aromatic molecules in zeolite NaY. 1. constrained reaction coordinate dynamics, J. Phys. Chem. B 101 (1997) 9476–9484. D. Plant, H. Jobic, P. Llewellyn, G. Maurin, Diffusion of CO2 in NaY and NaX Faujasite systems: quasi-elastic neutron scattering experiments and molecular dynamics simulations, Eur. Phys. J. Special Topics 141 (2007) 127–132. M.A. Granato, M. Jorge, T.H. Vlugt, A.E. Rodrigues, Diffusion of propane, propylene and isobutane in 13X zeolite by molecular dynamics, Chem. Eng. Sci. 65 (2010) 2656–2663. C. Laborde-Boutet, G. Joly, A. Nicolaos, M. Thomas, P. Magnoux, Selectivity of Thiophene/Toluene competitive adsorptions onto NaY and NaX zeolites, Ind. Eng. Chem. Res. 45 (2006) 6758–6764. S. Buttefey, A. Boutin, C. Mellot-Draznieks, A.H. Fuchs, A simple model for predicting the Na+ distribution in anhydrous NaY and NaX zeolites, J. Phys. Chem. B 105 (2001) 9569–9575. J. Zhang, N. Burke, S. Zhang, K. Liu, M. Pervukhina, Thermodynamic analysis of molecular simulations of CO2 and CH4 adsorption in FAU zeolites, Chem. Eng. Sci. 113 (2014) 54–61. H. Zheng, L. Zhao, J. Ji, J. Gao, C. Xu, F. Luck, Unraveling the adsorption mechanism of mono- and diaromatics in faujasite zeolite, ACS Appl. Mater. Interfaces 7 (2015) 10190–10200. Z. Sarbak, Desulfurization of ethanethiol over cadmium and mercury modified zeolite NaX, Appl. Catal. A-Gen. 147 (1996) 47–54. U. Hong, J. Karger, B. Hunger, N.N. Feoktistova, S.P. Zhdanov, In situ measurement of molecular diffusion during catalytic reaction by pulsed field gradient spectroscopy, J.Catal. 137 (1992) 243–251. H. Jobic, M. Bée, J. Caro, M. Bülow, J. Kärger, Molecular self-diffusion of methane in zeolite ZSM-5 by quasi-elastic neutron scattering and nuclear magnetic resonance pulsed field gradient technique, J. Chem. Soc. Faraday Trans. 1 85 (1989) 4201–4209. J. Karger, Measurement of diffusion in Zeolites—A never ending challenge? Adsorption 9 (2003) 29–35. D.F. Plant, G. Maurin, R.G. Bell, Diffusion of methanol in zeolite NaY: a molecular dynamics study, J. Phys. Chem. B 111 (2007) 2836–2844. H. Ramanan, S.M. Auerbach, M. Tsapatsis, Beyond lattice models of activated transport in zeolites: high-Temperature molecular dynamics of self-diffusion and cooperative diffusion of benzene in NaX, J. Phys. Chem. B 108 (2004) 17171–17178. V. Lachet, A. Boutin, B. Tavitian, A.H. Fuchs, Computational study of p-xylene/m-xylene mixtures adsorbed in NaY zeolite, J. Phys. Chem. B 102 (1998) 9224–9233. H. Klein, H. Fuess, Mobility of aromatic molecules in zeolite NaY by molecular dynamics simulation, J. Phys. Chem. B 100 (1996) 11101–11112.

M.H. Asl et al. / Journal of Molecular Graphics and Modelling 69 (2016) 61–71 [30] D.C. Rapaport, The Art of Molecular Dynamics Simulation, Cambridge University Press, 2004. [31] P. Demontis, E.S. Fois, G.B. Suffritti, S. Quartieri, Molecular dynamics studies on zeolites. 4. Diffusion of methane in silicalite, J. Phys. Chem. B 94 (1990) 4329–4334. [32] R.L. June, A.T. Bell, D.N. Theodorou, Molecular dynamics studies of butane and hexane in silicalite, J. Phys. Chem. 96 (1992) 1051–1060. [33] I. Daems, P. Leflaive, A. Me´ıthivier, J.F.M. Denayer, G.V. Baron, Size and packing related adsorption effects in the liquid phase adsorption of aromatics and alkenes on FAU type zeolites, Stud. Surf. Sci. Catal. 158 (2005) 1177–1184. [34] D. Keffer, H.T. Davis, A. McCormick, The effect of nanopore shape on the structure and isotherms of adsorbed fluids, Adsorption 2 (1996) 9–21. [35] F. Porcher, M. Souhassou, Y. Dusausoy, C. Lecomte, The crystal structure of a low-silica dehydrated NaX zeolite, Eur. J. Mineral. 11 (1999) 333–343. [36] W. Smith, T.R. Forester, DL POLY 2.0. A general-purpose parallel molecular dynamics simulation package, J. Mo. Graphics 14 (1996) 136–141. [37] P. Cosoli, M. Ferrone, S. Pricl, M. Fermeglia, Hydrogen sulfide removal from biogas by zeolite adsorption. Part II. MD simulations, Chem. Eng. J. 145 (2008) 93–99. [38] W.D. Cornell, P. Cieplak, C.I. Bayly, I.R. Gould, K.M. Merz, D.M. Ferguson, D.C. Spellmeyer, Thomas Fox, J.W. Caldwell, P.A. Kollman, A second generation force field for the simulation of proteins, nucleic acids, and organic molecules, J. Am. Chem. Soc. 117 (1995) 5179–5197. [39] C. Baehtz, H. Fuess, Naphthalene anthracene and 2,3-benzanthracene in zeolite Y, Phys. Chem. Chem. Phys. 4 (2002) 4543–4548. [40] L.W. Zhong, C. Cong, Y. Jian, Molecular dynamics simulation of self-diffusion coefficient and its relation with temperature using simple Lennard-Jones potential, Heat Trans.—Asian Res. 37 (2008) 86–93. [41] D.M. Ruthven, B.K. Kaul, Adsorption of aromatic hydrocarbons in NaX zeolite. 2. Kinetics, Ind. Eng. Chem. Res. 32 (1993) 2053–2057.

71

[42] E. Beerdsen, D. Dubbeldam, B. Smit, Loading dependence of the diffusion coefficient of methane in nanoporous materials, J. Phys. Chem. B 110 (2006) 22754–22772. [43] M.I. Falkovieh, A.A. Gureev, F. Aveisi, Y.u.A. Tltekov, V.T. Solodovnikova, Kinetics of Adsorption of Aromativ Hydrocarbons and Organic Sufur Compounds by Zeolite, UDC 665. 544: 541. 183.03, Plenum Publishing Corporation, 1982, pp. 13–15. [44] S.K. Walton, B.M. Abney, M.D. LeVan, CO2 adsorption in Y and X zeolites modified by alkali metal cation exchange, Microporous Mesoporous Mater. 91 (2006) 78–84. [45] C. Beauvais, A. Boutin, A.H. Fuchs, Adsorption of water in zeolite sodium-faujasite A molecular simulation study, C. R. Chim. 2005 (2005) 485–490. [46] M.A. Granato, T.J.H. Vlugt, A.E. Rodrigues, Molecular simulation of propane-Propylene binary adsorption equilibrium in zeolite 13X, Ind. Eng. Chem. Res. 46 (2007) 7239–7245. [47] S. Calero, D. Dubbeldam, R. Krishna, B. Smit, T.J.H. Vlugt, J.F.M. Denayer, J.A. Martens, T.L.M. Maesen, Understanding the role of sodium during adsorption: a force field for alkanes in sodium-Exchanged faujasites, J. Am. Chem. Soc. 124 (2004) 11377–11386. [48] L. Joos, J.A. Swisher, B. Smit, Molecular simulation study of the competitive adsorption of H2O and CO2 in zeolite 13X, Langmuir 29 (2013) 15936–15942. [49] A. Takahashi, F.H. Yang, R.T. Yang, New sorbents for desulfurization by Ð-Complexation: thiophene/Benzene adsorption, Ind. Eng. Chem. Res. 41 (2002) 2487–2496. [50] D.M. Ruthven, B.K. Kaul, Adsorption of aromatic hydrocarbons in NaX zeolite. 1. equilibrium, Ind. Eng. Chem. Res. 32 (1993) 2047–2052.