Quantum Effects on the Diffusivity of Hydrogen Isotopes in Zeolites

1 day ago - Our results, show at temperatures below 77K that H2 has a lower adsorption capacity and a lower diffusion coefficient than D2 and T2. Here...
0 downloads 0 Views 617KB Size
Subscriber access provided by Columbia University Libraries

C: Surfaces, Interfaces, Porous Materials, and Catalysis

Quantum Effects on the Diffusivity of Hydrogen Isotopes in Zeolites Jose Marcos Salazar, Michael Badawi, Bastien Radola, Matheiu Macaud, and Jean-Marc Simon J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.9b05090 • Publication Date (Web): 01 Sep 2019 Downloaded from pubs.acs.org on September 1, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Quantum Effects on the Diffusivity of Hydrogen Isotopes in Zeolites J. M. Salazar,∗,† M. Badawi,‡ B. Radola,† M. Macaud,¶ and J. M. Simon† †ICB-UMR 6303 CNRS Laboratoire Interdisciplinaire Carnot de Bourgogne Université de Bourgogne Franche-Comté 9, Av. A. Savary B.P. 47870 F-21078 Dijon Cedex, France. ‡ Laboratoire de Physique et Chimie Théoriques UMR CNRS 7019 Institut Jean Barriol FR2843 CNRS IUT de Moselle-Est, Département de Chimie Rue Victor Demange B.P. 80105, 57503 SAINT-AVOLD Cedex ¶ CEA, DAM, VALDUC, F-21120 Is-sur-Tille, France E-mail: [email protected]

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract The molecular sieving of H2 and its isotopes produced by nuclear plants is a long standing research where some adsorption processes are well identified. However, some mechanisms governing the diffusion of the proton and its isotopes within a faujasite NaX at low temperatures (40 -100)K still troublesome to characterize. Notwithstanding, an understanding of the processes governing adsorption diffusion within narrow pores is essential for the development of recycling procedures of gases produced by nuclear plants. At cryogenic temperatures it is well known that quantum effects are revealed and the heterogeneity of the guest structure plays an important role in the adsorption process. Here, we focus on the consequences of these two factors on molecular sieving and transport coefficients, based on molecular dynamics including the Feynman-Hibbs quantum approximation. Our results, show at temperatures below 77K that H2 has a lower adsorption capacity and a lower diffusion coefficient than D2 and T2. Here we give an original explanation of this diffusion inversion in terms of the activation energies. We show that this energy is greater for H2 by 30 % and 50% than for D2 and T2, respectively. Moreover, experimentally it has been shown that the pore heterogeneity of the faujasite NaX leads to an increase of the self diffusion coefficient with loading. For explaining this unexpected behavior, several authors have proposed either a scenario based on privileged adsorption sites or based on the residence time of molecules near the guest pore walls, followed by isotropic jumps. Here, we report an original analysis based on the mean square displacements revealing the presence of slow and fast molecular mobility regimes. The ensemble of our results provide useful physical information for the development of recycling procedures of gases produced by nuclear plants.

Introduction In a near future a source of energy may be provided by fusion reactors. 1 This new generation of reactors will produce and consume a considerable amount of gases composed of helium and H2 isotopes. These gases contain tritium compounds and for security reasons they must 2

ACS Paragon Plus Environment

Page 2 of 29

Page 3 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

be recycled in situ. In this context some technics available for sieving the proton and its isotopes are: diffusion through palladium membranes, 2 cryogenic distillation 3 or the chemical exchange method proposed by Girdler. 4 These costly technics have the disadvantage of having a small selectivity factor. Therefore a promising alternative is the use of sorbents, like: carbon nanotubes or MOF’s (metal organic frameworks). Unfortunately, these last two examples present the inconvenient formation of radioactive tritium compounds. 5 Nonetheless, the physicochemical stability of the well known zeolitic materials place them as an alternative for the physisorption of H2 and its isotopes at cryogenic temperatures. In these last materials, the space available for molecular motion of hydrogen isotopes is comparable with de Broglie wavelength, implying the presence of quantum effects. A theoretical treatment of these effects can be seen in the work of Beenakker et al. 6 where the authors proposed the concept of quantum molecular sieving originated by isotopic quantum fluctuations. Several experimental studies devoted to zeolites have shown that for temperatures below 100K D2 and T2 are preferentially adsorbed by the guest structure rather than H2 . The molecular dynamics (MD) studies of Kumar and Bhatia 7 have shown a reverse kinetic sieving effect where deuterium (D2 ) diffuses faster than hydrogen. At cryogenic temperatures classical arguments fail to explain this preferential adsorption. Notice that, at 40K the de Broglie wavelength for H2 (λ(H2 ) ≈ 2.18 Å) is greater than for D2 (λ(D2 ) ≈ 1.6 Å). In other words, H2 becomes more cumbersome than D2 as temperature decreases. When performing classical molecular dynamics simulations the quantum effects appearing at low temperatures can be represented by using the Feynman-Hibbs approach (FH). When the reduced mass and temperature are decreased this quantum correction is reflected on the Lennard Jones (LJ) potential as: an increase of the kinetic diameter and a reduction of the well depth. 8 A recent computational screening over 210 pure-silica zeolites reported an astonishing high selectivity of D2 against H2 (e.g., SD2 /H2 ∼ 105 for the BCT zeolite). 9 Clearly, these selectivities result from steric constraints 10 and are maximum for most narrow pores (e.g., BCT zeolite). However, the majority of zeolites industrially used (e.g., FAU, LTA, CHA and

3

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

MFI) have selectivities smaller than six. 11 For increasing the selectivity above this limiting value the role of cations needs to be studied systematically for an optimization of the zeolitic structures and the sieving process. In previous experimental studies for determining the transport coefficients of H2 in the faujasite NaX, by quasielastic neutron scattering (QENS) 7 and also by simulations, 12 have shown an increase of the self-diffusion coefficient of H2 with loading (the maximal capacity used by the authors was of 50%). For developing an explanation to this bizarre behavior some hypothesis have been formulated concerning the role of cations and the influence of the inhomogeneity of the guest pore in the adsorption process. The ideas developed for the diffusion of H2 from Kahn et al., 13 within a NaA zeolite (70-150K) and from Papadopoulos et al., 14 for a NaX zeolite are respectively: i) the presence of two populations of molecules: slow and fast, with the former undergoing isotropic jumps of mean length of 3.9 Å and a residence time between jumps of τ0 ≈ 10.8 ps. ii) the later authors, suggested the adsorption on discrete sites. In any case, there is a consensus to say that the inhomogeneity of the guest cage plays a non- neglieable role on the adsorption process. By looking at the faujasite NaX (Na86 Al86 Si106 O384 ) with cations disposed on sites II and III, we can realize with ease that the atomistic environment of the cations is not similar. The aim of the studies presented here is to enlarge the known scenario of the diffusivity of H2 and its isotopes, going from low to high loading on the faujasite NaX and from [40100]K. We will focus particularly on explaining the change of H2 and D2 diffusion behavior observed when decreasing the temperature, 7 based on the determination of activation energies. Moreover, we will present an original analysis based on the mean square displacement of molecules to investigate the assumptions made by Kahn et al. 13 and Papadopoulos et al. 14 Our paper will be organized as follows: after presenting the faujasite models and the simulation details, we will present the results obtained regarding the self-diffusion, the supercage heterogeneity and the activation energy of diffusion processes.

4

ACS Paragon Plus Environment

Page 4 of 29

Page 5 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Faujasite model and simulation method Faujasite NaX model The

simulated

zeolite

NaX

includes

86

Na

atoms

per

unit

cell

(Na86 Al86 Si106 O384 ) and has a cubic structure. 15 The structure is composed of supercages (diameter ≈ 13 Å) with a window opening of 7.4 Å and smaller sodalite cages (diameter ≈ 6.6 Å) with a window opening of 2.2 Å interconnected through hexagonal prisms (Fig.1). The position of the compensating cations Na+ inside the structure are located to inhibit the sodalite cages access for most sorbate molecules but water.

Molecular dynamics settings By MD we simulated eight elementary unit cells (Lx = 2a, Ly = 2a, Lz = 2a) of the faujasite NaX which has a crystal lattice parameter of a = 25.1 Å. 16 According to T. Bromley and coworkers 17 a flexible s tructure w as u sed. T his w as m otivated b y t he f act t hat molecular transport through a confining p orous e nvironment i s s trongly d ependent o n t he flexibility of the framework. Flexibility was imposed by considering bending and stretching of Si, Al, O and H atoms through the intramolecular potentials of the structure (Table 1). Also, we used three dimensional periodic boundary conditions. In all the performed simulations was applied a constant temperature T by using a Nosé-Hoover thermostat. The volume V and the number of atoms N are kept constant during the simulation (NVT). We set the time step, ∆t, to 5.0 ×10−4 ps and we simulate trajectories of 30 ns. For studying the loading effects on the transport properties we used several amounts of hydrogen (D2 or T2) ranging from 1 to 120 molecules per unit cell. The molecules of H2 , D2 and T2 are represented with the unified a tomic description (spherical model). The Na+ ions are linked to the atoms of the framework through: an electrostatic potential, 18 a 12-6 Lennard-Jones (LJ) potential (Table 2) and a Buckingham potential (Table 3). As mentioned in Sect. Introduction, for the faujasite NaX the difference 5

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 29

between pore diameter and the molecular hard core are comparable to the de Broglie wavelength of radial motion (λ = h/vradial ). In these conditions quantum effects are expected to arise. These non-classical effects and the use of cryogenic temperatures (e.g., 100K, 77K and 40K) advocate molecular sieving. 6,7 For considering quantum effects the path integral approach is the ideal method for this kind of studies. However, they have an overwhelming computational cost. Instead, we implemented in a MD code (developed in our group) the Feynman-Hibbs (FH) modified interaction potential: 8

VF H (r) =

6µ πβ¯ h2

!3/2 Z

6µ 2 dRU (|r + R|)exp − R πβ¯h2

!

(1)

by an expansion of the integrand in powers of R up to the 4th order, namely: β¯h2 VFH (rij ) = VLJ (rij )+ 24µ where µ =

mi ∗mj mi +mj

d2 VLJ (rij ) 2 dVLJ (rij ) 1 + + 2 drij rij drij 2 !

is the reduced mass, β =

1 , kB T

β¯ h2 24µ

!2

d4 VLJ (rij ) 4 d3 VLJ (rij ) + 4 3 drij rij drij (2)

kB is the Boltzmann constant h ¯ the Planck

constant and (

VLJ (rij ) = 4.0ij

σij rij

!12

σij − rij

!6 )

,

(3)

the Lennard-Jones potential (Table.2). Where rij is the distance between atoms i and j, σij and ij are the LJ parameters between atoms. The use of the FH potential implies, a reduction of the potential depth well and an increase of the kinetic diameter when the reduced mass and temperature are decreased. Note that the expression for the quartic Feynamn-Hibbs effective potential in Eq. 2 is slightly different than the one used in recent simulation studies on the subject. 7,10,19 For consistency, we re-deduced the expression from the general formula given by Sesé. 20 The quantitative difference between the previous and the corrected formulae are found to be negligible. As was demonstrated in Ref., 7 in the case of confined fluids (which is the case when the NaX structure is above 50% of its loading capacity) the quartic approximation to the FH potential is more accurate and this is because

6

ACS Paragon Plus Environment

!

Page 7 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

the dependence on temperature and mass give a more accurate determination of the quantum corrections. The quadratic approximation is enough for a dilute number of H2 molecules. Here, is used Eq.2 for the interactions between H2-H2, H2-O and H2-Na+ and similarly for the equivalent interactions of D2 and T2. We used a cutoff radius for all intermolecular interaction potentials of 13 Å.

H2 -Na interaction energy Previous studies on the diffusion of H2 where the guest structure is the faujasite NaX have used a force-field for the H 2 -Na+ i nteraction 12 with a set parameters ( H2 −N a /kb = 2 65K and σH2 −N a = 2.8 Å) which render the potential more attractive when compared to those used in our studies (Table.2). In the former set of parameters (adjusted to reproduce the results obtained by QENS), imply that the Na+ -H2 interaction will be favored at the expenses of other interactions inside the pore. The consequence is a higher adsorption at the beginning of the process and a mobility reduction of absorbent molecules. The adsorption difference observed between the parameters used in Ref. 12 and our results, lead us to perform DFT simulations with the VASP 21 package. Our objective being to validate the parameters we used for the H2 -Na+ interaction. With these simulations we do not search to re-parametrize the interaction potential but to have an energetic criterion (i.e., adsorption energy) comparable to the experimentally obtained. For this, the H2 -Na+ interaction within a faujasite NaX have been analyzed with the PBE GGA functional 22 associated with several additive dispersion correction schemes implemented in the VASP package. 23 The PAW scheme 23,24 has been used to treat the electron-ion interactions and has proven its high accuracy and convenience in a considerable number of applications and systems. 21 The plane wave cutoff energy was set to 500 eV. Considering the large size of the unit cell, the Brillouin zone integration was performed with one Γ-point only. The Kohn-Sham equations were solved self-consistently until the energy difference between cycles becomes lower than 1E-8 eV. The atomic positions have been fully optimized until forces are smaller than 0.01 eV/Å per atom. 7

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 29

For computing with a high accuracy the interaction energies between molecules and the faujasite NaX, we have chosen to include for the missing London dispersion interactions in conventional Kohn-Sham DFT by using several vdW (van Der Waals) correction schemes. 25,26 These methods include pairwise additive corrections schemes of Grimme 27–29 where the C6 coefficients are tabulated (D2) 27 or corrected depending on the coordination number (D3). 28 In the D3-BJ version, the damping function has been improved. 29 Tkatchenko and Scheffler 30 have developed another approach (TS) based on the determination of C6 on fly, depending on the chemical environment. In this study, we will use the original TS scheme 31 as well as the TS/HI version including an iterative Hirshfeld partitioning developed and implemented in VASP by Bučko and co-workers. 32,33 Finally, the MBD/FI 34 goes beyond the simple pairwise correction by taking into account many-body interactions as in the original MBD method, 35 but includes correctly the atoms ionicity. In particular, these three methods are available in the VASP code thanks to recent implementations. 31–34,36 Once the total energy of the different systems is obtained with or without dispersion correction schemes, the interaction energies of H2 -NaX at 0 K are calculated by using the following formula: ∆E int = E(N aX−X) − EN aX − EX ,

(4)

with E(N aX) - energy of the empty faujasite, EX - energy of the isolated molecule in gaseous phase and E(N aX−X) -energy of the faujasite with the adsorbed molecules. The contribution of dispersion forces ∆E disp to the adsorption energy is computed similarly:

∆E disp = E disp (N aX−X) − E disp N aX − E disp X .

(5)

In this work, we exchange 22 Si atoms by 22 Al atoms inside this primitive cell and add 22 Na atoms distributed among the I (4Na), I0 (2 Na), II (8 Na) and III (8 Na) sites according to the known experimental distribution, 37 resulting in a Si/Al ratio of 1.18 (Fig.1) . Given that H2 or D2 isotopes are known to be adsorbed in the NaX supercages, we

8

ACS Paragon Plus Environment

Page 9 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

focused on site II and III where we placed H2 molecules in the neighborhood of the Na+ cations located in these sites. The choice of site III is motivated for calculating the evolution of the interaction energy of H2 as a function of the Na-H2 distance: since site III is more distant from the framework than site II, the evolution trend can be followed with ease and avoid screening appearing due to the presence of other cations. In Table.4 are given the results for site III. Worthwhile to mention, is that the interaction energies obtained for H2 on site II or site III are similar, indeed. The PBE+vdW computed values are in a narrow interval: from -11.5 kJ/mol for PBE+D3 to -15.7 kJ/mol for PBE+D2. PBE+D3 gives the closest value to the experimental one of -6 kJ/mol determined in Ref. 10 PBE+TS/HI gives also an interesting value of -11.6 kJ/mol, showing that both Grimme and Tkatchenko and Scheffler’s approaches provide similar results. The gap of 5 kJ/mol between our DFT values and experimental values may be originated by: (i) the temperature, which is 0K for VASP and 77K for experiments, (ii) the lack of a description of relativistic effects in our DFT calculations and the quantum effect of the molecular core and iii) the inhomogeneity of the supercage. In our calculations when a H2 approaches to either site II or site III the calculated energy includes the atoms surrounding the cations. Nonetheless, our findings with DFT dispersion-corrected methods give an interaction energy quantitatively close to that obtained experimentally and by MD simulations. This result, allow us to say that we have to our disposal a reliable set of parameters for the H2-Na+ interaction (Table.2).

Results In previous studies we developed a comparison of experimental adsorption isotherms with those obtained from MD simulations at cryogenic temperatures for the faujasite NaX. 10 Here, we focus on the dynamical aspects for elucidating the undergoing kinetic processes. This is done by comparaison between from QENS data 12 and MD simulations.

9

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 29

Self-diffusion coefficient The self-diffusion coefficients, Dself , are calculated at equilibrium conditions for the faujasite NaX by using the Einstein’s relation: D

E

lim [r(t) − r(0)]2 = 6N Dself t,

t→∞

(6)

where r(t) are the H2 (or isotopes) molecular positions in time (t). hXi represent the time average over X. The calculated Dself by using the trajectories obtained with our MD code (including the Feynman-Hibbs modified potential (Eq.2) and the parameters for the LJ-potential given in Tables. 2 and 3 are given in Fig.3. As expected, the figure shows an increase of Dself with temperature. At 100K and 77K, the Dself of the lighter compound are higher than those for D2 H2 the heavier. However, at 40K and low loading this behavior is inverted (Dself > Dself ). This

can be interpreted as a quantum effect and is in agreement with previous observations by Kumar and Bathia. 38 By increasing the supercage loading, the mean free path and diffusivity (Eq.6) of molecules is decreased as shown in Fig.3.a. The reported self-diffusion coefficients for H2 are in agreement with those obtained by neutron scattering by Kahn et al.,. 13 The behavioral scenario reported in Ref., 12 for the faujasite NaX at 100K by MD simulations, shows an increase of Dself by a factor 2 when passing from low to relatively high loading. This behavior, may be interpreted by a trapping of the first arriving guest molecules at the supercage walls and producing a screening of subsequent H2 -Na and H2 -O interactions. This last favoring the H2 - H2 interactions and resulting on an increase of Dself for H2 . It is worthwhile to notice in Fig.4 of 12 that the maximum number of particles per unit cell used in their simulations is 64 corresponding roughly to 50% of the maximum loading. In our simulations, we increased the loading above 50% giving a decrease of Dself (Fig.3.b). This last is produced by the increased packing of molecules inside the pore which reduce their mean free path.

10

ACS Paragon Plus Environment

Page 11 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Supercage heterogeneity The increase of diffusivity at low loading by the adsorption of H2 within a faujasite NaX is, undoubtedly, an intriguing issue. Some authors privileged the idea of adsorption sites 14 where molecules will be fixed and the rest of molecules free to move inside the c avity. While other authors 13 advanced the hypothesis of two populations of molecules: slow and fast. By slow meaning, particles that reside a period of time close to the must attractive regions inside the supercage (e.g., the pore walls) and by fast particles those moving freely inside the pore. The last authors estimated the residence time, between isotropic jumps of length of ≈ 3.9 Å , to τ0 ≈ 11 ps at 100K. All authors agree to state that the supercage heterogeneity is responsible for an increase of Dself . This issue is, from an industrial point of view, crucial for orienting the choice of both the guest structure and the role of the substitution cation for increasing the selectivity. We performed an analysis based on the particle distribution inside the supercage (α). The distributions obtained with the parameters given in Table.2 and Table.3 show at 100K a maximum at σ c ∼ 3.5 Å from the α-center, corresponding to a position close to the cage wall. At high loading we obtain the same maximum but with a spreading of H2 from α-center and σ c (Fig.4.a). At 40K and for low and intermediate loading most of the H2 molecules are situated near the supercages walls, while at high loading we observe a concentration of molecules around the α-center (Fig.4.b). Our analysis did not permit to highlight the presence of privileged adsorption sites within the faujasite NaX. We used the H2 trajectories obtained by MD simulation for calculating the statistical distribution of mean square displacement (msd) for every molecule at different temperatures and loadings. The analysis were performed by taking time intervals of 500 ps and 5000 ps and revealed a surprising result. Fig5a shows the case at 500ps for 40K, 77K and 100K at high loading (near saturation). In this figure we can distinguish below the msd of 2000 Å2, and for all the cases analyzed, two peaks corresponding to particles with different mobility. By performing the same analysis but with an average every 5000ps(Fig.5b), there are not 11

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 29

well defined peaks. Nonetheless, at 40K we observed at the top of the curve peaks which correspond, probably, to different residence times near the pore wall. This last is also observed at 77K and at low loading with an average of 500ps (Fig.5.c), while at 5000ps the presence of characteristic peaks has desappeard (Fig.5.d). A similar behavior to the case of 77K was observed at 100K, indeed. This analysis highlights the presence of slow and fast populations of particles governing, partially, the diffusivity inside the supercage. Obviously, the density of slow molecules is temperature dependent: higher the temperature lower the densities of slow particles. Similar results were obtained for D2 and T2. Our analysis shows that the pore heterogeneity leads to the presence of slow and fast molecules which at low loading result on an increase of the self diffusivity. This result is not in contradiction with the suggestion of Papadopoulos 14 of molecules spending extra time on privileged sites.

Activation Energy For a given loading the dependence of Dself with temperature is given by: Ea D(T ) = D0 exp − RT 



,

(7)

where D0 is a pre-exponential factor and R the ideal gas constant. From this equation, we can determine the activation energy Ea . This quantity energetically tell us the difficulty of a molecule for moving within its environment, roughly speaking: larger the molecule is, higher is the energy needed to diffuse. The ensemble of D0 and Ea calculated at different loadings are gathered in Tables.5 and 6. The statistical uncertainties are lower than 20%. The results in Fig.6 for different coverages show an increase of the activation energy with loading. Despite the statistical uncertainties, the trend for all species is essentially linear. By a linear extrapolation at zero coverage we obtain 0.32 kJ/mol for H2 , 0.15 kJ/mol for D2 and 0.12 kJ/mol for T2 . This means that, H2 becomes more cumbersome than D2 and T2 when

12

ACS Paragon Plus Environment

Page 13 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

temperature is decreased and is reflected by a higher mobility of the isotopes. This explains that at low temperatures (e.g., 40K) the self-diffusion coefficient of T2 or D2 are higher than H2 . It is important to highlight the increase of the activation energy with loading (Fig.6). Henceforth, and by Eq.7, a decrease of the self-diffusivity as depicted in Fig.3.a. For analyzing the cooperative diffusive behavior given by the Fick’s diffusion law, we compute the Maxwell-Stefan’s and the Fick’s diffusion coefficients. The former coefficient can be deduced from the following equation: D

E

lim [R(t) − R(0)]2 = 6N DM S t,

t→∞

with R(t) =

PN

i=1

(8)

ri (t). The Maxwell-Stefan’s (DM S ) and the Fick’s diffusion (DF )

coefficients are linked by: DF = DM S Γ,

(9)

where Γ is a thermodynamic factor 39 given by: 1 Γ= RT

∂µ ∂lnNa

!

= T,V


, > − < Na >2

(10)

where µ and Na are the chemical potential and the number of molecules adsorbed per unit cell, respectively . We can determine Γ−1 either from the adsorption isotherms or from the instantaneous fluctuations of Na . 40 By using the last approach for the MD trajectories, up to 800 ps, we obtained for low coverage that all curves tend to 1 and they are indistinguishable (Fig.7). By increasing the loading, the curves split from each other and tend to 0 (saturation). Contrary to a Langmuir’s isotherm, the obtained trend is not linear. For a given coverage, the value of Γ−1 increases with temperature and the mass of the isotope. However, at higher temperatures (where quantum corrections stop to be relevant) this is less visible. This can be explained by the fact that the kinetic radius of H2 is increased at low temperatures.

13

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

In Figs.8 are given the Maxwell-Stefan’s and Fick’s diffusion coefficients computed from the simulations performed with our own code. It is well-known that these coefficients, when calculated from the same trajectories, exhibit larger statistical uncertainties than selfdiffusion coefficients. For the same coverage DM S increases with temperature and decreases as the isotope mass increases. Nonetheless, it is rather constant with loading (Figs.8.b and 9). Due to the thermodynamic factor Fick’s diffusion coefficient tends to increase with loading (supposed to take an infinite value at saturation) in agreement with QENS experiments (Fig.8.a). Notice that Dself , DM S and DF tend to the same value at low loading, as expected, and illustrated in Fig.9 for H2 at 40K.

Conclusions A crucial issue when performing molecular dynamics simulations concerns the reliability of the force field u sed f or a g iven i nteratomic i nteraction. T he m odification of th e potential parameters is a subtle task and in some cases can lead to misleading interpretations, since two set of parameters can provide similar adsorption isotherms. In our studies for the sake of having a reliable set of force-field parameters, we developed an energetic criterion based on a threefold method comparison (i.e., experimental, DFT and classical molecular dynamics) of the interaction energy. 10 The three methods gave approximately the same value of -6 kJ/mol validating the parameters we used for the H2-Na+ interaction. The analysis of the H2 molecules distribution inside the supercages, for low loading, shows the presence of a peak located approximately at 3.5 Å from the supercage center where H2 molecules are distributed. By increasing the loading, molecules spread around this center (Fig.4) and favor the adsorbate-adsorbate interactions. Interestingly, our analysis based on the MSD of the trajectories of H2 show the presence of two populations of particles slow and fast in agreement with the experimental finding of Kahn et al.,. 13 It is not excluded that the presence of slow and fast particles be compatible with the existence of privileged adsorption

14

ACS Paragon Plus Environment

Page 14 of 29

Page 15 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

sites, 14 but for validating this ansatz an enhanced representation of the cage heterogeneity is unavoidable. This task requires an enhancement of the adsorbate-adsorbent interactions, since they determine the selectivity for a mixture of gases. We have shown that in the temperature range from 40K to 100K and by increasing the loading up to saturation, the self diffusion coefficient behaves as follows: i) An initial increase of Dself at low loading (Fig.3.b) - which can be explained by an enhanced adsorption of H2 or its isotopes at the walls of the NaX supercages (Fig.4), produced by the heterogeneity of the pore, and leaving some molecules free to move at the center of the pore which in average are those responsible of the increase of the average diffusivity. ii) By increasing the loading (even when some molecules are close to the pore walls) the average adsorbate-adsorbate interaction increases, this unavoidably leads to a reduction of molecules mean free path and reflected on a decrease of the average diffusivity (Fig.3). At 40K our studies show a faster diffusivity for D2 and T2 than for H2 in agreement with previous simulations 7 and experimental findings. 12 Our analysis based on the activation energies shows that this diffusion inversion can be explained by the fact that H2 becomes more cumbersome (i.e., mobility reduction in its environment) than its isotopes. This mobility reduction is originated by quantum effects which are accentuated as the temperature is decreases (Fig.6). Our approach gives a different point of view which can be validated from experimental measurements. For us, this sort of validation is essential since our objective is to obtain and elucidate physicochemical information useful for the development of molecular sieving protocols of H2 and its isotopes. An original issue of our studies has been to give an energetic interpretation of the isotopic quantum effects over the diffusion coefficients based on the activation energy of diffusion. The adsorption energy determines the particles density inside a porous structure while the activation energy allows to explain energetic aspects of molecular transport. By definition the activation energy quantifies t he e volution o f D

self w ith

t emperature. We h ave shown,

that by lowering the temperature the quantum effects are enhanced (as expected) and this is

15

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

reflected on the activation energies values of the different species which take different values. The interest of the results reported in this manuscript, is to provide information about the kinetic processes developing within a zeolitic structure. Moreover, they can orient the choice of a guest structure or to develop specific tailored structure for recycling the gases produced by nuclear plants.

Acknowledgement The authors thanks the data center of the Université de Bourgogne for giving access to their computing facilities. B.R. thanks the C.E.A. for financial support.

16

ACS Paragon Plus Environment

Page 16 of 29

Page 17 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

References (1) Salazar, J. M.; Simon, J.-M. Two prototype of this kind of reactors are ITER (International Thermonuclear Experimental Reactors) and DEMO (Demonstration Power Station). 2015. (2) Harris, W. Separation Of hydrogen Isotopes By diffusion Through Palladium. Proc. Natl.Acad. Sci. 1933, 19, 4278–4289. (3) Oh, H.; Hirsher, M. Quantum Sieving for Separation of Hydrogen Isotopes Using MOF’s. Eur. J. Inorg. Chem. 2016, 27, 4278–4289. (4) Miller, A.; Girdler, Heavy Water: A Manufacturers Guide for the Hydrogen Century. Can. Nucl. Soc. Bull. 1983, 22, 1–14. (5) Saeki, M. Influence of radiation damage on diffusivity of tritium in graphite. Int. J. Appl. Radiat. Isot. 1983, 24, 739–742. (6) Beenakker, J.; Borman, V.; Krylov, S. Molecular transport in subnanometer pores: zero-point energy, reduced dimensionality and quantum sieving. Chem. Phys. Letters 1995, 232, 379–382. (7) Kumar, A.; Jobic, H.; Bhatia, S. Quantum effects on adsorption and diffusion of hydrogen and deuterium in microporous materials. J. Phys. Chem. 2006, 110, 16666–16671. (8) Feynman, R.; Hibbs, A. Quantum Mechanics and Path Integrals.; Dover Publications, Inc. New York., 1965. (9) Perez-Carbajo, J.; Parra, B. J.; Ania, C.; Merkling, P.; Calero, S. Molecular Sieves for the Separation of Hydrogen Isotopes. ACS Appl. Mater. Interfaces 2019, 11, 18833– 18840.

17

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(10) Salazar, J.-M.; Lectez, S.; Gauvin, C.; Macaud, M.; Bellat, J.; Weber, G.; Bezverkhyy, I.; Simon, J.-M. Adsorption of hydrogen isotopes in the zeolite NaX: experiments and simulations. Int. J. Hyd. En. 2017, 42, 13099–13110. (11) Giraudet, M.; Bezverkhyy, I.; Weber, G.; Dirand, C.; Macaud, M.; Bellat, J.-P. D2/H2 adsorption selectivity on FAU zeolites at 77.4 K: Influence of Si/Al ratio and cationic composition. Microporous and Mesoporous Materials 2019, 270, 211–219. (12) Pantatosaki, E.; Papadopoulos, G. K.; Jobic, H.; Theodorou, D. Combined atomistic simulation and quasielastic neutron scattering study of the low-temperature dynamics of hydrogen and deuterium confined in NaX zeolite. J. Phys. Chem. B 2008, 112, 11708–11715. (13) Kahn, R.; Cohen De Lara, E.; Viennet, E. Diffusivity of the hydrogen molecule sorbed in NaA zeolite by a neutron scattering experiment. The Journal of Chemical Physics 1989, 91, 5097. (14) Papadopoulos, G.; Jobic, H.; Theodorou, D. Transport Diffusivity of N2 and CO2 in Silicalite: Coherent Quasielastic Neutron Scattering Measurements and Molecular Dynamics Simulations. The Journal of Chemical Physics B 2004, 108, 12748–12756. (15) Olson, D. The crystal structure of dehydrated NaX. Zeolites 1995, 15, 439–443. (16) Jeffroy, M.; Draghi, C.; Boutin, A. Molecular simulation of zeolite flexibility. Molecular Simulation 2014, 40, 6–15. (17) van den Berg, A.; Bromley, T.; Ramsahye, N.; Maschmeyer, T. Diffusion of Molecular Hydrogen through Porous Materials: The Importance of Framework Flexibility. J. Phys. Chem. B 2004, 108, 5088–5094. (18) Fennell, C.; Gezelter, J. Is the Ewald summation still necessary? Pairwise alterna-

18

ACS Paragon Plus Environment

Page 18 of 29

Page 19 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

tives to the accepted standard for long-range electrostatics. J. Phys. Chem. 2006, 124, 234104. (19) Kumar, A.; Jobic, H.; Suresh, K.; Bhatia, K. Quantum Effect Induced Kinetic Molecular Sieving of Hydrogen and Deuterium in Microporous Materials. Adsorption 2007, 13, 501–508. (20) Sesé, L. Feynman-Hibbs Potentials and Path Integrals for Quantum Lennard-Jones System: Theory and Monte Carlo Simulations. Molec. Phys. 1995, 85, 931–947. (21) Hafner, J. Ab-Initio Simulations of Materials Using VASP: Density-Functional Theory and Beyond. J. Comput. Chem. 2008, 29, 2044–2078. (22) Perdew, J.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys Rev Lett 1996, 77, 3865–3868. (23) Kresse, G.; Hafner, J. Ab initio molecular dynamics for liquid metals. Phys Rev B 1993, 47, 558–561. (24) Blöchl, P. Projector Augmented-Wave Method. Phys. Rev. B 1994, 50, 17953–17979. (25) Chebbi, M.; Chibani, S.; Paul, J.; Cantrel, L.; Badawi, M. Evaluation of volatile iodine trapping in presence of contaminants: A periodic DFT study on cation exchangedfaujasite. Microporous Mesoporous Mater 2017, 239, 111–122. (26) Hessou, E.; Kanhounnon, W.; Rocca, D.; Monnier, H.; Vallières, C.; Lebègue, S.; Badawi, M. Adsorption of NO, NO2, CO, H2O and CO2 over Isolated Monovalent Cations in Faujasite Zeolite: A Periodic DFT Investigation. Theor Chem Acc 2018, 137, 161. (27) Grimme, S. Semiempirical GGA-type density functional constructed with a long-range dispersion correction. J Comput Chem 2006, 27, 1787–1799.

19

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(28) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu. Journal of Chemical Physics 2010, 132, 154104. (29) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the Damping Function in Dispersion Corrected Density Functional Theory. J. Comput. Chem 2011, 32, 1456–1465. (30) Tkatchenko A, S. M. Accurate Molecular Van Der Waals Interactions from GroundState Electron Density and Free-Atom Reference Data. Phys Rev Lett 2009, 102, 073005. (31) Bučko, T.; Lebègue, S.; Hafner, J.; Ángyan, J. G. Tkatchenko-Scheffler van Der Waals Correction Method with and without Self-Consistent Screening Applied to Solids. Phys. Rev. B 2013, 87, 064110. (32) Bučko, T.; Lebègue, S.; Hafner, J.; Ángyan, J. G. Improved Density Dependent Correction for the Description of London Dispersion Forces. J Chem Theory Comput 2013, 9, 4293–4299. (33) Bučko, T.; Lebègue, S.; Ángyan, J. G.; Hafner, J. Extending the applicability of the Tkatchenko-Scheffler dispersion correction via iterative Hirshfeld partitioning. Journal of Chemical Physics 2014, 141, 034114. (34) Gould, T.; Lebègue, S.; Ángyan, J. G.; Bučko, T. Fractionally Ionic Approach to Polarizability and van der Waals Many-Body Dispersion Calculations. J Chem Theory Comput 2016, 12, 5920–5930. (35) Tkatchenko, A.; DiStasio, R.; Car, R.; Scheffler, M. Accurate and Efficient Method for Many-Body van der Waals Interactions. Physical Review Letters 2012, 108, 236402. (36) Bučko, T.; Hafner, J.; Lebègue, S.; Ángyan, J. G. Improved Description of the Structure

20

ACS Paragon Plus Environment

Page 20 of 29

Page 21 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

of Molecular and Layered Crystals: Ab Initio DFT Calculations with van Der Waals Corrections. J. Phys. Chem. A 2010, 114, 11814–11824. (37) Frising, T.; Leflaive, P. Extraframework Cation Distributions in X and Y Faujasite Zeolites: A Review. Microporous and Mesoporous Materials 2000, 114, 27–63. (38) Kumar, A.; Bhatia, S. Quantum Effect Induced Reverse Kinetic Molecular Sieving in Microporous Materials. Phys. Rev. Lett. 2005, 95, 245901. (39) Allen, M.; Tildesley, D. Computer Simulation of Liquids; Oxford University Press: New York, 1984. (40) Simon, J.-M. Non-equilibrium Thermodynamics Applied to Adsorption. Experimental Thermodynamics Volume X : Non-equilibrium Thermodynamics with Applications 2015, 178–203. (41) Smit, B.; Maesen, T. Molecular simulations of zeolites: adsorption, diffusion, and shape selectivity. Chemical Review 2008, 108, 4125–4184. (42) Uytterhoeven, L.; Dompas, D.; Mortier, W. J. Theoretical investigations on the interaction of benzene with faujasite. J. Chem. Soc. Faraday Trans. 1992, 88, 2753–2760.

TOC Graphic

21

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 29

Figures and tables Table 1: Expressions for the stretching and bending interaction potentials used for the flexible zeolite structure. Here, are given ks and req for stretching and kb and θeq for bending, 41 rij is the distance between bonded atoms i and j and θijk is the angle between atoms i, j and k. i-j Si-O Al-O i-j-k O-Si-O O-Al-O Si-O-Al Si-O-Si

Stretching ks /kB (×106 K Å−2 ) 1.39 8.12 Bending kb /kB (×105 K deg−2 ) 3.00 1.00 0.48 0.60

Vstretch (rij ) = k2s (rij − req )2 req (Å) 1.65 1.75 Vbend (θijk ) = k2b (θijk − θeq )2 θeq (deg) 110 110 145 145

Table 2: Lennard-Jones parameter sets used for the interaction between H2 -O, H2 -Na+ , H2 -H2 and Na+ -Na+ . The H2 molecules are considered neutral/nonpolar. Lennard-Jones 12-6 parameters  (K) σ (Å) i-j kb H2 -O 40.0 3.83 H2 -Na+ 34.0 3.80 H2 -H2 , 26.0 2.556 D2 -D2 and T2 -T2 Na+ -Na+ 50.270 2.586

Table 3: Parameters used in the MD simulations for the interaction Na+ -O with a Buckingham potential. (N a+ −O)





Buckingham potential :VBu

(rij ) = α exp − rβij −

α/kB (K) 61.1558 × 106

γ/kB ( K.Å6 ) 76.5898 × 104 42

βÅ 4.0518

22

ACS Paragon Plus Environment



γ 6 rij



Page 23 of 29 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Table 4: Interaction energies obtained for H2 as a function of the distance to site III, within the NaX supercage, by using different dispersion methods available in the VASP package. The methods PBE+D3 and PBE+TS/HI give the closest values to the experimental one. vdW correction schemes D2 D3 D3-BJ TS TS/HI MBD/FI

∆Eint (kJ/mol) -15.7 -11.5 - 13.2 -15.3 -11.6 -13.7

Table 5: Self-diffusion activation energy Ea for H2 and its isotopes at several loadings Self diffusion activation energy Molecules loading H2 80 0.31505 160 0.2885 240 0.3172 400 0.3585 560 0.4015 720 0.3382 1283 0.5905

Ea (kJ/mol) D2 0.2231 0.2486 0.2533 0.3038 0.3462 0.4025 0.5496

T2 0.1884 0.2068 0.2493 0.2959 0.3407 0.3856 0.5181

Table 6: D0 in Eq.7 for H2 and its isotopes at several loadings D0 (10−8 m2 /s) Molecules loading H2 D2 80 20.0956 14.3263 160 17.0321 13.4961 240 13.0553 12.3938 400 10.5055 10.3563 560 7.5776 8.3536 720 2.0252 7.0006 1283 3.6594 2.7765

23

ACS Paragon Plus Environment

T2 12.4809 12.089 11.5236 9.4281 7.5127 6.1092 2.7723

The Journal of Physical Chemistry

Figure 1: Schematic representation of the faujasite NaX. The elementary unit cell is composed of 8 supercages, 8 sodalites cages and 16 prismes D6R.Figure taken from Ref. 10

4 2 Na-H2 energy (kJ/mol)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 29

0 -2 -4 -6 -8 -10 -12 1.5

2

2.5

3

4 3.5 4.5 distance (Å)

5

5.5

6

Figure 2: The graph gives the interaction energy of H2 as a function of the distance to site III of the faujasite NaX (Fig.1). These energies are not due solely to the H2 -Na+ , but includes a contribution of the other atoms forming the supercage.

24

ACS Paragon Plus Environment

Page 25 of 29

14

H2 40 K H2 77 K H2 100 K D2 40 K D2 77 K D2 100 K T2 40 K T2 77 K T2 100 K

10

-8

2

Dself (10 m /s)

12

8 6 4 2 0 0

20

40 60 coverage (molec./u.c.)

80

100

Figure 3: Self-diffusion coefficients as a function of coverage at different temperatures for H2 and its isotopes obtained by MD simulations with our in-house code.

90

25 H2-α cages 100 K N=16 H2-α cages 100 K N=280 H2-α cages 100 K N=719 H2-α cages 100 K N=1260

20

H2-α cages 40 K N=16 H2-α cages 40 K N=280 H2-α cages 40 K N=719 H2-α cages 40 K N=1260

80 70 60 g(r)

15 g(r)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

10

50 40 30 20

5

10 0 0

a)

0,1

0,2

0,3

0,4

0,5 0,6 r (nm)

0,7

0,8

0,9

0 0

1

0,1

0,2

0,3

0,4

0,5 0,6 r (nm)

0,7

0,8

0,9

1

b)

Figure 4: Distribution of H2 molecules inside the supercage (α), respectively at a) 100K and b) 40K for the potential given in Table.2.

25

ACS Paragon Plus Environment

The Journal of Physical Chemistry

0.0006

6e-05 40K, 719 mol 77K, 719 mol 100K, 719 mol

0.0005

Histogram

Histogram

4e-05

0.0003

3e-05

0.0002

2e-05

0.0001

1e-05

0

40K, 719 mol 77K, 719 mol 100K, 719 mol

5e-05

0.0004

0

2000

4000

8000

6000

0

10000

0

10000

20000

2

MSD / (Å )

a)

0.0003

40000

50000

60000

3e-05 77K, 16 mol 77K, 130 mol 77K, 719 mol

0.00025

Histogram

2e-05

0.00015

1.5e-05

0.0001

1e-05

5e-05

5e-06

0

77K, 16 mol 77K, 130 mol 77K, 719 mol

2.5e-05

0.0002

c)

30000 2 MSD / (Å )

b)

Histogram

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 29

0

5000

10000

15000 2 MSD / (Å )

20000

25000

0

30000

0

50000

1e+05

1.5e+05 2 MSD / (Å )

2e+05

2.5e+05

3e+05

d)

Figure 5: Distribution of molecular mean square displacement of H2 for two characteristic time: 500 ps (a and c) and 5000 ps (b and d). a) and b) show the effect of the temperature at high loading and c) and d) the effect of the coverage.

26

ACS Paragon Plus Environment

0.6 0.5

H2 D2 T2

0.4 0.3 0.2 0.1 0 0

20

40 60 coverage (molec./u.c.)

80

100

Figure 6: The figures gives the activation energy, Ea , calculated from Eq. 7 as a function of the coverage with a statistical uncertainties below 20%.

1 H2 40 K H2 77 K H2 100 K D2 40 K D2 77 K D2 100 K T2 40 K T2 77 K T2 100 K

0.8

-1

0.6 Γ

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

activation energy (kJ/mol )

Page 27 of 29

0.4 0.2 0 0

50

100 150 coverage (molec./u.c.)

200

Figure 7: The inverse of the thermodynamic factor Γ for H2 plotted as a function of coverage at different temperatures. Indeed, similar results were obtained for D2 and T2 , no depicted in the graph.

27

ACS Paragon Plus Environment

The Journal of Physical Chemistry

-8

2

Fick diffusion coefficient (10 m /s)

60 50 40 30

H2 40 K H2 77 K H2 100 K D2 40 K D2 77 K D2 100 K T2 40 K T2 77 K T2 100 K

20 10 0 0

20

a)

40 60 coverage (molec./u.c.)

100

H2 40 K H2 77 K H2 100 K D2 40 K D2 77 K D2 100 K T2 40 K T2 77 K T2 100 K

-8

b)

80

20

2

Maxwell diffusion coefficient (10 m /s)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 29

15

10

5

0 0

20

40 60 80 coverage (molec./u.c.)

100

120

Figure 8: a) Fick’s diffusion and b) Maxwell-Stefan’s coefficients as a function of coverage per unit cell (molec/u.c.) and at different temperatures for H2 , D2 and T2 . The depicted results exhibit large dispersions and it is only possible to extract the general trends. The statistical uncertainties for the calculated coefficients is ≈ 20%.

28

ACS Paragon Plus Environment

Page 29 of 29

40 Dself DMS DF

35 30 25

-8

2 -1

Diffusion / (10 m .s )

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

20 15 10 5 0

0

10

20

30

40 50 60 -1 Coverage / (mol.uc )

70

80

90

100

Figure 9: Comparaison of the self, Fick and Maxwell-Stefan diffusion coefficient of H2 as the function of the coverage at 40 K.

29

ACS Paragon Plus Environment