Download Study of Azide Ion Behavior in Carbon Nanotubes: Interactions and Energy Levels and more Summaries Chemistry in PDF only on Docsity!
A Novel Intermolecular Potential to Describe the
Interaction Between the Azide Anion and Carbon
Nanotubes
Stefano Battagliaa,b,1, Stefano Evangelistia,2, Thierry Leiningera,3, Fernando Piranib,c, Noelia Faginas-Lagob,4,∗ aLaboratoire de Chimie et Physique Quantiques, IRSAMC, Universit´e de Toulouse et CNRS, 118 Route de Narbonne, F-31062 Toulouse Cedex, France bDipartimento di Chimica, Biologia e Biotecnologie, Universit`a degli Studi di Perugia, Via Elce di Sotto 8, I-06123 Perugia, Italy cIstituto di Scienze e Tecnologie Molecolari CNR (CNR-ISTM), Perugia - Italy
Abstract
In this contribution we propose a novel and accurate intermolecular potential that can be used for the simulation of the azide anion confined inside carbon nanotubes of arbitrary size. The peculiarity of our approach is to include an explicit term, modeling the induction attractive contributions from the negatively charged azide ion, that can be generalized to other ions confined in carbon nanotubes of different size and length. Through a series of accu- rate DLPNO-CCSD(T) calculations, we show that this potential reproduces the ab initio interaction energy to within a few kcal/mol. The potential is implemented in a molecular dynamics program, with which we carried out illustrative simulations to demonstrate the effectiveness of our approach. At last, the guidelines provided by this investigation can be applied to build up force fields for many neutral/ionic molecular species confined within carbon
∗Corresponding author Email address: noelia.faginaslago@unipg.it (Noelia Faginas-Lago) (^1) ORCID: 0000-0002-5082- (^2) ORCID: 0000-0001-8782-443X (^3) ORCID: 0000-0002-7373- (^4) ORCID: 0000-0002-4056-
Preprint submitted to Elsevier October 17, 2019
nanotubes; a crucial requirement to carry out molecular dynamics simula- tions under a variety of conditions.
Keywords: Azide Anion, Carbon Nanotubes, Induction, Intermolecular Potential
- Introduction
The hollow structure of carbon nanotubes (CNTs) and the low reactivity of the internal wall provide an ideal environment to confine a wide variety of systems in their cavity, making them suitable for different types of applica- tions. Since 1993, when Ajayan and Iijima [1] experimentally demonstrated the possibility to fill CNTs by capillary suction proving feasible the early theoretical predictions, a lot of work has been carried out in this domain such that a large variety of confined and composite systems have now been investigated; e.g. fullerenes and metallofullerenes[2, 3], a long list of differ- ent metals[4–10], energetic molecules[11–14], and many more. The limited space inside the cavity of CNTs can also induce new phases of materials, with the most notable example being certainly that of water. Novel “ice” phases arising from encapsulated water were first predicted from theoretical calculations[15] and then observed experimentally over the years[16–18].
In the search for alternative and environmentally friendly energy sources, a promising class of systems is that of all-nitrogen molecules[19, 20]. Account- ing for as much as 78% of the total, nitrogen is the most abundant species in Earth’s atmosphere and is therefore widely and easily accessible. It is virtu- ally only present as di-molecular N 2 , which is by far the most stable form and it is not harmful for the environment. The energy difference between single, double and triple covalent bond is such that polynitrogen molecules store an enormous amount of energy and for this reason pure nitrogen clusters are considered as high energy-density materials (HEDM)[19, 20]. The first and most famous polynitrogen molecule is the azide anion, N− 3 , known since the end of the nineteenth century[21]. It took more than a
of the azide ion in crystal or gas phase[34–42] as well as solvated[43–55], N− 3 was never studied in a confined environment. Considering that N− 3 often represents a precursor in the synthesis of energetic molecules, e.g. the afore- mentioned cg-N phase synthesis starts from sodium azide[33], or it appears as an intermediate structure during the dissociation reaction of larger nitrogen clusters, e.g. the barrierless dissociation of the N− 5 ion[26, 56], it is certainly of great interest to study its behavior in a constrained environment such as that provided by carbon nanotubes, and in particular to understand the type of interactions between the confining and the confined fragments. In this work, the confinement of the azide anion inside carbon nanotubes of different lengths and diameters are presented and a newly developed inter- molecular potential is reported. A first part of this manuscript is devoted to an ab initio study of this system, providing important insight into the ad- sorption of the polynitrogen molecule in the cavity, followed by a systematic investigation where the finite size effects are studied by high-level coupled- cluster calculations, providing accurate, reference energies. Building on the knowledge acquired from the static calculations and compar- ing to the reference ab initio data, a novel intermolecular potential modeling the non-covalent interaction between N− 3 and carbon nanotubes is then de- rived and presented. The potential is implemented in a local version of the molecular dynamics (MD) program DL POLY 4.08[57, 58] and the working equations of the potential are reported here. The implementation of this po- tential allows to accurately study the dynamics of the confined ion in carbon nanotubes of any size and type, and in principle solvated in any non-ionic liq- uid. The implementation of the potential is tested by performing simulation of the confined system which reproduce well the ab initio results. Moreover, we should emphasize that the methodology proposed here to characterize static and dynamical properties of the prototype N− 3 molecular ion confined in CNTs can be also generalized to build up the force fields (FF) for many other neutral/ionic molecular species encapsulated in CNTs. The availability
of FF given in analytical form is a crucial condition to carry out molecular dynamics simulations under a variety of conditions and it is therefore of gen- eral interest for the characterization of the behavior of prototype polyatomic molecular ions confined in a nanotube. The paper is structured as follows. In section 2 we present the ab initio results along with the computational details of the performed calculations. In section 3 we focus on the derivation of the intermolecular potential and the parameters entering its formulation. In section 4 we present the results of MD simulations carried out to test the potential and finally in section 5 we draw the conclusions of this study.
- Ab initio Calculations
In order to propose the new intermolecular potential, we first carried out a series of ab initio calculations. In this part, carbon nanotubes were treated as finite-size systems, with the addition of hydrogen atoms at both ends in order to fill the valence shell. Only the armchair CNT(5, 5) was considered, as this was previously reported to be the CNT showing strongest interaction with the azide anion[14]. Clearly, the interaction between the CNT and the azide strongly depends on the length of the former, a detail that is often not discussed in similar studies. Thus, in order to obtain converged values with respect to the length of the nanotube, CNTs with different number of units were used. We shall use for the rest of this article the notation Λ-CNT(5, 5) to label a nanotube composed by Λ units, where a single unit corresponds to a [10]cyclophenacene.
2.1. Geometry
The geometries of the nanotubes and the azide anion were optimized us- ing restricted Kohn-Sham DFT, employing the B97D3 exchange-correlation functional with Becke-Johnson damping[59, 60]. For both fragments, Dun- ning’s triple-ζ cc-pVTZ basis set[61] was used, with the addition of diffuse
distance and orientation were relaxed, thus locating the optimal adsorption distance. One can imagine that, given the highly periodic structure of the carbon nanotubes, the number of minima in the PES is extremely large, however many of those are expected to be quite close in energy. Remarkably, it turns out that in the particular case of an azide anion confined inside a CNT(5, 5), the geometry with the N− 3 molecule placed exactly in the center of the nanotube and aligned to its principal axis results to be a local mini- mum. The lowest minimum found with this procedure is nevertheless very similar, with the azide anion still nearly perfectly aligned to the principal axis, but displaced by approximately 0.5 ˚A towards one of the two openings. The energy difference between these two minima is negligible (less than 0. 6 kcal/mol) and thus, for the sake of simplicity, the geometry with the N− 3 ion placed in the center of the CNT was used for the calculation of Eint at higher levels of theory and for nanotubes longer than 9 units. As already mentioned before, it is important to study the dependence of the interaction energy with respect to the length of the nanotube. In particu- lar, given that the confined species is a negatively charged ion, we should expect a slower convergence of this quantity compared to uncharged systems due to the more diffuse electron density. To obtain highly accurate inter- action energies for systems as large as the carbon nanotubes studied here, approximations to canonical post-Hartree-Fock methods were employed. In particular, the energy were obtained using the density-fitted spin-component- scaled MP2 method (RI-SCS-MP2)[69] and the local coupled cluster method DLPNO-CCSD(T) as implemented in the ORCA 4 program package[70, 71]. The cc-pVTZ basis set was used where computationally feasible, otherwise the cc-pVDZ one. Diffuse functions (aug-) were added on the nitrogen atoms in both cases. In order to assess the basis set incompleteness error, for the nanotube with Λ = 3, we have carried out an RI-SCS-MP2 calculation using the larger (aug-)cc-pVQZ basis set and extrapolated the results to the com- plete basis set limit (CBS) using the triple- and quadrupole-ζ bases energies.
In particular, for each fragment, we have extrapolated the RI-SCS-MP2 cor- relation energy (i.e. not the total RI-SCS-MP2 energy), labeled and then computed the interaction energy using the CBS values. In order to do so, we have used the two-point extrapolation scheme by Halkier et al. [72], with the formula Ecorr(∞) = E(X)X
3 − E(Y )Y 3
X^3 − Y 3
where X and Y are the cardinal numbers of the basis set (X = 3 and Y = 4 for triple- and quadrupole-ζ, respectively) and E(X) and E(Y ) are the corresponding energies. The interaction energies obtained for the three basis sets (we show also the values obtained with the double-ζ basis) and the CBS limit are reported in Table 1. As we can see, going from triple- to quadrupole- ζ, the interaction energy changes by less than 1 kcal/mol, whereas going to the CBS limit, the difference with the triple-ζ is only about 1 kcal/mol. In contrast, we clearly see that the double-ζ results underestimate substantially the interaction strength. According to this analysis we can conclude that the errors introduced by the finite basis will be very likely less than those introduced by the classical approximation of the force field.
The results obtained with the different methods, using the optimal geome- tries of the fragments and with the azide anion placed exactly in the center of the CNT, are shown as a function of the nanotube length in Figure 1. For all CNTs long enough to confine the azide anion, the reported energy is neg- ative, meaning a favorable interaction between the nanotube and the guest species (note that the case Λ = 1 is excluded from Figure 1 because given that the CNT is extremely short, the azide anion is in practice not confined anymore and the resulting interaction energy is positive, corresponding to the azide trying to “escape” the nanobelt). From the plot, it is immediately clear how the interaction energy has a slow convergence with respect to the number of units Λ, a behavior that we shall investigate further. Coupled-cluster-level energies were obtained in conjunction with a double-ζ basis set only, whereas a few points at triple-ζ quality are shown for RI-SCS-
action (in classical terms) were identified for the system under study: van der Waals (vdW) interactions (these are always present), induction effects due to the charged N− 3 fragment and electrostatic effects between the polarized nanotube and the anion. We shall see how these different interactions can be accurately modeled and added together to form a complete description of all non-covalent interactions present in the system.
3.1. Model Potential Definition
Van der Waals interactions are modeled according to the Improved Lennard- Jones (ILJ) potential, a modified form of the canonical Lennard-Jones po- tential proposed a few years ago by Pirani et al. [74]. The general formula reads
Vilj (rij ) =
[
m n(rij ) − m
rm rij
)n(rij ) − (^) n(rn(rij^ ) ij )^ −^ m
rm rij
)m] (5)
where
n(rij ) = β + 4. 0
rij rm
The indices i and j refer to atoms of the system which do not belong to the same molecule. For each unique pair of interacting atom types, there are three parameters to be set. The value rm corresponds to the position of the minimum of the potential energy surface for the interaction of the two atoms labeled by i and j. The depth of that minmum is given by the value of and the last parameter is β, appearing in Equation (6), which is related to the hardness of the system and usually varies between 7 and 9[74]. The value of m entering Equation (5) as an exponent and a prefactor is set according to the type and charge of the interacting centers: for atom–atom interactions m is equal 6, for ion–atom interactions is equal 4 and for ion–ion is set to 1. This is therefore not a freely varying parameter, but a value which is set a priori. The ILJ potential has proven very effective in a number of applications,
in particular dealing with polar molecules such as water[75], water-carbon systems[76], ion-π interactions[77, 78], and in general gases interacting with carbon nanostructures such as graphene and carbon nanotubes[79–84]. Albeit the ability of the ILJ potential to accurately reproduce the potential energy surface of interacting ion–neutral species, it is most often the case that the ions in question are individual atoms; this fact allows for a pairwise treatment of the interaction as it is done, e.g., in the ILJ potential. However, the presence of arbitrary, non-spherically-symmetric molecular ions substan- tially complicate this picture. The latter are responsible for strong induction effects acting on the interacting species, which are non-additive in nature[85] and increase the complexity of the model potential.
The induction effects generated by the presence of N− 3 can be modeled according to a relatively simple formula which arises directly from the general formula of an induced dipole, i.e.
μ∗^ = αE (7)
where μ∗^ is the induced dipole, α is the polarizability and E is the external electric field. In the system we are considering, the external electric field is the one generated by the partial atomic charges of the azide, whereas the polarizability is that of the carbon nanotube. For the specific case of the electric field generated by three point charges, the induction potential, obtained by integration of Equation (7) is given by
Vind(ri 1 , ri 2 , ri 3 ) = −^12 αi
[
q 1 r^2 i 1 +^
q 2 r^2 i 2 +^
q 3 r i^23
] 2
where ri 1 , ri 2 and ri 3 are the distances between atom i of the nanotube and the three nitrogen atoms of the azide anion, q 1 , q 2 and q 3 are their partial atomic charges and αi is the polarizability per volume of atom i. Note that the choice to model the electric field with three partial charges, one for each
3.2. Model Potential Parameters
With this choice of potential, we are now set to test its accuracy by com- paring the interaction energies obtained with the DLPNO-CCSD(T) method with those obtained from Equation (10). In this study, we used a combina- tion of experimental and ab initio data, together with correlation formulas derived from extensive empirical studies[86–88] to obtain the parameters used in the potential function. The most important parameters are the atomic po- larizabilities, which are used to determine the and rm values for Vilj , as well as appearing directly in the induction potential. The nitrogen atomic polar- izability is obtained from the experimental total mean polarizability of the azide anion[89], ¯α = 4.65 ˚A^3 , which was distributed on the three atoms ac- cording to their partial charges, resulting in a value of 1.9 ˚A^3 for the external atoms (charges q 1 = q 3 = − 0 .56 atomic units (au)) and a value of 0.85 ˚A^3 for the internal one (charge q 2 = +0.12 au). The polarizability of the hydrogen atoms was taken from Ref. 90 and set equal to 0.380 ˚A^3 , whereas that of carbon atoms from Ref. 91, where we note that the actual value depends on the specific curvature (diameter) of the nanotube. The polarizability for carbon atoms, along with the calculated ILJ parameters according to Ref. 86 are reported in Table 3. Similarly, the parameters for hydrogen and the hydrogen–nitrogen interaction are listed in Table 4. Atomic partial charges for the carbon nanotube were obtained through a natural population anal- ysis (NPA)[92] of the electron density calculated at B97D3/cc-pVTZ level of theory. Instead, those for the nitrogen atoms in the the electrostatic and induction potentials were set equal to q 1 = q 3 = − 0 .56 au (external atoms) and q 2 = +0.12 au (internal atom) according to the ab initio values reported by Le Borgne et al. [89]. As a matter of comparison, an NPA of N− 3 obtained with B97D3/aug-cc-pVTZ (both geometry and electron density) results in partial charges q 1 = q 3 = − 0 .55 au and q 2 = +0.10 au.
3.3. Potential Energy Minimum and Profile
To assess the accuracy of the potential energy function we compare the interaction energy obtained using the ab initio methodology and the one obtained with the model potential for all nanotubes considered so far. In Figure 2 the ab initio and the classical potential energies are reported for β = 8 as a function of the nanotube length Λ. As can be seen, for the short- est nanotube, the classical potential largely underestimates the interaction strength, however, starting from 7 units, the agreement is excellent. The fact that different types of interactions are explicitly considered naturally allows for an energy decomposition analysis. The plot depicted in Figure 3 shows the same two curves as in Figure 2, however, the relative contribu- tions of the three components Vilj , Vind and Vels are also reported. It is interesting to see the different ranges of the types of interaction. The vdW contribution described by the ILJ potential account for only a fraction of the total interaction energy. Moreover, the convergence with respect to the length of the nanotube is essentially immediate. Despite the fact that vdW interactions are usually called long-range in the quantum chemistry context, because exchange interactions usually decay exponentially as function of the distance, in this framework they actually should be considered as short-range due to their r−^6 dependence. The induction potential appears to converge also quite quickly with respect to the nanotube length. In this case however, the dependence is given by the inverse fourth power of the distance. Here, the magnitude of this contribution is very surprising; the explicit modeling of this interaction seems therefore crucial in order to provide a good description of the potential well. At last, the longest-range interaction is given by the Coulomb potential, whose strength is inversely proportional to the distance between the atomic centers. From this component it is clear how the shape of the total classical potential represented in Figure 2 is dominated by this term. Since between 5 and 7 units, the improved Lennard-Jones and the in- duction potentials are essentially constant, the convergence behavior of the
to the minimum of the potential well obtained in the first case; i.e. the opti- mal adsorption distance of the azide anion confined in the carbon nanotubes considered here. The corresponding profiles are shown in Figure 6, where as before the plots are arranged on top of each other, starting from the top with the CNT(5, 5). In all cases, if initially placed outside, the azide anion is attracted towards the nanotube and depending on its relative position with respect to the CNT axis a local minimum is present. An interesting result is the flipping of electrostatic interaction, which is attractive when the azide is outside the CNT since the external hydrogen atoms are positively charged, and repulsive inside since there is excess negative charge on the carbons. The last case analyzed is the rotation inside the nanotube reported in Fig- ure 7. This last series of profiles shows very clearly the effects of spatial confinement provided by the enclosing CNT. In a CNT(5, 5) the ion can barely rotate and tremendous repulsive walls build around the minimum of the PES. On the other hand, already with the slightly larger CNT(7, 7) the (repulsive) effects are only mildly perceived by the ion. For the largest nan- otube, the potential is basically flat, suggesting that a single N− 3 ion is likely to freely rotate even when confined.
- Molecular Dynamics Simulations
The intermolecular potential presented in the previous section was imple- mented in a locally modified version of the DL POLY 4.08 program[57, 58]. In particular, analytical energy gradients were derived for both the ILJ and the induction potentials, whereas the Coulomb potential is already available as part of the original software. It should be noted that the induction po- tential, being non-additive in nature, is forcedly hard-coded for the specific case of an electric field generated by three point charges. Nevertheless, there is no particular assumption towards the azide anion: any fragment modeled by three point charges can make use of this potential.
In order to validate the implementation, a series of MD simulations was
carried out on the same systems for which the energy profiles were presented, i.e. the azide anion confined in a CNT(5, 5), a CNT(7, 7) and a CNT(9, 9). In this case however, the carbon nanotubes were considered with PBC and hence without addition of hydrogen atoms at the two ends. This removes the need of electrostatic interactions, since the nanotube is not polarized. The simulations were carried in the canonical ensemble (NVT), using the Nos´e-Hoover thermostat with a relaxation constant of 0.5 ps at a tempera- ture of 300 K. Orthorombic periodic boundary conditions were applied, with a simulation box size of 40 × 40 × 32 .0036 ˚A^3 , where the CNT principal axis was aligned to the shortest edge of the box and centered with respect to the other two dimensions. We carried out an equilibration phase of 2 ns and a production phase of 3 ns, using a timestep of 1 fs. Statistical sampling was recorded every 10 ps. To minimize the number of parameters in these sim- ulations, the two fragments were considered frozen for the entire simulation such that no intramolecular potential was necessary. However, the introduc- tion of flexibility in the nanotube and the azide anion are easily possible. For each nanotube we repeated the simulation 10 times, each time starting from a new, randomly generated position of the azide ion inside the nanotube.
In Figure 8 we report the average adsorption energy as a function of time. Although the fluctuations for individual runs are quite severe due to the rigid-body approximation and to the use of a single azide ion, the adsorp- tion energy averaged over the 10 runs started from different configurations is significantly smoother. The mean energy values (i.e. the temporal average of the adsorption energy) is − 62 .98, − 42 .32 and − 34 .77 (in kcal/mol) for CNT(5, 5), CNT(7, 7) and CNT(9, 9), respectively, and is in good agreement with respect to the static calculations (CNT(5, 5)). In particular, the ex- trapolated DLPNO-CCSD(T)/cc-pVTZ interaction energy was estimated to be − 65 .93 kcal/mol, thus within only a few kcal/mol to the mean energy obtained with the dynamics. The advantage of such an intermolecular po- tential is that the interaction energy for larger systems such as the CNT(7, 7)
carbon nanotubes. This model is based on a series of high-level ab initio cal- culations performed specifically for this work as well as on a previous study carried out in our group[14]. This particular potential has been implemented in a local version of the molecular dynamics program DL POLY 4.08[57, 58], with which we have demonstrated its performance with a series of illustrative simulations on the confined azide anion. On one hand, in the construction of this model, we carried out extensive ab initio calculations which complement the work done in Ref. 14 and in partic- ular we have addressed the very important issue of length dependence in the calculation of interaction energies for species confined in carbon nanotubes. When working with finite-size systems, it is crucial to investigate the behav- ior, in principle, of any molecular property as a function of the length in order to obtain quantitative results. On the other hand, we have found that MD simulations based on the potential proposed here are able to reproduce the interaction energies obtained with the DLPNO-CCSD(T) method to within a few kcal/mol, which is a small fraction of the total interaction energy. Such a model retains the computational efficiency typical of classical molecular dynamics simulations, providing however a very good accuracy compared to the (much) more expensive ab initio calculations. Several other advantages apply as well. First and foremost the access to the dynamical properties of the system under a variety of conditions: indeed, this potential can be used in conjunction with any type of statistical ensemble. In our particular ex- ample, we have computed the mean adsorption distance of the azide anion confined in nanotubes of different diameters at constant volume and tem- perature, with the surprising result that larger nanotubes may have a more stabilizing effect than smaller ones. Another advantage is the general nature of this potential that can be applied to carbon nanotubes of arbitrary size and dimension, and in the limit even for graphene. In particular, we expect it to be applicable to carbon nanotubes of different chiralities, provided that model parameters are available. The latter can however be obtained through
experiments, calculations and/or phenomenological methods[76–79]. Open as well as periodic boundary conditions can be handled by this potential, as the electrostatic interaction term can be included or not according to it. As of now, there are two major limitations of the proposed model. First, it has been constructed for single-wall carbon nanotubes only. Its applica- bility on multi-wall CNTs is therefore at the moment not considered, nor suggested. A second issue is the tradeoff for the accuracy obtained: the ex- plicit modeling of the induction term is of foremost importance in order to obtain such accurate results, however, due to its non-additive nature, this term has been explicitly coded. This means that at present, only induction effects generated by a source of three point charges can be included. Actually, the presence of an explicit term describing induction effects is the key aspect of this work. The latter certainly reduces the transferability of the model proposed, although we should note that by borrowing techniques from the field of machine learning, such as automatic differentiation, it is in principle possible to generalize our approach to handle system with any number of point charges as sources[93]. At last, and perhaps most importantly, this study provides guidelines that can be used to create force fields for other neutral and ionic molecular species confined within carbon nanotubes, that still maintain the analytical form nec- essary to efficiently carry out molecular dynamics simulations on systems of any dimension.
- Acknowledgments
The results included in this publication have received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Sk lodowska-Curie grant agreement no642294. S. E. and T. L. ac- knowledge the “Programme Investissements d’Avenir” ANR-11-IDEX-0002- 02, reference ANR-10-LABX-0037-NEXT for financial support. N. F.-L. thanks MIUR and the University of Perugia for the financial support of the
