Intermolecular Force Field Parameters Optimization for Computer Simulations of CH4 in ZIF-8

The differential evolution (DE) algorithm is applied for obtaining the optimized intermolecular interaction parameters between CH 4 and 2-methylimidazolate ([C 4 N 2 H 5 ]) using quantum binding energies of CH 4 -[C 4 N 2 H 5 ] complexes.The initial parameters and their upper/lower bounds are obtained from the general AMBER force field. The DE optimized and the AMBER parameters are then used in the molecular dynamics (MD) simulations of CH 4 molecules in the frameworks of ZIF-8. The results show that the DE parameters are better for representing the quantum interaction energies than the AMBER parameters. The dynamical and structural behaviors obtained fromMD simulations with both sets of parameters are also of notable differences.


Introduction
Methane (CH 4 ) is one of the greenhouse gases caused by farming, livestock, burning biofuel, and burning fuels like coal, oil, and natural gas.In order to reduce this gas in the atmosphere, efficient and cost saving methods are studied worldwide.The application of porous material for high capacity and cost-effective gas separation and storage is one intensively investigated way [1][2][3][4][5][6][7].Zeolitic imidazolate frameworks (ZIFs), a subclass of metal organic frameworks (MOFs), which are composed of positive transition metal ions such as Zn, Co, and Cu are linked by ditopic imidazolate ligands.Because of their high porosity, thermal and chemical stability, and the exact structure, they are used in the greenhouse gases separation and storage [1][2][3][4][5][8][9][10][11].ZIF-8 [8][9][10]12] is one of the ZIFs, which is composed of Zn 2+ ion linked with ditopic 2-methylimidazolate ([C 4 The knowledge of chemical processes such as diffusion, adsorption, and reaction at the molecular level can help to design and improve the properties of these materials for efficient use in gas separation and storage [1-5, 8, 10, 11].The key to understanding the behavior of the adsorption and diffusion of CH 4 in ZIF-8 is the intermolecular interaction.Since the quantum method needs large resources and time consumption, it is not suitable for computer simulation which is composed of several atoms.Therefore, the development of intermolecular interactions functions is important.These functions are generated by the experimental data or the computational quantum data.Moreover, the optimization method is used to find the suitable forms and parameters of such functions.The famous form of the intermolecular interactions is the Lennard-Jones function [5,8].In this paper, the parameters of Lennard-Jones are optimized based on ab initio using the heuristic method called the differential evolution (DE) in order to obtain the specific parameters for use in the molecular dynamic simulations of CH 4 in the ZIF-8.

Models and Calculations
The geometrical structure of [C 4 N 2 H 5 ] − is cut directly from the framework of single-crystal XRD data [9] and the geometric model of CH 4 molecule is taken from the interaction   2.
The interaction energy with the corrected basis sets superposition error Δ CC is represented by the following equation [15]: where   (),   (), and   () are the total energy of where   = 4   12  ,   = 4   6  ,   = (  +   )/2, and   = √     . and  used for computing the interactions are shown in Table 1.
Although the general AMBER force field is famous for use in molecular dynamics, it is not suitable for some specific models.Thus, we use differential evolution algorithm [16] with the interaction energies from quantum mechanics method with zero charges to optimize the specific intermolecular interaction parameters  and .The objective function is the following chi-square equation: where  is the index of configuration,  CC () is the counterpoise corrected interaction energy, () is the Lennard-Jones energy, () = 1/| CC () −   | is the weighted value with respect to   which is set very close to the minimum energy of  CC ,  =  −  is the degree of freedom,  is the total number of energies, and  is the total number of parameters.The optimized molecular interaction parameters and the original AMBER [8] parameters are used in computer simulations for studying the adsorption and diffusion of CH 4 in porous material ZIF-8 by using DL POLY program (version 2.20) [17,18].The flexible ZIF-8 framework is composed of eight (2 × 2 × 2) unit cells and loadings of CH 4 molecules per unit cell are 1, 2, 4, 6, and 8, respectively.The NVT ensemble is set at 300 K and time step of 1.0 fs.The simulations are equilibrated for 1.0 ns (1,000,000 steps), and further 1.0 ns trajectories data are collected at every 200 time steps for studying the structural and dynamical properties.
The site-site radial distribution functions (RDFs),   (), are a statistical analysis which is the probability of finding an atom  around a reference atom : where () is the number of  atoms in thickness of radius of circle  at distance ,  is the volume of box, and   is the number of all  atoms.A smaller particle will diffuse without direction.The selfdiffusion coefficient,   (m 2 /s), is computed by mean squared  displacement (MSD) method using the Einstein relation for a three-dimensional system [5,19]: where  is the total number of guest particles and   () and   ( 0 ) are the position vectors of the diffusing molecule in framework at time  and time origin  0 .

Results and Discussion
The intermolecular interaction parameters   and   between atom of CH 4 and atom of [C 4 N 2 H 5 ] − obtained from AMBER force field and DE method are presented in Table 2.
It is clear that the optimized parameters provide the value of chi-square of 0.159 which is much smaller than 35.052 obtained with AMBER parameters.Thus, the energies with DE parameters are close to the QM energies more than those obtained with the AMBER parameters (see Figure 3).
The adsorption and diffusion of CH 4 molecules in ZIF-8 framework are investigated by molecular dynamic simulations using the DE and AMBER parameters for comparison.The site-site RDFs between all atoms in framework to atoms H and C of CH 4 molecules are presented in Figure 4.The results of RDFs obtained by all simulations indicate that atoms H4 and HT are favored sites for adsorption.The parameters of DE yield the average minimum distances from atoms H and C of CH 4 molecules to atoms of the framework as 1.45 Å and 2.15 Å, which are shorter than those of 1.77 Å and 2.33 Å obtained with the original AMBER parameters.These obtained distances agree well with the interaction profiles (see Figure 3) where both minimum binding energies and minimum binding distances of QM and DE profiles are stronger and shorter than the AMBER profiles.The simulations with AMBER parameters show that CH 4 molecules are distributed within the cavities and less diffuse between the cavities, while the CH 4 molecules are distributed within the cavities more close to atoms H4 and HT and more diffuse between the cavities with the optimized DE parameters in the simulations (see Figure 5).The self-diffusion coefficients for loadings 1, 2, 4, 6, and 8 of CH 4 molecules per unit cell obtained by simulations are shown in Table 3.The   values obtained from AMBER parameters are in agreement with previous works using the same parameters set.For all loadings, the   values obtained by the DE parameters are approximately one magnitude higher than those obtained by AMBER parameters.The higher distributions and diffusion obtained from DE parameters can be explained by interaction profiles and the obtained RDFs.Since the CH 4 molecules are more close to the imidazolate ligands which are flexible, the movements of CH 4 molecules are additionally induced by the ligands, and locating around the cage windows increases the chance of cage-to-cage migration.By the way, these obtained   values agree well with some experimental values that are in a range in the order between 10 −11 and 10 −9 m 2 /s [2,8,[20][21][22].More comparative works between experiments and simulations are required for more discussions.

Conclusions
The molecular dynamics simulations are performed for investigating the structural and dynamical properties of CH 4 molecules in the ZIF-8 framework.A total of 1377 configurations of CH 4 -[C 4 N 2 H 5 ] − complexes are used for obtaining the binding energies with BSSE correction at a level of B2PLYPD/6-31G(d).Those energies are used as the data in parameters optimization.Since the DE method is based   The molecular distributions obtained by the simulations with AMBER parameters show that the CH 4 molecules are distributed within the cavities and there is little diffusion between the cavities, while there is more diffusion between the cavities when using DE parameters.Moreover, the selfdiffusion coefficients obtained by the optimized parameters are higher than those obtained by the AMBER parameters for all loadings.These indicate the effects of flexible imidazolate ligands on the movements of their bounded CH 4 molecules.
N 2 H 5 ] − ) ligands with a chemical formula of Zn[C 4 N 2 H 5 ] 2 .It is applied for the separation of CH 4 from natural gas and the reduction of CH 4 in the atmosphere.Structures of CH 4 molecule, ZIF-8 framework, and [C 4 N 2 H 5 ] − molecule are shown in Figure 1.

Figure 2 :
Figure 2: Some configurations of CH 4 -[C 4 N 2 H 5 ] − complexes.In each complex, a methane molecule is started from the initial point and moved out along an arrow line.

Figure 5 :
Figure 5: Diffusion of CH 4 molecules in ZIF-8 using (a) AMBER parameters and (b) DE parameters.

4 N 2 H 5 ] − complexes are illustrated in Figure
the energy of CH 4 molecule, and the energy of [C 4 N 2 H 5 ] − molecule with the basis sets of both molecules, respectively [10].The Lennard-Jones function is used to estimate the intermolecular interaction between CH 4 with zero charges and [C 4 N 2 H 5 ] − molecules:

Table 2 :
Parameters   and   (kcal/mol) obtained with AMBER force field and DE method.

Table 3 :
The self-diffusion coefficients (10 −9 m 2 /s) of CH 4 -ZIF-8 obtained with AMBER and DE parameters.The site-site RDFs between atoms of the framework to atoms H and C of CH 4 molecule indicate that atoms H4 and HT are favored sites of adsorption obtained by all simulations.The parameters of DE give shorter average minimum distances from atoms H and C of CH 4 to all atoms of the framework than those obtained with AMBER parameters.These obtained results show correspondence between the structural properties and the energies profiles.