Julien Morthomas
email: julien.morthomas@insa-lyon.fr

William Gonçalves

Michel Perez

Geneviève Foray

Christophe L Martin

Patrice Chantrenne

A novel method to predict the thermal conductivity of nanoporous materials from atomistic simulations

Keywords: silica aerogels, thermal conductivity, nanoscale, molecular dynamics

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Silica aerogels are highly porous materials (more than 90% porosity) with pore size distribution of the order of a few nanometers. Such a size gives silica aerogels extremely low values of thermal conductivity, which is interesting if they are used as high-performance thermal insulators. However, these properties are obtained at the expense of their mechanical properties. Due to their nanostructure, the understanding and the prediction of their properties require atomic scale simulations. As silica aerogels are a porous heterogeneous material, the representative volume of the matter that has to be considered to predict its properties should contain several millions of atoms. With such a number of atoms, molecular dynamics simulations remain appropriate. However, the choice of an interatomic potential is challenging since the computational time has to be minimized.

Silica is a polar material, each atom thus carries an equivalent charge, inducing long range Coulombic interactions. Those have to be accounted for to calculate interatomic forces. Electrostatic forces have long-range interactions, which are costly in term of computation time. Using a BKS potential [START_REF] Van Beest | Force fields for silicas and aluminophosphates based on ab initio calculations[END_REF], Rajappa et al. [START_REF] Rajappa | A molecular dynamics study of ambient and high pressure phases of silica: Structure and enthalpy variation with molar volume[END_REF] have shown that the contribution of long-range Coulomb interactions to the total energy is lower for the amorphous phase than for the crystalline phase at the same density. This result is coherent with the previous work done of Carré et al. [3] who proposed a BKS potential with truncated coulomb interactions to successfully predict the static and dynamic properties of bulk silica. Recently, Gonçalves et al. [START_REF] Gonçalves | Molecular dynamics simulations of silica surface properties with truncated Coulomb interactions[END_REF] investigated the properties of (SiO2)n clusters and silica aggregates. For small aggregates of silica, the truncated BKS potential leads to the same 1 structural properties and surface energies as compared to the original BKS potential. This means that significant reduction of computation time (evaluated to a gain factor of 3000 on the CPU time per atom per time-step) is possible even when the ratio between the surface to the volume ratio of material increases drastically as for silica aerogels. Gonçalves et al. [START_REF] Gonçalves | Elasticity and strength of silica aerogels: A molecular dynamics study on large volumes[END_REF][START_REF] Gonçalves | Nanocompression of secondary particles of silica aerogel[END_REF] further used the truncated BKS potential to generate realistic silica aerogel structures. They have shown that the representative volume that has to be used to predict reliable mechanical properties of such material is of the order of 80 3 nm 3 . This is much larger than the maximum volume (20 3 nm 3 ) previously reached to predict mechanical behavior of silica aerogel from molecular dynamics simulations [START_REF] Lei | Determination of the Young's modulus of silica aerogels -an analyticalnumerical approach[END_REF][START_REF] Rivas Murillo | Structure and mechanical properties of silica aerogels and xerogels modeled by molecular dynamics simulation[END_REF][START_REF] Yeo | Silica Aerogels: A Review of Molecular Dynamics Modelling and Characterization of the Structural, Thermal, and Mechanical Properties[END_REF].

Heat transfer in silica aerogel is a combination of radiative and conductive heat transfer, both of which are coupled. Coquart et al. [START_REF] Coquard | Modelling of the conductive heat transfer through nano-structured porous silica materials[END_REF] and Wei et al. [START_REF] Wei | Thermal conductivities study on silica aerogel and its composite insulation materials[END_REF] predicted the equivalent thermal conductivity of aerogels, solving the coupled radiative transfer and conduction heat transfer equations. Coquart et al. [START_REF] Coquard | Modelling of the conductive heat transfer through nano-structured porous silica materials[END_REF] assumed that aerogels are made of a net of connected silica spheres, which diameter varies between 5 and 10 nm and they accounted for a variation of the thermal conductivity with the sphere diameter. This is not coherent with experimental observation [START_REF] Roiban | Advanced three dimensional characterization of silicabased ultraporous materials[END_REF] showing that silica aerogels may be considered as (SiO2)n clusters connected via nanometric SiO2 ligaments. Several authors have already used molecular dynamics simulations to predict the thermal conductivity of silica aerogel [START_REF] Yeo | Silica Aerogels: A Review of Molecular Dynamics Modelling and Characterization of the Structural, Thermal, and Mechanical Properties[END_REF][START_REF] Coquil | Molecular dynamics study of the thermal conductivity of amorphous nanoporous silica[END_REF][START_REF] Ng | A molecular dynamics study of the thermal conductivity of nanoporous silica aerogel, obtained through negative pressure rupturing[END_REF][START_REF] Yeo | Enhanced thermal characterization of silica aerogels through molecular dynamics simulation[END_REF][START_REF] Huang | Thermal conductivity of amorphous and crystalline thin films by molecular dynamics simulation[END_REF]. In these studies, due to computation time, the volume of the systems was smaller than 20 3 nm 3 and the minimum density is higher than the typical density of silica aerogels (around 70 -250 .

). Thus these system are not large enough to accurately describe the typical pore size distribution (centered on 10nm) observed experimentally [START_REF] Roiban | Advanced three dimensional characterization of silicabased ultraporous materials[END_REF][START_REF] Wang | Effect of preparation parameters on the texture of SiO2 aerogels[END_REF][START_REF]Aerogels handbook[END_REF]. Truncated BKS potential proposed by Carré et al. [START_REF] Carré | Amorphous silica modeled with truncated and screened Coulomb interactions: A molecular dynamics simulation study[END_REF] could open new avenue to calculate silica aerogel thermal conductivity at a reasonable cost. Although if this interatomic potential leads to appropriate surface properties [START_REF] Gonçalves | Molecular dynamics simulations of silica surface properties with truncated Coulomb interactions[END_REF], nanoporous silica aerogel structures and pore size distribution [START_REF] Gonçalves | Elasticity and strength of silica aerogels: A molecular dynamics study on large volumes[END_REF][START_REF] Gonçalves | Nanocompression of secondary particles of silica aerogel[END_REF], it does not warrant that the value of the thermal conductivity is reliable. In addition, in order to calculate the thermal conductivity of the silica aerogels the same conclusion could be made as for the mechanical properties, where it is necessary to study a sufficiently representative volume close to 100 3 .

Two methods based on molecular dynamics simulations are widely used to determine the thermal conductivity from molecular dynamics simulations [START_REF] Volz | Microscale and Nanoscale Heat Transfer[END_REF]: Equilibrium Molecular Dynamics (EMD) and Non-Equilibrium Molecular Dynamics (NEMD). In EMD simulations, the instantaneous variation of the heat flux in a system at equilibrium is calculated; the thermal conductivity is then extracted thanks to the Green Kubo Formulae [START_REF] Allen | Computer Simulation of Liquids[END_REF][START_REF] Frenkel | Chapter 4 -Molecular Dynamics Simulations[END_REF][START_REF] Ladd | Lattice thermal conductivity: A comparison of molecular dynamics and anharmonic lattice dynamics[END_REF][START_REF] Zwanzig | Time-Correlation Functions and Transport Coefficients in Statistical Mechanics[END_REF]. NEMD simulations are equivalent to well-known guarded hot plate experiments. It consists in simulating heat transfer in one direction of space by imposing a temperature difference between a hot zone and a cold zone. The technique to maintain the constant temperature in the hot and cold zones depends on the simulation strategy. Either the heat flux between the hot and cold zone is imposed [START_REF] Lukes | Molecular Dynamics Study of Solid Thin-Film Thermal Conductivity[END_REF][START_REF] Ikeshoji | Non-equilibrium molecular dynamics calculation of heat conduction in liquid and through liquid-gas interface[END_REF][START_REF] Jund | Molecular-dynamics calculation of the thermal conductivity of vitreous silica[END_REF], or the hot and cold zones are thermostated [START_REF] Maruyama | A molecular dynamics simulation of heat conduction in finite length SWNTs[END_REF][START_REF] Daly | Calculation of the thermal conductivity of superlattices by molecular dynamics simulation[END_REF][START_REF] Oligschleger | Simulation of thermal conductivity and heat transport in solids[END_REF][START_REF] Abs Da Cruz | Modified embeddedatom method interatomic potential and interfacial thermal conductance of Si-Cu systems: A molecular dynamics study[END_REF]. It has been shown by several authors that the two methods give consistent results [START_REF] Schelling | Comparison of atomic-level simulation methods for computing thermal conductivity[END_REF][START_REF] Poetzsch | Interplay of disorder and anharmonicity in heat conduction: Molecular-dynamics study[END_REF][START_REF] Landry | Complex superlattice unit cell designs for reduced thermal conductivity[END_REF][START_REF] Mahajan | Estimating thermal conductivity of amorphous silica nanoparticles and nanowires using molecular dynamics simulations[END_REF]. EMD exhibits larger uncertainties than NEMD since it requires a larger number of time steps to calculate an accurate correlation function of the instantaneous heat flux. Moreover, it appears that when the system is not homogeneous, then significant differences may appear [START_REF] Termentzidis | Cross-plane thermal conductivity of superlattices with rough interfaces using equilibrium and non-equilibrium molecular dynamics[END_REF]. This is because the EMD method is based on the Green Kubo formula, which is only valid for homogeneous systems. In this study the NEMD method, detailed in section II.C.1., is used to determine the thermal conductivity of silica aerogels.

The duration of NEMD simulations (number of time steps) and thus the computation time is governed by the characteristic time, , to reach the steady state for heat transfer. To evaluate , we consider a solid system made of a continuous material. The system is defined by its volume, , which a characteristic length, , a specific heat, , a density, and its thermal conductivity, . If the material is opaque (no radiative heat transfer) and if is much larger than the phonon mean free path of the material , then the Fourier's law applies to describe heat transfer within this system. If the system, initially at a temperature , is placed in a heat source at temperature then the temperature of the system will tend towards the temperature of the medium. The duration of the transient state is proportional to the characteristic time = 2 / .

For the silica aerogel systems studied here, the density varies between 2255

. 3 (dense silica) and 250 . 3 (silica aerogel which porosity is about 90%), and the characteristic length between 20 and 90

. At ambient temperature, the specific heat is equal to 720 .

1 . 1 . However, Molecular Dynamics (MD) simulations are based on the Newton equation and there is no quantification of energy, so whatever the temperature level, the specific heat of each atom is 3 (with the Boltzmann constant), which is equivalent to a specific heat of 1240 .

1 . 1 . The thermal conductivity is approximately 1.5 1 . 1 for dense silica and we expect that it decreases for silica aerogels down to values of the order of 0.05 . 1 . 1 [START_REF] Coquil | Molecular dynamics study of the thermal conductivity of amorphous nanoporous silica[END_REF][START_REF] Ng | A molecular dynamics study of the thermal conductivity of nanoporous silica aerogel, obtained through negative pressure rupturing[END_REF][START_REF] Yeo | Enhanced thermal characterization of silica aerogels through molecular dynamics simulation[END_REF][START_REF] Huang | Thermal conductivity of amorphous and crystalline thin films by molecular dynamics simulation[END_REF]. Using these orders of magnitude, we may expect simulation durations should be performed within 10 to 100 . This would lead to prohibitive computation times that would prevent a large number of systems of being tested especially those close to 100 3 in size.

A new simulation strategy to optimize the computational time and extract the value of the thermal conductivity from NEMD simulations with large volumes of nanoporous silica has been developed and is described in section II. This method is used to study the minimum size of the volume of silica aerogel that has to be considered to predict the thermal conductivity of silica aerogel which density is 250

. 3 . It is also used to assess the relationship between the silica aerogel density and its thermal conductivity. Results are compared to previous theoretical predictions and experimental values in section III.

Models and methods

Simulation strategy

Considering that silica aerogel has an amorphous structure, the phonon mean free path is of the same order of magnitude as the atomic distances and heat transfer through silica is diffused. Thus, during the transient state, the temperature profile in the system obtained with NEMD simulations should be the same than in an equivalent system made of a continuous material of the same geometry and thermal properties than the system described at the atomic system. To optimize the computational time, a three steps simulation strategy to calculate the thermal conductivity of silica aerogels has been developed here: (i) The evolution of temperature profiles during the transient state is calculated by NEMD from silica aerogel simulation boxes. (ii) The evolution of the temperature profile is also calculated by solving the one dimension second order heat transfer differential equation (HTDE) based on a guess of thermal conductivity. (iii) The evolutions of these two temperature profiles are then compared in order to extract a new value for thermal conductivity by identification using a least mean square algorithm. Steps (ii) and (iii) are thus repeated until a satisfactory convergence for thermal conductivity is achieved. Details on these different steps are given below.

(i) We start from cubic simulation boxes of nanoporous silica of volume 3 with a prescribed density and equilibrated at 300 . A given energy quantity is removed from the cold zone and given to the hot zone at each time step using a velocity rescaled technic, see Fig. 1a. As periodic conditions are used, to avoid direct heat transfer from the hot to the cold zones, atoms have been fixed between these two zones. Thus, heat transfer may only occur between the two zones through the silica aerogel. The temperature profiles , are stored in regularly spaced slices at positions (the middle of each slice ) of 4Å thickness and at different times . The value of the energy quantity and the number of time steps has to be adapted for each system size and density to ensure both a large enough temperature gradient compared to temperature fluctuations and a good convergence.

(ii) The temperature profiles , of a 1D system of length initially at temperature , 0 = 300 (see Fig. 1b) is calculated at the same locations and times than for the NEMD simulation, solving the second order heat transfer differential equation (HTDE): , = ,

To reproduce the same condition of heat transfer than for the MD simulations, the heat flux = !/ " ∆ , with $ the time step used in MD simulation, is imposed on both sides of the system:

= %0, "&, ∀ : , = )
At the first time this step, a guess of the thermal conductivity is used ranging between 1.5 We observed that the temperature field calculated with MD simulations is submitted to fluctuations,

8

~: / with the number of atoms and the temperature level. This is due to the fact that the number of atoms is relatively small. These fluctuations inhibit the functional to converge towards zero. The variation of the thermal conductivity between two iterations should decrease towards zero if there were no temperature fluctuations. Due to temperature fluctuations, when the iteration number is large enough, the thermal conductivity oscillates around a mean value with a relative variation of less than 1% which is the chosen convergence criterion. The quality of the results may be checked by looking at the mean deviation ; between the temperature calculated using the HTDE and the one obtained from MD simulation. The mean deviation is divided by the mean value of the temperature fluctuation 8 1 for the interpretation:

; = 1 8 1 . / + +< -. , / 0 -1 -. , / 0< 3 4 /56 3 7
.56

where and are the number of slices and the number of consecutive time steps used for the summation average at the time , respectively. After convergence, should be almost equal unity. Values of much smaller than unity means that the numerical solution of the HTDE gives almost the same results than MD simulations despite the temperature fluctuations, which is not physical. This might happen if the number of times used for the comparison is too small. Values of much larger than unity means that the numerical solution of HTDE differs significantly from MD simulations. In that case, the two models might differ physically. In step (i) and (ii), the initial temperature is chosen equal to 300 . However, as already stated, MD only solves the classical equation of Newton and no quantification of energy is accounted for. So, the comparison of MD results for the thermal conductivity is theoretically reliable (if the interatomic potential is accurate) for temperatures higher than the Debye temperature. In our case, our aim is to compare thermal conductivity values for different system sizes and densities. Qualitatively, the choice of the temperature range is not important. Choosing a higher temperature range would have resulted in higher temperature fluctuations. The ambient temperature offers a good compromise between the temperature gradient and temperature fluctuations.

Interatomic potential and boxes preparation

In this study, all MD simulations are performed using the Wolf BKS potential (van Beest, Kramer and van Santen), proposed by Carré et al. [START_REF] Carré | Amorphous silica modeled with truncated and screened Coulomb interactions: A molecular dynamics simulation study[END_REF][START_REF] Shcheblanov | Detailed analysis of plastic shear in the Raman spectra of SiO 2 glass[END_REF], which mainly differs from the original BKS potential by the introduction of a cut-off on Coulombic long-range interactions (shifted Wolf method [START_REF] Wolf | Exact method for the simulation of Coulombic systems by spherically truncated, pairwise r-1 summation[END_REF]) and ensuring the continuity of both the potential and the force at the cutoff radius as proposed by Fennell et al. [START_REF] Fennell | Is the Ewald summation still necessary? Pairwise alternatives to the accepted standard for long-range electrostatics[END_REF]. Fennell et al. [START_REF] Fennell | Is the Ewald summation still necessary? Pairwise alternatives to the accepted standard for long-range electrostatics[END_REF] have shown very good agreement between this method and the full Ewald summation used for the original BKS potential. In a previous study [START_REF] Gonçalves | Molecular dynamics simulations of silica surface properties with truncated Coulomb interactions[END_REF], we confirmed the transferability of the potential of Carré et al. [START_REF] Carré | Amorphous silica modeled with truncated and screened Coulomb interactions: A molecular dynamics simulation study[END_REF] for reproducing silica surface properties. The gain in computation time obtained by the truncation of Coulombic interactions opens the possibility to generate sufficiently large volumes of silica aerogel to account for the pore size distribution observed experimentally [START_REF] Gonçalves | Elasticity and strength of silica aerogels: A molecular dynamics study on large volumes[END_REF][START_REF] Gonçalves | Nanocompression of secondary particles of silica aerogel[END_REF]. Further, Carré et al. [START_REF] Carré | Amorphous silica modeled with truncated and screened Coulomb interactions: A molecular dynamics simulation study[END_REF] showed that the vibrational density of state is well reproduced in comparison with the original BKS. Moreover, considering the large surface to volume ratio of highly porous materials such as silica aerogels, it is critical that the interatomic potential be able to reproduce the structure and energy of amorphous silica surfaces. Both truncated and original BKS potentials have been extensively used to study the surface properties of amorphous silica [START_REF] Gonçalves | Molecular dynamics simulations of silica surface properties with truncated Coulomb interactions[END_REF][START_REF] Rarivomanantsoa | Classical molecular dynamics simulations of amorphous silica surfaces[END_REF][START_REF] Roder | Structure and dynamics of amorphous silica surfaces[END_REF]. Those studies suggest that the BKS potential is a reliable choice to investigate thermal properties of highly porous materials.

The Wolf BKS potential, proposed by Carré et al. [START_REF] Carré | Amorphous silica modeled with truncated and screened Coulomb interactions: A molecular dynamics simulation study[END_REF], writes:

Φ >? @ ABC D = E > E ? F G H D I H D + KL >? F M N OP - Q >? D R -SL >? F M T,UV N OP - Q >? D W, X R YZ I X D
where

G H D = S 1 D - 1 D W,H Y + D -D W,H D W,H , I H D = F [ \ ] -M M T,\ 0 ] I X D = F [ UV ] -M M T,UV 0 ]
^, _ stand for Silicon or Oxygen and F is the elementary charge. Effective charge values, E > , E ? , and parameters, L >? , >? and Q >? can be found in [START_REF] Van Beest | Force fields for silicas and aluminophosphates based on ab initio calculations[END_REF]. The expression G H D introduces a finite distance cut-off D W,H on the coulombic interactions and it ensures the continuity of both the potential and the forces at D = D W,H [START_REF] Wolf | Exact method for the simulation of Coulombic systems by spherically truncated, pairwise r-1 summation[END_REF][START_REF] Fennell | Is the Ewald summation still necessary? Pairwise alternatives to the accepted standard for long-range electrostatics[END_REF]. I H D and I X D functions are introduced to smooth interactions at the cut-off distance D W, X and D W,H . Here, D W,H = 10.17L ` and D W, X = 5.5L ` and the widths of the smoothing function are a H = a X = 0.5L ` [START_REF] Shcheblanov | Detailed analysis of plastic shear in the Raman spectra of SiO 2 glass[END_REF]. The original BKS potential, which pertains to the Coulomb-Buckingham type, does not tend to positive infinity as D tends to zero. Thus, atoms may approach too close from each other if the temperature or the pressure of the system are too high. This may induce chaotic dynamics. A solution proposed by Shcheblanov et al. [START_REF] Shcheblanov | Detailed analysis of plastic shear in the Raman spectra of SiO 2 glass[END_REF], and adopted here, is to add a strong repulsive term for short range interactions: where ; >? , ! >? and * >? parameters are listed in [START_REF] Shcheblanov | Detailed analysis of plastic shear in the Raman spectra of SiO 2 glass[END_REF].

Simulations are performed with the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) [START_REF] Plimpton | Fast Parallel Algorithms for Short-Range Molecular Dynamics[END_REF]. The integration of the equations of motion is achieved using a velocity-Verlet algorithm and a 0.25 timestep needed for stability. Temperature and pressure are respectively maintained by a Langevin thermostat and a Berendsen barostat.

Dense amorphous silica samples are prepared from melting of β-cristobalite cubic structures of 52728 atoms at 5000 K. The resulting silica liquid is then quenched to 300 K at 4.7 × 10 12 . 1 . Simulations boxes are equilibrated during 200 in the NPT ensemble at 300 K and zero stress. This box is then replicated in the three directions of space in order to get the required system size. It is then equilibrated during 1 in the NPT ensemble at 300 K and zero stress. This procedure is similar to the "melt-quench-duplicate" procedure described in [START_REF] Zhong | Deformation behavior of metallic glasses with shear band like atomic structure: a molecular dynamics study[END_REF][START_REF] Zhong | On the critical thickness for non-localized to localized plastic flow transition in metallic glasses: A molecular dynamics study[END_REF][START_REF] Zhou | Non-localized deformation in metallic alloys with amorphous structure[END_REF] and used to create "big samples" of metallic glasses. Periodic boundary conditions are used in all directions. The final amorphous silica box has a density of 2255 ± 1

. 3 .

Silica aerogel simulations boxes are generated with Kieffer's method [START_REF] Kieffer | Generation of fractal structures by negative pressure rupturing of SiO 2 glass[END_REF] by isostatically and instantaneously stretching by approximately 10% by steps a dense amorphous silica box in all directions.

Between each stretch, a relaxation stage in the NPT ensemble is applied during 50 i to maintain a 300 K temperature and a zero pressure. The stretching/relaxation stages are repeated until the desired density is reached. More details can be found in Gonçalves et al. [START_REF] Gonçalves | Elasticity and strength of silica aerogels: A molecular dynamics study on large volumes[END_REF]. With this method and the truncated BKS interatomic potential, the pore size distribution is around 10

for a = 250 .

. This pore size is comparable to the experimental values [START_REF] Roiban | Advanced three dimensional characterization of silicabased ultraporous materials[END_REF][START_REF] Wang | Effect of preparation parameters on the texture of SiO2 aerogels[END_REF][START_REF]Aerogels handbook[END_REF] as for commercial silica aerogel obtained from sol-gel process [START_REF] Chal | Nanostructured silica used in super-insulation materials (SIM), hygrothermal ageing followed by sorption characterizations[END_REF] and leading to monolith like aerogel particles or composite made of aerogels.

Results

Amorphous dense silica

The thermal conductivity of dense amorphous silica is calculated by two methods: first by using classical NEMD simulations and second by relying on the simulation strategy described in section II.A. The aim is to validate our simulation strategy for a simple and wellknown case of dense amorphous silica. Formally, in the classical NEMD simulation, the average temperature gradient / k of the linear temperature profile is determined for the steady state conditions. The thermal conductivity is calculated from the Fourier's law:

= - ) / k
As the steady state has to be reached, the system size should not be too large to minimize the computation time. However, the value of the thermal conductivity is not reliable if the system size is too small. Thus, several simulations for different system sizes are run. The inverse of the thermal conductivity is plotted against the inverse of the system size to extrapolate the value of the thermal conductivity for an infinite system size. In our case, the optimum section of the system is 14.1 × 14.1 nm² and the system length varied between 3.5 to 28.3 nm.

The temperature level is imposed in the hot and cold source using a temperature rescale method, allowing the calculation of the heat flux exchanged between the sources. Figure 2 shows the temperature profile in the longest system. The duration of the simulation to reach the steady state is 1.3 ns. In steady state, the heat flux flowing between the hot and cold sources is equal to 3.3 .2 . The temperature profile between the sources is approximatively linear which means that the phonon mean free path is much lower than the system size, as it should since the structure is amorphous. Figure 3 shows the inverse of the thermal conductivity against the inverse of the system size. The thermal conductivity increases when the system length increases. Compared to previous calculations [START_REF] Coquil | Molecular dynamics study of the thermal conductivity of amorphous nanoporous silica[END_REF][START_REF] Ng | A molecular dynamics study of the thermal conductivity of nanoporous silica aerogel, obtained through negative pressure rupturing[END_REF][START_REF] Yeo | Enhanced thermal characterization of silica aerogels through molecular dynamics simulation[END_REF] with the original BKS potential, our thermal conductivity values are larger. This might be due to the choice of the system section of the system which is larger than previous works. Indeed, the thermal conductivity also varies when the section of the system increases. Thus, we increased the section of our systems such that the value of the thermal conductivity does not change significantly. The thermal conductivity for an infinite length is obtained from the extrapolation of the linear fitting of the points towards 1/" = 0; it is equal to 2.4 ± 0.1 . 1 . 1 . Despite our values for each system length being larger than previous results, the extrapolated value for the infinite system length is of the same order as the one obtained with the original BKS interatomic potential. All these values are almost twice as large as the experimental values (1.4 . 1 . 1 ). Indeed, Yeo et al. [START_REF] Yeo | Enhanced thermal characterization of silica aerogels through molecular dynamics simulation[END_REF] also calculated the thermal conductivity of dense silica aerogel using the Tersoff interatomic potential proposed by Munetoh et al. [START_REF] Munetoh | Interatomic potential for Si-O systems using Tersoff parameterization[END_REF]. They obtained a value of the thermal conductivity equal to 1.19 . 1 . 1 . Actually, the choice of an interatomic potential to study thermal properties is not that straight. First MD simulations are classical simulations that do not account for the energy quantification of the phonon. Second, the interatomic potential should reproduce the harmonic (dispersion curves) and anharmonic behavior (phonon relaxation time) of the material. So, whatever the interatomic potential, MD simulations may result in over or under-estimation of the thermal conductivity for temperatures lower than the Debye temperature. So, as already discussed for the choice of the interatomic potential, the optimization of computational time and the ability of the BKS interatomic potential to reproduce the surface and structural properties of silica aerogel confirm our interatomic potential choice to investigate silica aerogel thermal properties. Since the bulk thermal conductivity of silica aerogel obtained with the BKS interatomic potential is not correct, only relative variations of the thermal conductivity of silica aerogel from the bulk value will be considered.

The new methodology for the thermal conductivity prediction is used for comparison with these first results. The system is a dense amorphous silica prepared as mention in section II.B. To test the ability of the method to tackle with large systems, the volume of amorphous silica is set to 37

, containing approximately 3.3 million of atoms. The initial temperature is equal to 300 K. The system is submitted to a heat flux of 1.18 IG.

during 0.78 on the opposite sides of the sample in the x direction. With this value the temperature variation between the hot and cold sources is around 18 K, which is one order of magnitude larger than the temperature fluctuations (0.43 K). Figure 4 shows the variation of the thermal conductivity against the simulation duration of the heat transfer. The thermal conductivity tends to a limit value equal to 2.35 . 1 . 1 when the simulation duration increases. This value is very close to the value obtained from the classical NEMD method. Figure 5 compares the temperature profiles obtained from the MD and HTDE simulations near the hot and cold sources. The quality of the diffusive model in dense silica is good: we find that the mean deviation ; is almost equal to 1. This validates the new methodology to determine the thermal conductivity from MD simulation for amorphous material. It is thus used for silica aerogel in the remaining of the paper.

As explained in section 2.2, all big dense amorphous silica systems are constructed by replications of a small dense amorphous silica configuration which are then relaxed during 1 at 300 K and zero stress before heat transfer simulation. After replications, the samples can be annealed at high temperature in order to eliminate possible artifacts of periodicity from replication [START_REF] Cao | Structural processes that initiate shear localization in metallic glass[END_REF]. In order to prove that our procedure of replication without annealing have no influence on the thermal conductivity, the dense amorphous silica system of dimension 14.1 × 14.1 × 28.3

obtained by replication has been annealed at 1000 o for 0.5 then quenched to 300 K at 4.7 × 10 6 o. 6 and relaxed for 0.5 . The thermal conductivity of this annealed system obtained by NEMD simulation is 2.31 G.

6 o 6 which is very close to the one previously obtained without annealing 2.32 G.

6 o 6 . This confirms the reliability of our simulation strategy.

Silica aerogel

In previous studies, due to the used of the original BKS interatomic potential, the maximum system size that was considered was 20 nm, which is of the order of magnitude of the pore characteristic size in silica aerogels. Moreover, the minimum value of the density (300 / 3 in Yeo et al. [START_REF] Yeo | Enhanced thermal characterization of silica aerogels through molecular dynamics simulation[END_REF]) is larger than the one investigate here (250 / 3 ). Thanks to the use of the truncated BKS interatomic potential and the new simulation strategy to determine the thermal conductivity larger systems may be considered. Our first aim is to study the influence of the system size on the thermal conductivity of silica aerogel to look for the system size that would lead to the thermal conductivity of the bulk silica aerogel. Thus, using the Kieffer's method [START_REF] Kieffer | Generation of fractal structures by negative pressure rupturing of SiO 2 glass[END_REF], seven aerogel volumes, all with a density equal to = 250

. 3 , with different volumes are generated, = 30 3 , 40 3 , 50 3 , 60 3 , 70 3 , 80 3 and 90 3 3 is imposed on the opposite sides of the sample. This procedure is repeated on all three directions ( , s or k) in order to induce a 1D temperature field variation.

Results are detailed for the largest system 90 3 3 first. As for the dense silica, the thermal conductivity has been identified as a function of the simulation duration (figure 6). The thermal conductivity also tends towards a limit value after 2.6 , a time smaller than the evaluated characteristic time to reach the permanent regime larger than 10 . However, in this case, the mean temperature deviation ; is not constant and increases when the simulation duration increases (figure 6) and becomes larger than the temperature fluctuation of the MD simulation (0.53 K). The comparison between temperatures obtained from MD and HTDE simulations on figure 7 suggests that the two models (MD and HTDE) differ slightly. This disagreement might be explained by the small variation of the potential energy of the system during MD simulation due to the relaxation of the system. This effect cannot be accurately described in the HTDE but actually has a negligible influence on the value of the thermal conductivity (see discussion in appendix). Therefore, for all the systems considered hereafter, the value of the thermal conductivity is determined with the same methodology than for the dense silica. All the values of the thermal conductivities obtained for all the system sizes and directions of heat transfer are gathered on figure 8. For smaller systems, the values are quite different in the x, y and z directions while for large systems, this difference decreases and almost vanishes. As the characteristic size of the pores is 10 then the anisotropy of the thermal conductivity for the smaller systems is due to the anisotropy of the geometry of the system. This anisotropy originates from the pore coalescence during the negative pressure cycles applied for the aerogel's generation process due to the very low target density. For larger systems, the geometry becomes isotropic, so does the thermal conductivity. This result is coherent with the one already shown by Gonçalves et al. [START_REF] Gonçalves | Elasticity and strength of silica aerogels: A molecular dynamics study on large volumes[END_REF] for the elastic modulus.

These predicted values of the thermal conductivity of silica aerogel might be quite different from the real one since, as discussed in the previous section, the prediction for dense amorphous silica is an overestimation of the experimental value. Moreover, our silica aerogel systems are quite different from real systems that have already been experimentally characterized. Real systems exhibit a two-scale organization: the nanoscale that has the same geometric characteristic than our silica aerogel and a larger scale due to the forming process of the final insulators. However, in previous studies, it has been shown that the macroscopic properties actually depend on the nanoscale. This is the case for the variation of the elastic modulus and thermal conductivity which depends on the system density [START_REF] Neugebauer | Thermal conductivity and characterization of compacted, granular silica aerogel[END_REF][START_REF] Nocentini | Hygro-thermal properties of silica aerogel blankets dried using microwave heating for building thermal insulation[END_REF]. Gonçalves et al. [START_REF] Gonçalves | Elasticity and strength of silica aerogels: A molecular dynamics study on large volumes[END_REF] has also shown that the elastic modulus power law dependence has the same exponent value as modeling or experimental studies [START_REF] Woignier | Mechanical Properties and Brittle Behavior of Silica Aerogels[END_REF][START_REF] Wang | Mixing and packing of binary hydrophobic silica aerogels[END_REF][START_REF] Ma | Mechanical structure-property relationship of aerogels[END_REF]. So, the aim is to predict the density variation of the thermal conductivity of silica aerogel, to be compared with the experimental one. , with respectively volume t = 90 , 74 3 , 70 3 , 68 3 , 64 3 , 60 3 , 56 3 , 53 3 , 48 3 , 42 3 and 37 3 3 , are generated. The number of atoms is almost the same for all these systems. The absolute value of the heat flux imposed on each side of the simulation boxes are respectively: p = 4.09, 3.81 3.82, 3.81, 3.90, 39.6, 40.9, 42.5, 46.2, 53.4 and 118.2 × 10 r G.

. The heat flux for systems with high density is higher in order to prevent that the temperature fluctuations are of the same order as the temperature difference imposed on each side of the simulation box. The variation of the .K -1 )

x direction y direction z direction thermal conductivity is then plotted as a function of the system density (figure 9). For each density, the values of the thermal conductivities in the , s and k directions are almost the same, indicating that for each density, the system size is large enough to be considered homogeneous. Also, the thermal conductivity clearly varies as a power law function of the density: = u

Our coefficient v = 1.77 is much larger than the one predicted by Ng and Yeo [START_REF] Ng | A molecular dynamics study of the thermal conductivity of nanoporous silica aerogel, obtained through negative pressure rupturing[END_REF] v = 1.01 , for a density range from = 320 to 990 . , but their system size was almost five time smaller than ours. However, note that the system size has to be large enough to ensure that the system is representative of bulk silica aerogel. Fricke et al. [START_REF] Hummer | Proc. 2nd Int. Symp. Aerogels 24[END_REF] derived a model for the thermal conductivity of aerogels. Based on this model and the sound velocity measurements in aerogels Gross et al. [START_REF] Gross | Ultrasonic velocity measurements in silica, carbon and organic aerogels[END_REF] and Hrubesh et al. [START_REF] Hrubesh | Thermal properties of organic and inorganic aerogels[END_REF] calculated a value of b = 1.88. Thus, our value of b is in good agreement with this prediction. This coefficient b has also already been experimentally determined by Fricke et al. [START_REF] Fricke | Optimization of monolithic silica aerogel insulants[END_REF] and Jain et al. [START_REF] Jain | Processing dependent thermal conductivity of nanoporous silica xerogel films[END_REF]. The reported experimental values of b range between 1 and 1.65. In these two papers, it has been highlighted that b strongly depends on the elaboration process.

Another empirical model has been proposed by Emmerling et al. [START_REF] Emmerling | Scaling properties and structure of aerogels[END_REF]. It is based on the fractal dimension w x of silica aerogels. As w x is around 2 for the density range of interest here (around 200

. 3 ) they determined a value b = 1.5, which seems to be in good agreement with experimental results [START_REF] Scheuerpflug | Thermal properties of silica aerogels between 1.4 and 330 K[END_REF]. However, for our silica aerogel systems, w x varies when the density varies [START_REF] Gonçalves | Elasticity and strength of silica aerogels: A molecular dynamics study on large volumes[END_REF], the value of b should also varies with the density. To summarize, the relative variation of the thermal conductivity as a function of density in our simulated silica aerogel is in good agreement with previous results. However, the comparison with experimental results is quite difficult since the interatomic potential does not allow for the real complexity of chemical interactions and chemical diversity to be reproduced.

Summary

The prediction of silica aerogel thermal conductivity is of great interest for material engineering. Due to the nanostructure of silica aerogel, the atomic description of the material is mandatory for the conduction contribution to the thermal conductivity. Molecular Dynamics simulation is an appropriate tool but the drawback of this technique is the large computational time linked to the choice of the interatomic potential and also due to the simulation strategy to predict the thermal conductivity.

In this paper, the computational time is first significantly reduced thanks to the use of the interatomic potential from Carré et al. [START_REF] Carré | Amorphous silica modeled with truncated and screened Coulomb interactions: A molecular dynamics simulation study[END_REF] which has already been evaluated for MD simulation of silica aerogel mechanical properties [START_REF] Gonçalves | Elasticity and strength of silica aerogels: A molecular dynamics study on large volumes[END_REF][START_REF] Gonçalves | Nanocompression of secondary particles of silica aerogel[END_REF]. Secondly, we took advantage of the diffusive nature of heat transfer in silica aerogel to propose a new methodology for the thermal conductivity prediction. It is based on the comparison of the transient temperature field given by the MD simulation and by the solution of the heat transfer diffusive equation. This new methodology makes it possible to reduce the computing time from above 10 with classical NEMD to 2.6 for our large systems of almost 90 3 3 .

The method has been validated by comparing the prediction of the thermal conductivity of dense silica using classical NEMD and the proposed methodology. It is then used to predict the thermal conductivity of silica aerogels. The thermal conductivity of silica aerogels has been determined as a function of the system size for a density of 250 . 3 . For such a density, the average pore size is 10 nm. The anisotropy of the thermal conductivity due to the presence of pores in the aerogel decreases when the system size increases and the system may be considered isotropic when it reaches 90 nm. This is the first time such a system size has been considered for heat transfer simulation. The thermal conductivity of silica aerogel has also been calculated as a function of the aerogel density. As expected, the thermal conductivity exhibits a power law dependence with the density. The power factor is in good agreement with the experimental results.

This new methodology to predict the thermal conductivity might be used to study materials in which conduction heat transfer is diffusive but which require the atomic description of large volume. The results are shown on figure A4 and compared to the one obtained without the source term. Accounting for the source term do not change the value of the thermal conductivity, but it decreases significantly the mean deviation ;. However, ; is still larger than 1. To our point of view, this is due to the fact that in our model, the source term is homogeneous while in MD simulation, this source term is due to local relaxation. This is however not possible to describe the spatial variation of the source term. So, there will always be a model difference between the HTDE and the MD simulation. The figure A5 confirms this assumption: the new temperature histories near the heat sources with MD simulation and the solution of HTDE with the source term exhibit a better agreement near the hot source and it seems that the agreement is less good near the cold source (see figure 7 for the same comparison without the source term) : due to the local nature of the source term it is not possible to have a perfect agreement between the two solutions in the whole system. As a conclusion, it is clear that the value of the thermal conductivity does not change if the term source is considered or not in the HTDE (several tests have been done for other system sizes and density), thus the results for the thermal conductivity are all obtained solving equation 1 of the paper without the source term.

  Using a least mean square algorithm, a new value of the thermal conductivity is calculated by minimizing the functional * = + + , -. , / 0 -1 -. ,

Fig. 1 .

 1 Fig.1. Configuration for a) the NEMD simulations and for b) the classical 1D transient simulation of heat transfer in a homogeneous continuous medium.

Fig. 2 .

 2 Fig.2. Temperature profile in steady state in a dense amorphous silica which section is 14.1×14.1 nm².

Fig. 3 .

 3 Fig.3. Inverse of the thermal conductivity function of the inverse of the system size for dense silica.

Fig. 4 .

 4 Fig.4. Identified thermal conductivity of dense silica function of the simulation duration.

Fig. 5 .

 5 Fig.5. Temperature variation in bulk silica during heat transfer simulations, comparing MD and HTDE results for thermal conductivity identification.

  . A heat flux of p = 8.34, 6.57, 4.77, 4.38, 3.82, 4.08, and 4.09 × 10 r G.

Fig. 6 .

 6 Fig.6. Thermal conductivity (point) and adimensional mean deviation ; (square) versus the simulation duration.

Fig. 7 .

 7 Fig.7. Temperature profile obtained from MD and HTDE for the (90 ) system.

Fig. 8 .

 8 Fig.8. Thermal conductivity values obtained for different system size. The density of silica aerogel is 250 . . Eleven aerogel systems of different densities = 250, 280, 320, 375, 450, 530, 660, 780, 1050, 1590 and 2255 ., with respectively volume t = 90 , 743 , 703 , 683 , 64 3 , 603 , 563 , 533 , 483 , 423 and 373 3 , are generated. The number of atoms is almost the same for all these systems. The absolute value of the heat flux imposed on each side of the simulation boxes are

Fig. 9 .

 9 Fig.9. Thermal conductivity of silica aerogel function of the density.

Fig. A1 .

 A1 Fig.A1. Mean temperature of the system obtained from the solution of the HTDE with the source term compared to MD mean temperature.

Fig. A2 .

 A2 Fig.A2. Potential energy variation of the system during MD heat transfer simulation.

Fig. A3 .

 A3 Fig.A3. Principle of the energy variation during the transient NEMD simulation.

Fig. A5 .

 A5 Fig.A5. Temperature history near the hot and cold sources. Comparison between MD results and the HTDE with source term.

Fig. A4 .

 A4 Fig.A4. Thermal conductivity and mean deviation ; as a function of the simulation duration without (points and filled squares) or with (circles and empty squares) source term.

Ng et al.[START_REF] Ng | A molecular dynamics study of the thermal conductivity of nanoporous silica aerogel, obtained through negative pressure rupturing[END_REF], S=3.58 2 nm[START_REF] Rajappa | A molecular dynamics study of ambient and high pressure phases of silica: Structure and enthalpy variation with molar volume[END_REF] Ng et al.[START_REF] Ng | A molecular dynamics study of the thermal conductivity of nanoporous silica aerogel, obtained through negative pressure rupturing[END_REF], S=2.[START_REF] Yeo | Enhanced thermal characterization of silica aerogels through molecular dynamics simulation[END_REF]

Acknowledgements

This work was performed using HPC resources from the FLMSN, "Fédération Lyonnaise de Modélisation et Sciences

Numériques", partner of EQUIPEX EQUIP@MESO. Funding for this project was provided by a grant from la Région Rhône-Alpes Auvergne (AURA).

Appendix

Model difference between HTDE and MD simulations

For dense silica, the temperature fields calculated solving the HTDE (equations 1 and 2) and Molecular Dynamics simulations match quite well (see figure 5). Thus, the identification of the thermal conductivity is reliable. For silica aerogel, there is a difference between the HTDE solution and the MD simulation: this is shown on figure 7 and moreover the mean deviation ; is larger than 1 (figure 6). Thus, one has to discuss the reason of this difference. To illustrate the discussion, the largest system (90 3 3 ) with the lowest density ( 250. 3 ) is considered. First, it is important to understand the origin of the difference between the temperature fields, since it might help to modify the model used to simulate heat transfer in the aerogel considered as a continuous medium.

The main point is that the mean temperature of the system (so its kinetic energy) during the MD varies between 300.8K and 302.1K (figure A1). This is not expected since during heat transfer simulation the same quantity of energy is given to the hot source and taken from the cold source of the system (Fig. 1a); the mean temperature of the system must remain constant during all the simulation. We also noted that during the MD simulation, the total energy is constant, which is coherent with the boundary condition for heat transfer simulation while the potential energy decreases (Figure A2) and the kinetic energy (so its temperature) increases. Thus, we think that the scenario that leads to the temperature increase is the following (Figure A3): during the simulation, the system experiences relaxations, giving rise to a decrease of potential energy by the value ! M corresponding to a virtual state of the system which correspond to the relaxed system. As the total energy of the system is constant and due to the equipartition theorem, the energy ! M due to the relaxation is given back to the system in the virtual state, half of it is transformed into kinetic energy (which is equivalent to a temperature rise) and the other part into potential energy. So, between the initial state and the final state of the simulation, the kinetic energy increases while the potential energy decreases. For the heat transfer model, this is equivalent to add a term source in the HTDE: , = , +

Using the time variation of the temperature and potential energy, it is possible to identified a source term function of time. To identify the thermal conductivity, the only change to do in the methodology is to add the source term in the HTDE.


S Cantrel

M Lebègue

Badawi

H Jabraoui

E P Hessou

S Chibani

L Cantrel

S Lebègue

M Badawi

Adsorption of volatile organic and iodine compounds over silver-exchanged mordenites A comparative periodic DFT study for several silver loadings

de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

During a core meltdown accident in a nuclear power plant, the iodine species CH3I and I2 may be released to the atmosphere, resulting in negative health effects directly related to radiation [START_REF] González | CHERNOBYL VIS-À-VIS THE NUCLEAR FUTURE: AN INTERNATIONAL PERSPECTIVE[END_REF] due to their relatively easy assimilation by thyroid. Many efforts have been made to confine these species and limit their spread in the environment. In other words, the effective capture of iodine compounds is of great importance in terms of public safety. Nowadays, a common strategy to trap gaseous iodinecontaining compounds involves the iodine precipitation by using porous solids [START_REF] Huve | Porous sorbents for the capture of radioactive iodine compounds: a review[END_REF][START_REF] Riley | Materials and processes for the effective capture and immobilization of radioiodine: A review[END_REF] such as polymers [START_REF] Rachuri | Two-Dimensional Coordination Polymers Comprising Mixed Tripodal Ligands for Selective Colorimetric Detection of Water and Iodine Capture[END_REF], nanocomposites [START_REF] Ma | Highly Efficient Iodine Capture by Layered Double Hydroxides Intercalated with Polysulfides[END_REF], metal-containing systems [START_REF] Ponce-Vargas | Tiara-like Complexes acting as Iodine Encapsulating Agents: The Role of M•••I Interactions in [M(μ-SCH2CO2Me)2]8⊂I2 (M = Ni, Pd, Pt) Inclusion Compounds[END_REF], and zeolites [START_REF] Kitamura | Silverimpregnated alumina for removal of radioactive methyl iodide[END_REF][START_REF] Li | Capture of organic iodides from nuclear waste by metal-organic framework-based molecular traps[END_REF][START_REF] Jolley | Desorption of organics from silver zeolite[END_REF][START_REF] Sava | Capture of Volatile Iodine, a Gaseous Fission Product, by Zeolitic Imidazolate Framework-8[END_REF][START_REF] Nandanwar | Capture of harmful radioactive contaminants from off-gas stream using porous solid sorbents for clean environment -A review[END_REF]. The latter exhibit appealing features such as a high removal efficiency, simple equipment design, and a low maintenance cost.

Among the most efficient porous sorbent materials, silver-exchanged zeolites present the unique ability to strongly capture iodine species by the formation of AgI precipitates [START_REF] Maeck | Highly efficient inorganic adsorber for airborne iodine species (silver zeolithe development studies)[END_REF][START_REF] Rastunov | Evaluation of the sorbent layer thickness in iodine filters[END_REF]. Indeed, Chapman and coworkers [START_REF] Chapman | Radioactive Iodine Capture in Silver-Containing Mordenites through Nanoscale Silver Iodide Formation[END_REF] employed a differential approach for the analysis of the pair distribution function (PDF) of an I2-treated Ag-MOR, and they reported the structure and distribution of AgI correlated to the iodine trapping. To explore deep this phenomenon, Azambre and Chebbi [START_REF] Azambre | Evaluation of Silver Zeolites Sorbents Toward Their Ability to Promote Stable CH3I Storage as AgI Precipitates[END_REF] studied different silver exchange zeolites (Faujasite X and Y, Mordenite, *BEA, MFI, FER) to probe the relationships between chemical, structural and textural parameters, as well as their iodine trapping efficiencies (CH3I). They observed the formation of AgI for all silver zeolites as well as the fact that the amount of AgI precipitates is directly related to the silver content. Moreover, Chapman et al [START_REF] Chapman | Radioactive Iodine Capture in Silver-Containing Mordenites through Nanoscale Silver Iodide Formation[END_REF] reported that Ag-MORs represent an important safe pathway for radioactive iodine capture since they present the ability of using pore-blockage to further trap the iodine for long-term storage. All these characteristics make Ag-MORs used for iodine trapping under mild operating conditions in nuclear fuel reprocessing facilities [START_REF] Chapman | Radioactive Iodine Capture in Silver-Containing Mordenites through Nanoscale Silver Iodide Formation[END_REF][START_REF] Bruffey | Capture of Elemental and Organic Iodine from Dilute Gas Streams by Silver-exchanged Mordenite[END_REF][START_REF] Vance | X-ray studies of iodine sorption in some silver zeolites[END_REF]. However, their disadvantages are the high cost of silver and the possible inhibiting effect of gaseous species on the adsorption of the iodine species [START_REF] Jolley | Desorption of organics from silver zeolite[END_REF][START_REF] Chibani | A DFT investigation of the adsorption of iodine compounds and water in H-, Na-, Ag-, and Cu-mordenite[END_REF][START_REF] Chibani | Impact of the Si/Al ratio on the selective capture of iodine compounds in silver-mordenite: a periodic DFT study[END_REF][START_REF] Chebbi | Evaluation of volatile iodine trapping in presence of contaminants: A periodic DFT study on cation exchangedfaujasite[END_REF].

Among the potential inhibiting species formed during nuclear severe accidents [START_REF] Nenoff | Silver-mordenite for radiologic gas capture from complex streams: Dual catalytic CH3I decomposition and I confinement[END_REF][START_REF] Soelberg | Deep Bed Iodine Sorbent Testing FY 2011 Report[END_REF], we can find some VOCs such as aromatics (benzene, 1.3-dimethylbenzene), ketone (2-propanone), alcohols (ethanol, methanol, 2-propanol), and alkane (methane, cyclohexane) as suggested by Jolley et al [START_REF] Jolley | Desorption of organics from silver zeolite[END_REF] The possible source of organic molecules are numerous in nuclear power plants such as paints containing solvents [START_REF] Wren | Dissolution of organic solvents from painted surfaces into water[END_REF], polymers (electric cables) and oils. However to our knowledge, there is no experimental data concerning the possible inhibiting effect of these molecules on the trapping of iodine on silver zeolites. Indeed, some theoretical works have been recently published [START_REF] Chibani | A DFT investigation of the adsorption of iodine compounds and water in H-, Na-, Ag-, and Cu-mordenite[END_REF][START_REF] Chibani | Impact of the Si/Al ratio on the selective capture of iodine compounds in silver-mordenite: a periodic DFT study[END_REF][START_REF] Chebbi | Evaluation of volatile iodine trapping in presence of contaminants: A periodic DFT study on cation exchangedfaujasite[END_REF], confirming that other contaminants such as CO or H2O can inhibit the adsorption of the iodine species. In this context, Chebbi and coworkers [START_REF] Chebbi | Evaluation of volatile iodine trapping in presence of contaminants: A periodic DFT study on cation exchangedfaujasite[END_REF] reported the inhibiting effect of CO on the I2 adsorption on Ag-exchanged faujasite (Ag-FAU), while Chibani evidenced the inhibiting effects of CO and H2O on I2 and CH3I trapping by Ag-MOR [START_REF] Chibani | A DFT investigation of the adsorption of iodine compounds and water in H-, Na-, Ag-, and Cu-mordenite[END_REF]. Remarkably, the inhibiting effect of the former contaminants could be considerably limited when the Si/Al ratio decreases [START_REF] Chibani | Impact of the Si/Al ratio on the selective capture of iodine compounds in silver-mordenite: a periodic DFT study[END_REF].

The objective of this work is to study if some organic contaminants that exist in nuclear installations are capable to inhibit the capture of iodine species. To this aim, we have used DFT to assess the inhibiting effects of the following VOCs: methane (CH4), cyclohexane (C6H12), benzene (C6H6), 1,3dimethylbenzene (C8H10), methanol (CH3OH), ethanol (C2H5OH), propan-2-ol (C3H7OH), and propanone (C3H6O), on the iodine trapping into Ag-MOR. Three Si/Al ratios for the Ag-MOR have been considered on the basis of experimental data: Si/Al = 5, 11, and 47 [START_REF] Nenoff | Silver-mordenite for radiologic gas capture from complex streams: Dual catalytic CH3I decomposition and I confinement[END_REF].

The paper is organized as follows: first, we detail our computational procedure, then, the obtained results with Ag-MOR for several Si/Al ratios are presented and discussed with a focus on the potential inhibiting effect of the VOCs on the trapping of iodine species, and finally we outline the main conclusions of our study.

Computational details

Structural models

The considered mordenite framework displays three different types of cavities: the main channel (MC), the side channel (SC), and the side pocket (SP). The MC is obtained from 12-membered rings stacked in a parallel array along the c-axis. The SP is accessible from the MC through 8-membered rings (more details are available in Ref 19). The MC has a large aperture to accommodate molecules of aromatic sizes [START_REF] Liu | A THEORETICAL STUDY ON ADSORPTION OF THIOPHENE, BENZENE, AND HEXANE IN Na-, Ag-, AND Cu[END_REF][START_REF] Benco | Adsorption of H2O, NH3 and C6H6 on alkali metal cations in internal surface of mordenite and in external surface of smectite: a DFT study[END_REF] as well as the adsorption of iodine molecules is favored in this channel compared to the SC and SP [START_REF] Chibani | A DFT investigation of the adsorption of iodine compounds and water in H-, Na-, Ag-, and Cu-mordenite[END_REF][START_REF] Chibani | Performance of CuII-, PbII-, and HgII-Exchanged Mordenite in the Adsorption of I2, ICH3, H2O, CO, ClCH3, and Cl2: A Density Functional Study[END_REF].

To obtain Ag-MORs, some silicon atoms are substituted by aluminum atoms to create a negative charge, allowing to add the silver cation. In the literature, the pure siliceous mordenite conventional cell contains 144 atoms (48 Si and 96 O) [START_REF] Alberti | The crystal structure refinement of a natural mordenite[END_REF], exhibiting an orthorhombic structure (space group:

Cmcm) with lattice parameters: a = 18.094 Å, b = 20.516 Å and c = 7.524 Å. Here our Ag-MORs structure are based on the structure optimized by Bucko et al [START_REF] Bučko | The role of spatial constraints and entropy in the adsorption and transformation of hydrocarbons catalyzed by zeolites[END_REF] with the primitive monoclinic cell doubled in the c-direction, its lattice parameters are estimated at a=b=13.805 Å, c = 15.212 Å and

= β
82.81°. Based on previous works [START_REF] Chibani | A DFT investigation of the adsorption of iodine compounds and water in H-, Na-, Ag-, and Cu-mordenite[END_REF][START_REF] Chibani | Impact of the Si/Al ratio on the selective capture of iodine compounds in silver-mordenite: a periodic DFT study[END_REF] which systematically analyzed the distribution of Al and Ag atoms inside Ag-MOR, we have selected the most energetically stable site of the MC, namely T1-E according to the Mortier's [START_REF] Mortier | Compilation of Extra Framework Sites in Zeolites (Structure Commission of the International Zeolite Association)[END_REF] classification. This location is the most frequently occupied by cations [START_REF] Mortier | Compilation of Extra Framework Sites in Zeolites (Structure Commission of the International Zeolite Association)[END_REF], and has been also chosen in other theoretical investigations [START_REF] Bučko | Adsorption and vibrational spectroscopy of ammonia at mordenite: Ab initio study[END_REF][START_REF] Bučko | Adsorption and Vibrational Spectroscopy of CO on Mordenite: Ab initio Density-Functional Study[END_REF][START_REF] Benco | Adsorption of H2O, NH3 and C6H6 on alkali metal cations in internal surface of mordenite and in external surface of smectite: a DFT study[END_REF].

Periodic density functional theory calculations

The adsorption of iodine and organic compounds onto Ag-zeolites is investigated using the Vienna Ab-initio Simulation Package (VASP) [START_REF] Kresse | Ab initio molecular dynamics for liquid metals[END_REF]. These calculations employed the PBE functional [START_REF] Perdew | Generalized Gradient Approximation Made Simple[END_REF] and the PAW method developed by Blöchl [START_REF] Blöchl | Projector augmented-wave method[END_REF] and later adapted by Kresse and Joubert [START_REF] Kresse | From ultrasoft pseudopotentials to the projector augmented-wave method[END_REF]. The plane wave cutoff energy has been set to 450 eV and a Gaussian smearing σ = 0.1 eV is used to help the convergence. Atomic positions are fully optimized up to all forces are less than 0.02 eV / Å per atom.

The Kohn-Sham equations were solved self-consistently until the energy difference between cycles take a value smaller than 10 -6 eV. The Brillouin-zone sampling was restricted to the Γ-point. The considered dispersion correction scheme in our calculations is the Tkatchenko-Scheffler scheme with iterative Hirshfeld partitioning [START_REF] Tkatchenko | Accurate Molecular Van Der Waals Interactions from Ground-State Electron Density and Free-Atom Reference Data[END_REF][START_REF] Bučko | Improved Density Dependent Correction for the Description of London Dispersion Forces[END_REF] as implemented in VASP by Bucko et al [START_REF] Bučko | Improved Density Dependent Correction for the Description of London Dispersion Forces[END_REF][START_REF] Bučko | Extending the applicability of the Tkatchenko-Scheffler dispersion correction via iterative Hirshfeld partitioning[END_REF] So far, this method shows promising results to describe adsorption phenomena onto zeolites [START_REF] Chibani | A DFT investigation of the adsorption of iodine compounds and water in H-, Na-, Ag-, and Cu-mordenite[END_REF][START_REF] Chibani | Impact of the Si/Al ratio on the selective capture of iodine compounds in silver-mordenite: a periodic DFT study[END_REF][START_REF] Bučko | Improved Density Dependent Correction for the Description of London Dispersion Forces[END_REF][START_REF] Bučko | Extending the applicability of the Tkatchenko-Scheffler dispersion correction via iterative Hirshfeld partitioning[END_REF].

The adsorption phenomena can be evaluated using the interaction energy of zeolite-adsorbate, which can be expressed using the combination of the following three terms:

∆E = E -E -E (1)
Where:

-Ezeolite: the energy of the clean MOR;

-EX: the energy of the isolated molecule in gaseous phase;

-Ezeolite-X: the energy of MOR with adsorbed molecules; Also, the dispersion energy contribution ∆E to the interaction energy is obtained as:

∆E = E -E -E (2)
To further improve our understanding of the nature of the bonding interactions, we have determined the difference in electronic density (Δρ) and the difference in Bader charge (ΔQ). To visualize the charge density difference (Δρ) (using VESTA [START_REF] Momma | VESTA: a three-dimensional visualization system for electronic and structural analysis[END_REF]) upon adsorption of molecules into MOR, we have calculated the charge densities of the clean MOR (ρzeolite), of the isolated molecule in the gaseous phase (ρx), and of the MOR with adsorbed molecules (ρzeolite-X.). Then to get (Δρ), we have used the following expression:

Δρ= ρzeolite-X -ρzeolite -ρX (3)
The Bader approach consists in partionning the electronic density of a system in different volumes which when integrated over space, can provide the charge of an atom [START_REF] Bader | Atoms in molecules: a quantum theory[END_REF][START_REF] Tang | A grid-based Bader analysis algorithm without lattice bias[END_REF][START_REF] Sanville | Improved grid-based algorithm for Bader charge allocation[END_REF][START_REF] Henkelman | A fast and robust algorithm for Bader decomposition of charge density[END_REF]. To get the difference in electronic charge ΔQ, we have computed the Bader charge for each atom (Z) in the zeolite-X, in the clean zeolite and in the isolated molecules.

ΔQ= Q Atom Z in zeolite-X -Q Atom Z in zeolite or X (4)
Where:

Q Atom Z in zeolite-X : the Bader charge for all the atoms in the MOR with adsorbed molecules.

Q Atom Z in zeolite or X : the Bader charge for all the atoms in the clean zeolite or in the isolated molecule in gaseous phase.

Results and discussions

Herein, we report the potential impact of certain species of VOCs on the removal of iodine released in nuclear waste. We have classified the investigated VOCs into two categories as follows:

-Hydrocarbons: benzene (C6H6), 1,3-dimethylbenzene (C8H10), cyclohexane (C6H12), and methane (CH4).

-Oxygenated hydrocarbons: methanol (CH3OH), ethanol (C2H5OH), propan-2-ol (C3H7OH), and propanone (C3H6O).

Our computed interaction energies of the iodine molecules and VOCs with Ag-MORs (Si/Al = 5, 11, and 47) are presented in Table 1. We also verified the dissociative adsorption and then confirmed that there was no dissociation observed in the studied molecules (see supporting information, the case of alcohol molecules).

Iodine species

Regarding the interactions between iodine species and Ag-MOR as a function of Si/Al ratio, we found that they increase as the Si/Al ratio decreases, especially in the case of I2. Indeed, as indicated in Let us now explain these trends: for the adsorption CH3I, the atom of I interacts with 2 atoms of Ag located on the same side in the main channel for Si/Al = 11 and 5 (see Fig. 1), which leads to similar interaction energy for both Si/Al ratio. The I-CH3 distance does not evolve upon the adsorption of CH3I whatever the Si/Al ratio, suggesting that the dissociation of iodomethane would be more difficult than the one of molecular iodine under these conditions. Indeed, in the case of I2, we can find that its atoms can interact with more silver atoms (4 Ag atoms for Si/Al = 11 and at least 6 Ag atoms for Si/Al = 5 see Fig. 1), which leads to a significant increase in its total interaction energy (50 kJ/mol for Si/Al =11 and 120 kJ/mol for Si/Al =5) compared to its adsorption into Ag-MOR for Si/Al = 47. Besides, the long-range dispersion contributions of both I2 and CH3I in to their total interaction energies present almost the same values. So, the difference can come from the DFT contributions to the overall interaction energies: the corresponding values for I2 are always higher than those of CH3I irrespective of the Si/Al ratio.

The increase of the total interaction energy with the decrease of the Si/Al ratio can be understood by studying the electron density (Δρ) and Bader charge (ΔQ) differences induced by the adsorption. In this vein, we found that the interactions between I and Ag atoms present a significant electron transfer, where a high Ag content (Si/Al = 5) results in a large charge transfer between the Ag atoms and the iodine atoms of I2 (See Fig. 1). For the ratio 5, our calculated Bader charges (see, Table SI-2), show that the iodine atoms bear the largest positive charge (+0.15e) and the Ag atoms bonded directly to the iodine bear a negative charge. In other word, iodine atoms pick up electrons during the adsorption process. For all the rest of the atoms, the Bader charge are almost unaffected by the adsorption of I2.

However, upon adsorption of CH3I, the iodine atoms always carry a negative charge of -0.24e (see, Table SI-2) and -2.21e for Ag-MOR (Si/Al=5) and Ag-MOR (Si/Al=47), respectively. Overall, these findings indicate that the I-I bond in I2 becomes weak upon adsorption, which could promote the formation of AgI precipitates as explained by Bucko et al [START_REF] Bučko | Dissociative iodomethane adsorption on Ag-MOR and the formation of AgI clusters: an ab initio molecular dynamics study[END_REF].

Hydrocarbons

Aromatic molecules

We found that the interaction energies of aromatic compounds increase when the Si/Al ratio decreases.

Indeed, the total interaction energies of benzene/ show that about half of the interaction between benzene or 1,3-dimethylbenzene and Ag-MORs is due to dispersion forces. This contribution can be correlated with π-cation interactions, which are ubiquitous in systems involving aromatic compounds and metallic ions as has been confirmed by different theoretical [START_REF] Witte | Vacuum level alignment at organic/metal junctions: "Cushion" effect and the interface dipole[END_REF][START_REF] Bilić | Adsorption of Benzene on Copper, Silver, and Gold Surfaces[END_REF][START_REF] Chwee | Adsorption studies of C6H6 on Cu (111), Ag (111), and Au (111) within dispersion corrected density functional theory[END_REF] and experimental [START_REF] Tsai | Platinum metal surface chemistry of benzene and toluene[END_REF][START_REF] Grassian | Vibrational electron energy loss spectroscopic study of benzene, toluene, and pyridine adsorbed on palladium(111) at 180 K[END_REF] investigations. Generally, aromatic molecule adsorptions on cations are particularly interesting in view of the participation of π-electrons in the total interaction. Benzene and its derivatives can form a variety of charge-transfer complexes depending on their ability to donate electrons or to accept electrons [START_REF] Mulliken | Molecular complexes: a lecture and reprint volume[END_REF].

In this regard, Kozyra et al [START_REF] Kozyra | From Electron Density Flow Towards Activation: Benzene Interacting with Cu(I) and Ag(I) Sites in ZSM-5. DFT Modeling[END_REF] have studied the π-cation interactions between benzene and Ag sites in ZSM-5, where a large charge transfer from benzene to Ag occurs. This is in line with our own results as a significant electron transfer between benzene or 1,3-dimethylbenzene and Ag-MOR is observed for all the tested Si/Al ratios and especially with the ratio 5, see Fig. 3. For the latter ratio, benzene and 1,3-dimethylbenzene exhibit for most composed atoms positive Delta Q (Table SI-2), highlighting their important interaction with the MOR structure.

To discuss the adsorption of both aromatic molecules into Ag-MORs in term of geometric parameters and with π-cation interactions [START_REF] Munakata | Supramolecular Silver(I) Complexes with Highly Strained Polycyclic Aromatic Compounds[END_REF][START_REF] Munakata | Construction and Conductivity of W-Type Sandwich Silver(I) Polymers with Pyrene and Perylene[END_REF][START_REF] Griffith | Metal ion-aromatic complexes. XX. Preparation and molecular structure of anthracenetetrakis(silver perchlorate) monohydrate[END_REF], we have reported in Table 2 the distances between the silver atoms and the center of the aromatic rings. For Si/Al = 47, a single Ag atom interacts with the aromatics rings of benzene and 1,3-dimethylbenzene at a distance of 2.66 Å and 2.74 Å, respectively.

Here, the 1,3-dimethylbenzene molecule is more adsorbed than benzene by around 27 kJ/mol. For Si/Al = 11, the adsorption of 1,3-dimethylbenzene into Ag-MOR is achieved through the interaction between its aromatic ring and one Ag atom at a distance of 2.68 Å. In contrast, for the benzene adsorption into Ag-MOR for Si/Al = 11, there are 2 Ag atoms bonded to the aromatic ring at distances of 2.89Å and 3.17Å, which show that the adsorption of 1,3-dimethylbenzene is stronger than with benzene by around 34 kJ/mol.

In the case of Ag-MOR for Si/Al=5, 1,3-dimethylbenzene and benzene interact with two and three silver atoms, respectively. The reported distances are 2.73 Å and 2.79 Å for 1,3-dimethylbenzene and 3.09 Å, 3.13 Å and 3.27 Å for benzene (see Table 2) . In this case, the 1,3-dimethylbenzene is more adsorbed than the benzene molecule by around 43 kJ/mol. Fig. 3 is showing the electron transfers occurring upon adsorption of benzene and 1,3-dimethylbenzene in Ag-MORs (see Fig. 3).

Cyclohexane

Regarding the influence of the Si/Al ratio on cyclohexane-Ag interactions, the total interaction energy of cyclohexane is equal to -173.2 kJ/mol for Si/Al = 5 while it is -153.9 kJ/mol for Si/Al = 11, and -120.3 kJ/mol for Si/Al = 47. Here, we found that cyclohexane having no π-electrons is mainly adsorbed by dispersion forces whatever the Si/Al ratio. This result is in a good agreement with reported by Suda et al [START_REF] Suda | Interaction of benzene, cyclohexene, and cyclohexane with the surface of titanium dioxide (rutile)[END_REF], who have experimentally investigated the adsorption of benzene, cyclohexene, and cyclohexane on the surface of titanium dioxide.

As reported by Barthomeuf et al [START_REF] Barthomeuf | Adsorption of benzene and cyclohexane on faujasite-type zeolites. Part 2.-Adsorption site efficiency and zeolite field influence at high coverage[END_REF], cyclohexane and benzene have the same number of carbon atoms and similar polarizabilities, but differ in that benzene can give rise to specific interactions between its aromatic ring and cations in zeolites. In our work, we found that benzene is more adsorbed than cyclohexane (which have no π-electrons) by around 46 kJ/mol in Ag-MOR (Si/Al=5), 31 kJ/mol in Ag-MOR (Si/Al=11), and 36 kJ/mol in Ag-MOR (Si/Al=47). So, the benzene π-electrons contribute to the increase of total interaction energy. While zeolites have localized charges that render their surface heterogeneous [START_REF] Aristov | Dependence of adsorption of krypton and xenon by crystals of zeolite LiX and NaX on pressure and temperature[END_REF], it can be assumed that the cyclohexane-silver interactions are weak enough to permit the adsorption of cyclohexane on an apparently homogeneous surface.

Methane

The total interaction energy of methane is -46.6 kJ/mol for Ag-MOR (Si/Al = 5), -40.5 kJ/mol for Ag-MOR (Si/Al = 11), and -45.5 kJ/mol for Ag-MOR (Si/Al = 47). We can remark that the total interaction energy has almost the same values for all used Si/Al ratios in Ag-MOR. Also, the dispersion energy increases to -25.1 kJ/mol for Si/Al = 5 against -21.8 kJ / mol for Si/Al = 11, and -17.0 kJ/mol for Si/Al = 47. Therefore at low Ag content in Ag-MOR, the adsorption of methane is mainly directed by short range interactions (more than 62% of the total interaction energy). In contrast, the methane adsorption is mainly achieved through dispersion energy (long-range interaction) when the silver content of the zeolite is high (54%).

Oxygenated hydrocarbons

The effect of the Ag content in Ag-MOR on the oxygenated molecules adsorption is also studied. We choose as case studies the alcohols: methanol, ethanol, propan-2-ol, and propanone. The interaction energies of methanol/ethanol/propan-2-ol goes from -117.6/-146.4/-164.0 kJ/mol (for Si/Al = 5 to -106.5/-125.3/-147.7 kJ/mol for Si/Al = 47, and it increases with the carbon atoms number (See Fig. 5),

suggesting that there exists a linear function relating the interaction energy variation with respect to the carbon atoms number (nC) of the studied alcohols, which can be expressed as follows:

ΔEint(nC) = a(nC) + b (3)
Where "a" and "b" represents the contribution of dispersion and DFT energies to the total interaction energy, respectively. To combine the different contribution reported in Table 1 and in Eq.3, we can identify the different energetic terms contributing to the adsorption of each alcohol. Table 3 shows the dispersion and DFT contributions of the carbon and oxygen atoms to the adsorption of the investigated alcohols. Thus, these findings indicate that the high adsorption of alcohols is directly related to the carbon chain. Indeed, going from one alcohol to another by adding a carbon, one can estimate an increase of the interaction energy to be ~20 kJ/mol.

The comparison of the C-O distance of alcohols before and after their adsorption enable us to evaluate the occurrence of C-O bond breaking upon adsorption and therefore the possibility of formation of undesirable byproducts that can affect the adsorption of iodine species. We found that this distance in methanol/ethanol/propan-2-ol varies from to 1.44 Å (gaseous phase, for all alcohol molecules) to about1.45 Å/1.45 Å/1.48Å (adsorbed state) for every Si/Al ratio (See Table 3 and Fig. SI-2). Therefore the propanol is expected to react upon adsorption much more than methanol or ethanol.

The adsorption of propanone in Ag-MOR prefers to be via its oxygen atoms. So, as Si/Al ratio decreases in Ag-MOR, the total interaction energy raises to -158.4 kJ/mol (for Si/Al = 5) compared to -142.7 kJ/mol (for Si/Al = 11) and -146.9 kJ/mol (for Si/Al = 47). On the other hand, we found that propan-2-ol is slightly more adsorbed than propanone by around 6.1 kJ/mol for Si/Al = 5, 16.8 kJ/mol for Si/Al = 11, and 0.8 kJ/mol for Si/Al = 47.

3.4. Effect of the Si/Al on the selective trapping of iodine species towards VOCs

For Ag-MOR with different Si/Al ratios, we have compared the adsorptions of the iodine species and the investigated VOCs as reported in Figure 2 and Figure 4. For each ratio, we show the most adsorbed VOCs and the less adsorbed iodine compound (I2 for Si/Al = 47, CH3I for Si/Al = 5 and 11)

as follows:

• For Si/Al = 47, the most absorbed VOCs are 1,3-dimethylbenzene/benzene/propan-2ol/propane and they present higher total interaction energies than I2 by around ~ 50/23/15/14 kJ/mol, respectively.

• For Si/Al=11, 1,3-dimethylbenzene/benzene/cyclohexane/propan-2-ol are more adsorbed than CH3I by around ~ 68/33/2/7 kJ/mol, respectively.

• For Si/Al = 5, the most strongly adsorbed VOCs are 1,3-dimethylbenzene/benzene, showing a total interaction energy higher than that of CH3I by around ~ 37/30 kJ/mol.

In general, aromatic molecules potentially present the strongest inhibiting effects on the iodine removal compared to other VOCs. More precisely, benzene, 1,3-dimethylbenzene and propan-2-ol are expected to inhibit the adsorption of I2 and CH3I.

Moreover, we conclude that in Ag-MOR with a low Si/Al ratio, larger molecules such as I2, benzene, 1,3-dimethylbenzene and propan-2-ol have total interaction energies that are 5-35% higher than in Ag-MOR with a high Si/Al ratio. This increase can be correlated with the increase in the number of cations interacting with those molecules. This behavior is clearly confirmed by the electronic density map (see Figs 1 and3). In the opposite case, small molecules such as methane and CH3I interact with almost the same number of cations, even though the Si/Al ratio decreases, resulting in the neglected effects of increasing cations on total interaction energies of these molecules.

the adsorption of I2 and CH3I. Finally, the inhibiting effects of 1,3-dimethylbenzene and benzene on I2 trapping are also limited for a low Si/Al ratio, but these effects are still present in the trapping of CH3I.

Conclusion

Using periodic DFT calculations, we have investigated the interaction between Ag-MORs and different type of gaseous molecules which can be present in the atmosphere of a containment environment generated from a severe nuclear accident, with the aim to determine the influence of the Si/Al ratio on the adsorption selectivity of iodine compounds towards VOCs.

We have found that the adsorption of iodine species is strongly influenced by the Si/Al ratio variation in Ag-MOR. Indeed, the interaction energies of I2 and CH3I with Ag-MORs strongly increase when the silver content increases, i.e. the Si/Al ratio decreases. So, both CH3I and I2 prefer to be absorbed by Ag by their iodine I atom with at a distance I-Ag of 2.6 Å, which is accompanied by a large electron transfer between Ag atoms and I atoms. On the other hand, when the iodine molecule is absorbed in Ag-MOR for Si/Al = 5 and 11, the I-I bond is elongated. Thus, these results indicate that the I-I bond becomes weak upon adsorption and the formation of AgI precipitates could be favored.

Regarding the VOCs, the cyclohexane, benzene, 1,3-dimethylbenzene and propan-2-ol are more adsorbed than iodine species in Ag-MOR for Si/Al=47, and 11, which can play an inhibiting effect on the iodine trapping. For the considered aromatic molecules, benzene and 1,3-dimethylbenzene prefer to be absorbed on Ag + through π-electrons at a distance of Ag-aromatic center ranging from 2.65 Å to 3.16 Å. In the case of the alcohols, the investigated molecules are adsorbed via their oxygen atom at a distance of 2.2Å. This issue can be overcome by decreasing of the Si/Al ratio. However, even for a low Si/Al content in Ag-MOR, the problem of the inhibiting effects of benzene and 1,3dimethylbenzene on the trapping of CH3I is still present. Besides, methane exhibits a weaker

Fig. 1 .

 1 Fig.1. Difference in electronic density (Δρ) induced by adsorption in various Ag-MORs of I2 (top) and CH3I (bottom). The blue (yellow) zones indicate density increase (decrease).

Fig. 2 .Fig. 3 .

 23 Fig.2. Interaction energies ΔEint (kJ/mol) of the investigated hydrocarbons and iodine species adsorbed on Ag-MOR with three Si/Al ratios of 5, 11 and 47.

Fig. 4 .Fig. 5 .

 45 Fig.4. Interaction energies ΔEint (kJ/mol) of the investigated oxygenated and iodine species adsorbed in Ag-MOR with three Si/Al ratios of 5, 11 and 47.

Table 1 ,

 1

	length. Concerning Ag-MOR with Si/Al = 47, each iodine atom of I2 or of CH3I interacts with an Ag
	atom. Here CH3I is more adsorbed than I2 by around 30 kJ/mol. However, for Si/Al = 11 and 5, the
	molecule of I2 is more adsorbed than CH3I by an amount of about 30 kJ/mol and 109 kJ/mol,
	respectively.

the total interaction energies of I2/CH3I are -132.6/-152.1 kJ/mol for Si/Al = 47, -183.8/-151.8 kJ/mol for Si/Al = 11, and -297.6/-188.5 kJ/mol for Si/Al = 5, in line with an increase of the I-I bond

  51% for Si/Al = 5, 42/51% for Si/Al = 11, and 47/50% for Si/Al = 47. We note that these molecules have almost the same contribution regardless of the Si/Al ratio (about 50%). These results

1,3-dimethylbenzene are -155.8/-183.3 kJ/mol in the case of Si/Al = 47, -185.1/-219.9 kJ/mol for Si/Al = 11 and -219.4/-254.3 kJ/mol for Si/Al = 5 (See Figure 2 and Table

1

). The long-range dispersion contributions of benzene/1,3-dimethylbenzene are 41/

Table 1

 1 Computed total interaction energies ∆E [kJ/mol] of methane (CH4), cyclohexane (C6H12), benzene

	(C6H6), 1,3-dimethylbenzene (C8H10), methanol (CH3OH), ethanol (C2H5OH), propan-2-ol (C3H7OH) and
	propanone (C3H6O) with Ag-MORs. The contributions of dispersion energies to the interaction energies
				ΔEint (ΔEdis) in kJ/mol
		Si/Al ratio	5	11	47
		I2	-297.6 (-56.0)	-183.8 (-42.9)	-132.6 (-46.5)
	Iodine species	CH3I	-188.5 (-50.9) -151.8 (-41.4)	-152.1 (-36.2)
		C8H10	-262.1 (-34.9) -219.9 (-114.1)	-183.3 (-92.2)
	Hydrocarbons	C6H6	-219.4 (-92.1) -185.1 (-78.9)	-155.8 (-67.7)
		C6H12	-173.2(-111.3) -153.9 (-97.3)	-120.3 (-81.8)
		CH4	-46.6 (-25.1)	-40.5 (-21.8)	-45.7 (-17.0)

∆Edis[kJ/mol] are reported in this table in parentheses. We note that the site used to adsorb these molecules is found to be in the main channel of mordenite structure.

Table 2

 2 Distances between Ag cations and aromatic centers of benzene and 1,3-dimethylbenzene in Ag-MOR

	Molecules	Si/Al ratio	π-cation bond (Å)
		5	3.09 ; 3.13 ; 3.27 ; 5.25 ; 5.26 ; 5.42 ; 5.68 ; 5.73
	Benzene	11	2.89 ; 3.17 ; 5.13 ; 5.36
		47	2.66
		5	2.73 ; 2.79 ; 3.79 ; 4.48 ; 4.91 ; 6.00 ; 6.21 ; 6.40
	1,3-dimethylbenzene	11	2.68 ; 3.87 ; 4.227 ; 4.84
		47	2.74

Table 3

 3 DFT and dispersion contributions of oxygen and carbon atoms to the alcohol adsorptions into Ag-

	MORs.
		Si/Al	CH3OH	C2H5OH	C3H7OH
	ΔEint (kJ/mol)	5	-117.6	-146.4	-164.5
		11	-111,0	-133.7	-159.5
		47	-106.5	-125.3	-147.7
	ΔEdisp (kJ/mol)	5	-32.3	-51.4	-74.7
		11	-29.3	-45.8	-64.4
		47	-24.4	-39.5	-53.6
	ΔEDFT kJ/mol	5	-85.3	--95.0	-89.8
		11	-81.7	-87.9	-95.1
		47	-82.1	-85.8	-94.1
	ΔEdisp (CHx) in	5	-23.4	-46.9	-70.3
	kJ/mol	11	-24.3	-48.5	-72.8
		47	-20.6	-41.2	-61.9
	ΔEdisp (OH) in	5	-8.9	-4.5	-4.4
	kJ/mol	11	-5.0	+2.7	+8.4
		47	-3.8	+1.7	+8.3

Table 4 C

 4 -O distance (Å) in the gas state and the adsorbed state of methanol, ethanol and propan-2-ol.

			Adsorbed state in Ag-MOR (Si/Al ratio)
		Gaseous state
			Ag-MOR (5)	Ag-MOR (11)	Ag-MOR (47)
	CH3OH	1.43	1.45	1.45	1.45
	C2H5OH	1.44	1.45	1.46	1.46
	C3H7OH	1.44	1.49	1.48	1.48

Acknowledgments

This work was performed in the frame of the French research program ANR-11-RSNR-0013-01 called MiRE (Mitigation of Iodine Releases in the Environment). We thank The PMMS (Pôle Messin de Modélisation et de Simulation) and GENCI-CCRT/CINES (Grant No. x2018-A0030810169) for providing us computer time.

inhibiting effect on iodine trapping.

Our findings show that benzene, 1,3-dimethylbenzene and propan-2-ol (in a lesser extent) are the VOCs present in the containment atmosphere which can mostly affect the capture of iodine species.

We have highlighted that Ag-MORs with a low Si/Al ratio should be used to ensure the trapping of iodine in presence of VOCs. Remarkably, at low Si/Al ratio, the complete adsorption of I2 and CH3I could be achieved in the presence of alcohols, whereas only the adsorption of CH3I is affected by the presence of aromatic hydrocarbons.

Associated content

Supporting Information. Contains configurations of alcohols adsorbed in Ag-MOR (Si/Al=5=) and computed Bader charges for several systems.