Prediction of the infrared spectrum of water on graphene substrate using hybrid classical/quantum simulation

Surficial water can have an important influence on physical, chemical, and biological properties. An understanding of surficial water is thus crucial for technological advancement. Spectroscopic methods, such as IR spectroscopy, are among the most powerful tools for investigating the structure of surficial water. IR spectroscopy has been used to measure the vibration of water clusters in order to clarify the complex hydrogen bond networks of water molecules.^1,2) Molecular simulations are also useful for elucidating the microscopic structure of surficial water molecules.^3,4) Simulation results can be used to complement the results of experiments and facilitate microscopic interpretations of the behavior of surficial water.

In this paper, we calculate the IR spectrum of surficial water on a graphene substrate using a hybrid method that combines classical molecular dynamics (MD) and density functional theory (DFT) analyses. Graphene is the simplest substrate for exploring the IR spectrum of surficial water because it is atomically flat and nonpolar. The combination of classical MD and DFT analyses was successfully applied to estimate the IR spectrum of bulk water in our previous paper.⁵⁾ In this method, a DFT simulation (frequency analysis) is performed using locational information obtained from a classical MD simulation. A time-dependent IR spectrum can be obtained by integrating multiple trials of the frequency analysis. The calculated spectra were found to be in good agreement with the experimental IR spectra.⁵⁾

Water on a graphene substrate forms a multi-layer structure.^6–8) In the present study, we prepared two models: one with a single layer of water molecules in direct contact with the graphene (1L model), and the other with an additional layer not in direct contact with the graphene (2L model). The thickness of the lower layer was 0.5 nm (defined as the first layer), and that of the upper layer was 0.2 nm (defined as the second layer).⁶⁾ In this paper, the water in the first layer is denoted as "surficial water". The purpose of this study is to find distinctive features in the IR spectrum of surficial water by comparing it to that of bulk water and clarifying the origin of these features.

First, classical MD simulations were conducted for 20 H₂O molecules (1L model) and 35 H₂O molecules (2L model) on a graphene layer consisting of 60 carbon atoms. The rectangular simulation boxes were 1.28 nm × 1.23 nm × 3.36 nm with periodic boundary conditions. The classical equation of atomic motion was numerically solved by the velocity Verlet algorithm with a time step of 1.0 fs. The IR spectrum is influenced by the water clusters, and the structures of the water clusters are partially dependent on the diffusion coefficient. Thus, the MD temperature was set to 250 K with a Nose thermostat,⁹⁾ at which the water diffusion coefficient in these MD simulations corresponds to the experimental value obtained at 300 K.¹⁰⁾ The long-range Coulomb interactions were calculated using the Ewald method.¹¹⁾ The classical MD simulations were conducted using MS Forcite Plus in Materials Studio 8.0¹²⁾ with COMPASS II¹³⁾ (Condensed Phase Optimized Molecular Potentials for Atomic Simulation Studies). The charge distributions for hydrogen and oxygen in the water molecules were set to −0.41 and 0.82 eC, respectively, and the carbon in graphene was set to have no charge. DFT simulations were continuously executed for the boxes obtained in the classical MD simulations with periodic boundary conditions. The geometric coordination was minimized with the oxygen atoms of water molecules fixed, and the IR spectra were calculated via diagonalization of the dynamical matrix (Hessian matrix) of water molecules. Hence, the structure of graphene was assumed to be fixed. The DFT calculations were iterated 20 times for simulation boxes extracted every 200 fs from the classical MD simulations. The geometrical optimization and frequency analysis of water molecules were performed using DMol3 in Materials Studio 8.0 using the Perdew–Wang 1991 generalized gradient approximation¹⁴⁾ (PW91) with a triple numerical basis with polarization functions¹⁵⁾ (TNP).

Figure 1 shows the IR spectra of the first layer in the 1L and 2L models and in bulk water⁵⁾ calculated using the hybrid simulations, and Table I shows the wavenumber of each mode and the intensity ratio of the stretching mode to the bending mode for water molecules. The total IR spectra are composed of the derivations from the graphene substrate and water in the first and second layers. In order to estimate the IR spectra of surficial water, we extracted the normalized intensity for each vibration mode that originated from water molecules in the first layer. The IR spectra of water in the first layer are what remain after subtracting the IR intensity associated with the graphene and water in the second layer from the total IR spectrum. All spectra were smoothed by the following Gaussian function,

$$\begin{eqnarray}&&S(x)=\displaystyle \frac{{A}_{i}}{\sigma \sqrt{2\pi }}\exp \left(-\displaystyle \frac{{\left(x-{x}_{i}\right)}^{2}}{2{\sigma }^{2}}\right)\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&\sigma =\displaystyle \frac{{\rm{FWHM}}}{2\sqrt{2{\rm{ln}}2}}\end{eqnarray} \tag{ 2 }$$

where A_i is the intensity of the line spectra, x_i is the wavenumber of the line spectra, and FWHM is the full width at half maximum, which is set to 60 cm⁻¹ in this study.

Zoom In Zoom Out Reset image size

Three obvious bands can be observed for the 1L, 2L and bulk models in Fig. 1, which correspond to symmetric and asymmetric stretching vibrations (3000–3800 cm⁻¹), intramolecular bending vibrations (1500–1700 cm⁻¹) of water molecules, and librational motion (0–1000 cm⁻¹).¹⁶⁾ These bands appear for both the surficial water and bulk water. The wavenumber for the bending mode of water on graphene is very similar to that for the bulk water. The difference between them is no more than 7 cm⁻¹. The intensity ratios for the surficial water are also close to that for bulk water. On the other hand, the stretching mode for surficial water is shifted to a higher wavenumber, by approximately 20 cm⁻¹, compared to that for the bulk water. In addition, the small bands around 2950 cm⁻¹ appear only for the surficial water.

The stretching band for water depends on the length of the OH bonds in the water molecules.^17,18) The longer the OH bond becomes, the lower the wavenumber for the stretching mode becomes. Furthermore, the OH bond length is closely related to the strength of the hydrogen bond. Figure 2 shows IR spectra of representative small water clusters on graphene (purple curves) and isolated water clusters (green curves). The aerial structure of each cluster is shown on the right-hand side of Fig. 2. The cluster names were taken from Ref. 5 and Ref. 19. Some previous studies^20,21) have reported that the water clusters are adsorbed on a graphene substrate. The geometrical optimization and frequency analysis were conducted separately from the main calculations using DMol3 in Materials Studio 8.0 with PW91/TNP. The harmonic bands for free OH bonds associated with monomers in this case are over 3800 cm⁻¹ for both water clusters and graphene and isolated clusters. The length of these free OH bonds is 0.967 Å. The hydrogen-bonded OH bonds in the dimers produce harmonic bands at 3500–3700 cm⁻¹. The length of these bonds is 0.975 Å. The wavenumber associated with the non-hydrogen-bonded OH in the dimer, whose bond lengths are the same as those in the monomer, is close to the wavenumber for the monomer. In general, the harmonic bands for hydrogen-bonded OH shift to lower wavenumbers as the number of water molecules in a water cluster increases because the hydrogen bonds grow stronger.²²⁾

Zoom In Zoom Out Reset image size

Figures 3(a)–3(c) show the distributions of OH bond lengths after DFT geometrical optimization for bulk water, surficial water in the 1L model, and surficial water in the 2L model, respectively. The most frequent length of OH bonds is 0.985 Å for the bulk water model, and 0.975 Å for surficial water in the 1L and 2L models. Thus, the harmonic modes for surficial water shift to higher wavenumbers compared to those for the bulk water. The change in the OH bond lengths can be explained by the specific surficial water cluster structures. In bulk water, the major water cluster type is a tetrahedron (three-dimensional four-coordinated), whereas on graphene, the water cluster types are the two-dimensional tetramer, pentamer, and hexamer.⁶⁾

Zoom In Zoom Out Reset image size

Figure 4 shows the number distribution of hydrogen bonds. In this paper, the number of hydrogen bonds are defined as the number of donors from one water molecule. If both OH bonds in a water molecule interacted with another water molecule via a hydrogen bond, "2 bonds" was recorded. The label "0 bonds" was recorded for the monomer of a water molecule. For hydrogen bonds the O–O distance between two water molecules should be less than 3.60 Å, the O–H distance should be less than 2.45 Å, and the H–O–⋯–O angle should be under 30^◦.²³⁾ According to Fig. 4, the number of monomers in bulk water is only slightly higher than that in surficial water, whereas the number of "1 bond" structures is higher and that of "2 bonds" structures is lower in surficial water compared to those in bulk water. The water molecules on the surface form polygonal clusters; however, not all OH bonds make hydrogen bonds. In the 1L model, there are no water molecules above the first layer, and thus upward-facing OH bonds become free bonds. In addition, the main reason why the number of free OH bonds is higher in surficial water, which is also applicable to the 2L model, is the existence of a surface for which hydrogen bonds with downward-facing OH bonds are not permitted. In other words, the quasi-two-dimensional structure of water clusters results in a higher number of free OH bonds.

Zoom In Zoom Out Reset image size

We now focus on the OH bonds in Fig. 3 whose length is more than 1.01 Å. Small peaks appear in the histograms for the 1L [Fig. 3(b)] and 2L [Fig. 3(c)] models at a length of about 1.01 Å; these peaks do not appear in the histogram for bulk water (Fig. 3(a)). Figure 3(d) shows the histograms in Figs. 3(a)–3(c) smoothed using Eqs. (1) and (2) with FWHM = 0.004 Å. A small peak appears at a length of about 1.01 Å in the 1L and 2L models, which is in good agreement with the small sub-band in the IR spectra of surficial water [Fig. 2]. Vibrational bands of less than 3000 cm⁻¹ mostly originated from OH bonds with a length of more than 1.01 Å as determined by the atom eigenvectors. The long OH bonds indicate strong hydrogen bonds. Water molecules on graphene are densely packed (approximately three times denser than those in bulk water⁶⁾); therefore, some water molecules become good donors or acceptors, depending on the surrounding environment.

We predicted the IR spectrum of surficial water on a graphene substrate using a hybrid method that combines classical MD and DFT,⁵⁾ and compared the results with the calculated IR spectrum of bulk water. The bending modes for surficial and bulk water were similar. For surficial water, the stretching mode shifted to a higher wavenumber than that for the bulk water and a small sub-peak appeared at around 2950 cm⁻¹. The former is caused by free OH bonds and the latter is caused by long OH bonds with a length of more than 1.01 Å. The differences between the calculated IR spectra for bulk and surficial water were explained based on the simulation results in this study; however, the IR spectrum of surficial water on graphene has not yet been obtained experimentally. We should compare the calculated spectrum to the experimental one once it is available.

Acknowledgments