Konstantin S Smirnov
email: konstantin.smirnov@univ-lille.fr

Structure and sum-frequency generation spectra of water on neutral hydroxylated silica surfaces †

Structure and sum-frequency generation spectra of water on neutral hydroxylated silica surfaces

Introduction

Silicon dioxide is one of the most abundant minerals in the Earth's crust and it plays a key role in many natural geological and atmospheric processes. Silica is indispensable for numerous applications, such as glass industry, heterogeneous catalysis, manufacturing electronic and optical devices, and it also finds the usage in biomedical applications and pharmaceutical industry, e.g. as implants and drug delivery systems. 1,2 Silica/water interfaces are ubiquitous and play an important role in all the above mentioned areas. A detailed knowledge of structure of silica/water interface and of dynamics of species in the interfacial region is of a great value for the understanding of interfacial phenomena and underlying microscopic mechanisms.

A common problem in obtaining the information on surface processes is that the surface-to-bulk ratio is generally too small and therefore, the measured signal is dominated, as a rule, by a bulk contribution. With this respect, surface-sensitive and surface-specific experimental techniques, such as atomic force microscopy, 3,4 X-ray reflectivity 5,[START_REF] De Vries | Surface atomic structure of KDP crystals in aqueous solution: An explanation of the growth shape[END_REF] and nonlinear vibrational sum-frequency generation (VSFG) and second-harmonic generation spectroscopies [START_REF] Ong | Polarization of water molecules at a charged interface: second harmonic studies of the silica/water interface[END_REF][START_REF] Du | Vibrational spectra of water molecules at quartz/water interfaces[END_REF][START_REF] Bain | Sum-frequency vibrational spectroscopy of the solid/liquid interface[END_REF][START_REF] Tian | Recent progress on sum-frequency spectroscopy[END_REF] are of great help to get insight into the behavior of interfacial region. [START_REF] Zaera | Probing liquid/solid interfaces at the molecular level[END_REF][START_REF] Björneholm | Water at interfaces[END_REF] The surface specificity of the VSFG and second-harmonic generation methods arises from the broken inversion symmetry at the surface that results in a spectral response due to the second-order nonlinear polarizability χ (2) . [START_REF] Boyd | Nonlinear Optics[END_REF][START_REF] Shen | Fundamentals of Sum-Frequency Spectroscopy[END_REF] Among these methods, the VSFG technique is of particular interest as measuring χ (2) as a function of frequency ω yields the molecular-level information due to vibrational resonances characteristic of species present at the interface. Moreover, the phase-sensitive VSFG spectroscopy allows to obtain a direct information about the orientation of molecules on the surface from the sign of spectral features. [START_REF] Ostroverkhov | New information on water interfacial structure revealed by phase-sensitive surface spectroscopy[END_REF][START_REF] Myalitsin | Water structure at the buried silica/aqueous interface studied by heterodyne-detected vibrational sum-frequency generation[END_REF][START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF] A very recent review by Backus et al. [START_REF] Backus | The mineral/water interface probed with nonlinear optical spectroscopy[END_REF] summarizes the outcome of experimental studies of mineral/water interfaces, and silica/water ones in particular, by the VSFG and second-harmonic generation techniques.

Atomic level information about solid/water interfaces can be obtained by molecular modeling and among different approaches, molecular dynamics (MD) simulations are of particular relevance. The computations can be performed either with a set of effective interatomic potentials (classical MD) or by computing energies and forces with an electronic-structure method on-the-fly (ab initio MD), and the length-and time-scales of MD simulations nicely fit those probed by spectroscopic techniques. The behavior of water on silica surfaces was the subject of numerous modeling studies. The review by Rimola et al. [START_REF] Rimola | Silica surface features and their role in the adsorption of biomolecules: Computational modeling and experiments[END_REF] discusses the outcome of some of these works, while examples of more recent investigations can be found in ref. 20-31. Results of these studies have revealed a complexity of the structural organization of silica/water interfaces that was shown to depend on a number of factors such as surface crystallinity, morphology and charge, presence of surface functional groups and composition of aqueous phase. A number of modeling works on silica/water interfaces attempted to relate their structure of to the nonlinear spectra. [START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF][START_REF] Dewan | Structure of water at charged interfaces: A molecular dynamics study[END_REF][START_REF] Joutsuka | Effects of third-order susceptibility in sum frequency generation spectra: a molecular dynamics study in liquid water[END_REF][START_REF] Tuladhar | Ions tune interfacial water structure and modulate hydrophobic interactions at silica surfaces[END_REF][START_REF] Isaienko | Hydrophobicity of hydroxylated amorphous fused silica surfaces[END_REF][START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF][START_REF] Pezzotti | What the Diffuse Layer (DL) Reveals in Non-Linear SFG Spectroscopy[END_REF][START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF][START_REF] Chen | Molecular dynamics study of the electric double layer and nonlinear spectroscopy at the amorphous silica-water interface[END_REF][START_REF] Lyu | Potassium hydroxide concentration-dependent water structure on the quartz surface studied by combining sum-frequency generation (SFG) spectroscopy and molecular simulations[END_REF][START_REF] Smirnov | Structure and sum-frequency generation spectra of water on uncharged Q 4 silica surfaces: a molecular dynamics study[END_REF] As at normal pH conditions, silica surface possesses a negative charge due to a partial deprotonation of surface OH groups, the computational works considered charged surfaces, [START_REF] Dewan | Structure of water at charged interfaces: A molecular dynamics study[END_REF][START_REF] Joutsuka | Effects of third-order susceptibility in sum frequency generation spectra: a molecular dynamics study in liquid water[END_REF][START_REF] Tuladhar | Ions tune interfacial water structure and modulate hydrophobic interactions at silica surfaces[END_REF][START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF][START_REF] Pezzotti | What the Diffuse Layer (DL) Reveals in Non-Linear SFG Spectroscopy[END_REF][START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF][START_REF] Lyu | Potassium hydroxide concentration-dependent water structure on the quartz surface studied by combining sum-frequency generation (SFG) spectroscopy and molecular simulations[END_REF] addressed the effect of concentration and nature of charge-balancing ions on the spectra. [START_REF] Dewan | Structure of water at charged interfaces: A molecular dynamics study[END_REF][START_REF] Tuladhar | Ions tune interfacial water structure and modulate hydrophobic interactions at silica surfaces[END_REF][START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF] An attention has also been paid to a contribution of the third-order susceptibility and of surface electric field to the SFG signal. [START_REF] Joutsuka | Effects of third-order susceptibility in sum frequency generation spectra: a molecular dynamics study in liquid water[END_REF][START_REF] Tuladhar | Ions tune interfacial water structure and modulate hydrophobic interactions at silica surfaces[END_REF][START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF][START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF][START_REF] Chen | Molecular dynamics study of the electric double layer and nonlinear spectroscopy at the amorphous silica-water interface[END_REF][START_REF] Lyu | Potassium hydroxide concentration-dependent water structure on the quartz surface studied by combining sum-frequency generation (SFG) spectroscopy and molecular simulations[END_REF] The interpretation of VSFG experiments needs to take into account the heterogeneity of interfacial region and to separate contributions arising from different parts of aqueous phase and from the solid surface. Structure of solid/liquid interfaces is commonly described with electrical double layer or Gouy-Chapman-Stern model. [START_REF] Sparnaay | The Electrical Double Layer[END_REF][START_REF] Schmickler | Electrochemical theory: Double layer[END_REF] In this model, the interfacial region consists of a compact layer of species in the immediate proximity to the surface, called Stern layer, and a more distant diffuse layer (DL) loosely associated with the surface. Generally speaking, both the layers can contribute to the SFG signal and the knowledge of their spectral signatures is essential for the molecular-level description of aqueous interface. Experimentally, however, the different contributions to the spectrum are nearly inextricable. Making use of the double-layer theory, Shen et al. [START_REF] Wen | Unveiling microscopic structures of charged water interfaces by surface-specific vibrational spectroscopy[END_REF] suggested a scheme allowing to subtract the DL contribution and then successfully applied their method to obtain the phase-sensitive VSFG spectra of water layer next to a prototype charged surface at different pH values of the liquid phase.

Recently, Gaigeot and co-workers proposed a definition of different regions at solid/water interfaces based on the analysis of MD simulation results. [START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF] For this purpose, the authors used three descriptors, namely the water density, number of intermolecular hydrogen bonds, and the strength and directionality of H-bond network as a function of distance from surface. Molecules in a direct contact with the surface were assigned to a bonded interfacial layer (BIL), where the interaction with the surface was found to reduce the number of intermolecular hydrogen bonds and strongly alters the H-bond network. A more distant region, where water molecules have already local environment characteristic of bulk but the isotropy of H-bonds network is still absent, was associated with the diffuse layer. The authors argued that DL exists only at charged surfaces as it is arises from H 2 O molecules re-oriented by the surface electric field. Making use of the approach, subsequent studies by the same group explained the VSFG spectra of BIL and DL of water on surfaces of α-quartz and amorphous silica. [START_REF] Pezzotti | What the Diffuse Layer (DL) Reveals in Non-Linear SFG Spectroscopy[END_REF][START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF] Thus, characteristics of BIL spectrum were correlated with structure of water in the layer, surface crystallinity, and with the strength and directionality of interactions. [START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF][START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF] The VSFG spectrum of DL was found to arise from the third-order nonlinear susceptibility of bulk and the intensity and sign of spectral band are related to the strength and direction of the surface electric field. [START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF] Furthermore, the electric double layer formation, surface charge, and nominal pH of aqueous phase could also be inferred from the spectrum of diffuse layer, according to the authors. [START_REF] Pezzotti | What the Diffuse Layer (DL) Reveals in Non-Linear SFG Spectroscopy[END_REF] The behavior of water on charged silica surfaces are of great interest from both fundamental and technological viewpoints. Analysis of simulation results for such systems is, however, complicated by the existence of localized charges due to the deprotonation of silanol groups, resulting surface electric field and charge-balancing ions present in the aqueous phase. High computational expenses of ab initio MD method significantly limit the size of the system, range of possible ion concentrations and the time-scale accessible in the simulations. With this respect, classical MD simulations of mineral/water interfaces permit to investigate much larger systems for longer period of time, and yield phase-sensitive VSFG spectra in a good agreement with both experiments and ab initio MD. [START_REF] Smirnov | Structure and sum-frequency generation spectra of water on uncharged Q 4 silica surfaces: a molecular dynamics study[END_REF][START_REF] Kroutil | Phase-sensitive vibrational SFG spectra from simple classical force field molecular dynamics simulations[END_REF] Characteristics of neutral silica/water interfaces are important reference data and their analysis allows to get a revealing insight into the behavior of the more realistic systems. The purpose of the present work is to study the relationship between the structure of water on neutral hydroxylated silica surfaces and the VSFG spectrum of the interface as a function of the surface and liquid phase composition. The study is sought to answer the following questions:

(i) how do the surface morphology and density of surface OH groups influence the structure of water layer on neutral hydroxylated silica surfaces, (ii) how does specific structural organization of the interfacial water reflect itself in the nonlinear vibrational SFG spectrum, (iii) how does the chemical composition of aqueous phase affect water structure and impact the nonlinear spectra of the interface.

For this purpose, structural organization and nonlinear vibrational sum-frequency generation spectra of two series of silica/water interfaces were investigated by classical molecular dynamics simulations. The first series consisted of neutral crystalline and amorphous silica surfaces in contact with a neat water (NW), while the surfaces were exposed to an electrolyte solution (ES) in the second series of simulations. The latter series also aims at representing real experimental conditions where the aqueous phase always contains a certain number of solvated ions necessary to adjust the solution pH and to keep the ionic strength at a desired value. The next section presents silica surfaces used in the study, discusses simulation setup and effective interaction potentials, describes computed characteristics of the systems, and provides detail on the MD simulations.

Models

Models of silica surface.

In order to investigate the influence of surface characteristics on the structure and nonlinear spectra of silica/water interface, the work considered four models of silica surface with different morphology and density of surface hydroxy groups. A first model is the (20 2) surface of α-cristobalite called "CRB" throughout the paper, Fig. 1a. Coordinates of atoms in the CRB model were taken from the silica surface model database by Heinz et al. [START_REF] Emami | Force field and a surface model database for silica to simulate interfacial properties in atomic resolution[END_REF] The surface has terminal silanol groups with the density OH = 4.7 OH/nm 2 , and it exhibits two types of hydroxyls that differ in the position of oxygen atoms with respect to the surface. Outer hydroxy groups labelled "OH (OUT)" in Fig. 1a, are by 0.48 Å more lifted than inner "OH (IN)" groups. The second and third models were constructed using the structure of α-quartz 44 and correspond to the (101) and (001) surfaces. In what follows, these two models are respectively referred to as "Q101" and "Q001" and they are shown in Fig. 1b and Fig. 1c. The Q101 surface has vicinal silanols with the density OH = 5.6 OH/nm 2 and, similarly to the CRB surface, it features outer and inner OH groups with the former lying by 0.76 Å above the latter, Fig. 1b. The Q001 surface is characterized by the highest silanol density of 9.6 OH/nm 2 among the surfaces studied in the work and it exhibits rows of geminal OH groups. A half of these OH groups typically forms H-bonds with hydroxyls of a neighbor row and their OH bonds lie nearly parallel to the surface, these in-plane silanols are marked as "OH (IN)" in Fig. 1c. The remaining silanols orient their OH bonds outwards the surface and they are labelled as "OH (OUT)" in Fig. 1c. In addition to the three crystalline surfaces, a surface of amorphous silica with the silanol surface density of 4.7 OH/nm 2 was taken from the database by Heinz et al. [START_REF] Emami | Force field and a surface model database for silica to simulate interfacial properties in atomic resolution[END_REF] This surface is hereafter called "AS surface". All the surfaces are neutral, no deprotonation of OH groups is considered and the simulations correspond to the point of zero charge conditions, pH ∼ 2 -4. Nevertheless, the model surfaces can also be representative of water/silica interface at a higher pH (up to pH ∼ 6 -6.5). At this pH value, the surface charge is less than -0.15 |e -|/nm 2 (ref. 45-48) that gives one negatively charged SiO -group per 6.5 nm 2 at most, i.e. large surface areas are neutral.

Simulation setup.

A layout of MD simulation box used in the work is depicted in Fig. 2. In the layout, a slab of water molecules is enclosed between two silica surfaces (slit-pore setup). The Cartesian z-axis is directed perpendicular to the surfaces and the periodic boundary conditions are applied to the system in three directions. Note that the simulation cell in the setup possesses the inversion symmetry. In order to test the influence of the layout on the computed characteristics, the CRB/NW interface was also modeled in a layout where only one side of water slab was in contact with the surface, whereas the second was exposed to the vacuum (slab setup). Results obtained for the system with the two simulation setups are compared in ESI † .The comparison shows that the characteristics of the system are independent on the setup, while the use of the slit-pore layout permits to enhance the statistical reliability of results.

Liquid phase.

The slab of water was obtained in the following way. First, a system of bulk liquid water was carefully equilibrated in an MD simulation for ambient temperature. Then, a slab of molecules of appropriate dimensions was cut out from the final configuration of the bulk simulation and inserted between the two surfaces of silica layer. The interlayer spacing ∆z (Fig. 2) was adjusted so that the water density in the central part of the slab was equal to the density of bulk liquid water. The thickness of silica layer and ∆z were greater than 2R C , where R C is the cut-off radius of short-range interactions (vide infra), that excluded the direct interaction of molecules adsorbed on the two surfaces across either the liquid phase or the solid. In order to obtain the 1 M NaCl electrolyte solution, the appropriate number of water molecules in the corresponding silica/NW system was at random replaced by sodium and chlorine ions, so that the ions were at least 5 Å distant. The system was then equilibrated in a few ns MD run without adjustment of the interlayer spacing ∆z. Parameters of MD simulation box and characteristics of the systems are reported in Table 1.

Potential models.

Water molecules were described with a flexible SPCFw model by Wu et al. [START_REF] Wu | Flexible simple point-charge water model with improved liquid-state properties[END_REF] with the intramolecular part taken from ref. 38. Nonbonded interactions in the systems were modeled by the sum of Lennard-Jones (12-6) and Coulomb potentials

E (ij) NB (r) = D (ij) 0 R (ij) 0 r 12 -2 R (ij) 0 r 6 + q i q j r (1)
where i and j denote two atoms separated by a distance r. Parameters R 0 , D 0 and q for the intermolecular interactions were those of the SPC water model, [START_REF] Wu | Flexible simple point-charge water model with improved liquid-state properties[END_REF][START_REF] Berendsen | Interaction models for water in relation to protein hydration[END_REF] parameters for the silica atoms were taken from the INTERFACE force field. [START_REF] Emami | Force field and a surface model database for silica to simulate interfacial properties in atomic resolution[END_REF][START_REF] Heinz | Thermodynamically consistent force fields for the assembly of inorganic, organic, and biological nanostructures: The interface force field[END_REF] A recent modeling study of water on a neutral Q 4 silica surface has shown that the force field correctly reproduces the hydrophobic character of the surface. [START_REF] Smirnov | Structure and sum-frequency generation spectra of water on uncharged Q 4 silica surfaces: a molecular dynamics study[END_REF] Parameters for the Na + and Cl -ions were taken from ref. 52. The Lorentz-Berthelot combination rules were used to obtain parameters for unlike atoms. Conventional 1,2 and 1,3 exclusion rules were applied to the nonbonded interactions.

The energies O-H bond-stretching and Si-O-H angle-bending coordinates of the silica surface OH groups were described with a Morse potential and a harmonic angle-bending potential, respectively. The energy parameter of the Morse potential was taken to be equal to that of the water model [START_REF] Smirnov | Structure and sum-frequency generation spectra of water on uncharged Q 4 silica surfaces: a molecular dynamics study[END_REF] and the equilibrium OH bond length was set to 1 Å. The remaining parameters of the potentials were adjusted to mimic the OH bond-stretching vibrational frequency, the Si-O-H angle and angle-bending vibrational frequency of the silica surface OH groups. [START_REF] Tsyganenko | Study of deformation vibrations of hydroxyl groups of a silicon dioxide surface with respect to spectra of combined frequencies[END_REF][START_REF] Morrow | Surface vibrational modes of silanol groups on silica[END_REF][START_REF] Carteret | Vibrational properties of silanol group: From alkylsilanol to small silica cluster: Effects of silicon substituents[END_REF][START_REF] Bandura | Periodic density functional theory study of water adsorption on the α-quartz (101) surface[END_REF] All parameters of the interatomic potentials are listed in Table 2. a Ow and Hw denote oxygen and hydrogen atoms of water; see Fig. 1a for the notation of other atoms.

Computations

3.1 Characteristics of silica/water interface.

The silica/water interfaces were characterized by a number of quantities. First, profiles of relative density ρ * (z) = ρ(z)/ρ 0 , where ρ 0 is the density of bulk liquid water, and of orientational order parameter S 2 (z) for the molecular dipole and the HH vector along the z-axis were computed. [START_REF] Smirnov | A molecular dynamics study of the interaction of water with the external surface of silicalite-1[END_REF][START_REF] Smirnov | Structure and sum-frequency generation spectra of water on uncharged Q 4 silica surfaces: a molecular dynamics study[END_REF] The S 2 (z) profile for a vector is defined as

S 2 (z) = 1 N z N z ∑ i=1 P 2 (cos θ i ) , (2)
where cos θ i is the cosine between the vector and the z-axis, P 2 () is the second-order Legendre polynomial N z is a number of molecules in with CoM z-coordinate from z to z + ∆z. S 2 = -0.5 corresponds to a mean vector orientation perpendicular to the z-axis, S 2 = 1 points to the parallel-to-axis orientation, while S 2 = 0 indicates a random orientation of the vector with respect to the axis. Following Gaigeot et al., [START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF] the network of intermolecular hydrogen bonds was described by a joint probability density distribution Q(r, u, Z ) of the orientation of H-bond donor (D) → acceptor (A) vector and of the DA distance

Q(r, u, Z ) = ∑ i∈Z δ(r -R (i) DA ) • δ(u -cos ψ (i) ) , (3)
where

R (i)
DA and cos ψ (i) are the DA distance and the cosine of angle between the DA vector and the z-axis, respectively, Fig. 3. The sum in (3) runs over donor atoms i with z-coordinate within a region Z. Hydrogen bonds were identified using the geometric criteria. [START_REF] Prada-Gracia | The quest for self-consistency in hydrogen bond definitions[END_REF] Specifically, an H-bond was determined when the DA distance R DA was smaller than 3.3 Å (position of the 1 st minimum in the O-O pair radial distribution function) and the D-H•••A angle was greater than 140 • . The same criteria were used for the analysis of number of hydrogen bonds in different Z regions and for different DA pairs. For an isotropic (centrosymmetric) medium the two-dimensional Q(r, u, Z ) map is characterized by a homogeneous probability density without a preferred DA vector orientation.

In addition to the S 2 (z) order parameter profiles, the orientation of water molecules with respect to the z-axis was characterized by a joint probability density distribution P(u 1 , u 2 , Z )

P(u 1 , u 2 , Z ) = 1 N Z ∑ i∈Z δ(u 1 -cos φ (i) 1 ) • δ(u 2 -cos φ (i) 2 ) (4)
with cos φ

(i)
k being the cosine of angle between the k-th OH bond vector (k = 1, 2) of molecule i and the z-axis, Fig. 3. The normalization factor N Z is given by 1 -1

P(u 1 , u 2 , Z ) du 1 du 2 = 1. (5)
In what follows, the paper uses a relative probability density P(u 1 , u 2 , Z ) = P(u 1 , u 2 , Z )/P(u 1 , u 2 , ∞), where P(u 1 , u 2 , ∞) stands for the probability density distribution (4) in an isotropic medium (taken as bulk liquid water). [START_REF] Willard | Instantaneous liquid interfaces[END_REF] 3.2 Nonlinear spectra calculation.

Nonlinear sum-frequency generation spectrum of water slab was obtained with the time correlation function approach. [START_REF] Morita | A theoretical analysis of the sum frequency generation spectrum of the water surface. II. Time-dependent approach[END_REF] According to the approach, the pqr component of the resonant part (R) of the frequency-dependent nonlinear second-order susceptibility tensor χ (2) (ω) is given by

χ (2),R pqr (ω) = iω k B T ∞ 0 dt e iωt M r (0) • A pq (t) , (6)
where M r and A pq stand for the r and pq elements of the system dipole M and polarizability A, respectively. The computation of the quantities is discussed in ESI † . The work focuses on a spectrum obtained as an average of the xxz and yyz elements of the χ (2),R tensor. This χ

(2),R ⊥ spectrum, where the and ⊥ symbols respectively designate x and y, and z components, can be related to the SFG signal recorded with the ssp polarization combination, that is the sand p-polarized visible and infrared pulses, respectively, and the s-polarized sum-frequency radiation.

Time correlation function in (6) included both intra-and intermolecular contributions with the intermolecular terms limited to neighbor molecules at distances lesser that 5.6 Å. A Hann apodization window was applied to the correlation function prior to performing the Laplace-Fourier transform and the resulting spectrum was smoothed by using a Gaussian filter with the parameter σ equal to 3 cm -1 . The spectra were normalized by the surface area (Table 1). Details on the spectra calculations can be found elsewhere. [START_REF] Smirnov | Structure and sum-frequency generation spectra of water on uncharged Q 4 silica surfaces: a molecular dynamics study[END_REF]

Molecular dynamics simulations.

MD simulations and analysis were performed with an in-house codes. The classical equations of motion were integrated using the velocity Verlet algorithm with the time-step of 0.5 fs. The real space cut-off radius R C for the van der Waals and electrostatic interactions was equal to 10 Å. The electrostatic interactions were handled using a damped shifted-force modification of Wolf method [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] with a damping parameter of 0.2 Å -1 . [START_REF] Fennell | Is the Ewald summation still necessary? Pairwise alternatives to the accepted standard for long-range electrostatics[END_REF] The discontinuity of the short-range energies and forces at the cut-off distance was treated with a shifted-force technique. In the simulations, Si and O atoms were kept fixed at their initial positions, while H atoms of hydroxy groups were allowed to move. To separate the vibrational dynamics of silica surface OH groups from that of water molecules, the mass of silanol H atoms was set to the mass of deuterium.

Reported results are the averages of 50 MD runs performed for temperature of 293 K. Each run consisted of 10 ps equilibration period in a canonical (NVT) ensemble that was followed by a 50 ps NVE stage in which the coordinates and velocities of atoms were saved each 4 fs during the last 40 ps. The final configuration of previous run was used as the starting configuration of next run, while the atomic velocities were anew chosen from the Maxwell-Boltzmann distribution. Temperature in the NVT simulations was controlled via a chain of Nosé-Hoover thermostats. Prior to starting the series of the production runs, the systems were equilibrated in an NVT simulation of at least 2 ns length.

Results and discussion

Structural characteristics of silica/water interfaces

Profiles of water density and S 2 order parameters. Profiles of the relative water density ρ * along the z-axis are presented in Fig. 4a. The profiles illustrate a layered structure of interfacial water with the presence of several layers on the crystalline surfaces. No evident relation between the layered water organization and the surface density of silanol groups can be inferred from the profile analysis. On the other hand, Figure 4a shows that a well defined pattern of layers near the crystalline surfaces blurs into a faded, almost featureless density profile near the surface of amorphous silica and therefore, the surface crystallinity strongly affects the structure of interfacial water. The presence of ions in the liquid phase does not notably influence the ρ * profiles and a slightly smaller density in the central part of water slab can be ascribed to a lesser number of water molecules in the liquid phase of the SiO 2 /ES systems. Figure 4b and Figure 4c display z-profiles of S 2 order parameters (2) for the dipole and HH vector of H 2 O molecules, respectively. At distances corresponding to the first peak in the density profile, negative values of the S 2 order parameter for the dipole indicate that molecules tend to orient their dipoles perpendicularly to the z-axis, whereas the positive S 2 values for the HH vector suggest a favored orientation of the vector towards the surface. The presence of ions in the liquid phase only slightly affects the S 2 dipole parameter while that for the HH vector shows a somewhat higher sensitivity to the composition of aqueous phase. Analysis of the z-profiles in Fig. 4 permits to identify different interfacial regions Z in the liquid phase. The first water layer, in the immediate proximity to the surface, corresponds to the bonded interfacial layer as defined in ref. 33,41. Accordingly, this region is labelled as BIL in Fig. 4a and the layer limit z BIL is placed at the z-coordinate of the first minimum in the density fluctuations. The fluctuations of ρ * decay at the distance of 6 -8 Å from z BIL for all the surfaces studied and, in a somewhat arbitrary way, the region from z BIL to z DL = z BIL + 7 Å is defined as the diffuse layer in the present work. At the z = z DL distance, the S 2 order parameters reach the zero value typical of the isotropic vector orientation. Beyond z DL one finds a bulk zone (BLK), where the structural descriptors attain values characteristic of bulk liquid water.

Ions distribution.

Figure 5 shows distribution of ions in the liquid phase of the Q101/ES and AS/ES systems. Because of their smaller size (and an effective negative surface charge, vide infra) the Na + ions approach the surface closer than large chlorine anions. One observes a slight increase of cations concentration in BIL compared to the DL and BLK regions with a surface density of ions equal to 0.20 and 0.25 Na + per nm 2 for the Q101/ES and AS/ES systems, respectively. Bouhadja and Skelton found no marked accumulation of sodium cations on the Q101 surface in 0.5 M NaCl solution, unless the surface is negatively charged. [START_REF] Bouhadja | Dynamical properties of water and ions at the quartz (101)-water interface at a range of solution conditions: A classical molecular dynamics study[END_REF] From the results reported in that paper for a neutral Q101 surface, one can estimate the surface density of sodium cations as ca. 0.21 nm -2 in a region including two first layers in the Ow density profile. This value well agrees with the Na + surface density obtained in the present study for the system. A density of 0.12 nm -2 for Na cations adsorbed on the Q101 surface was obtained by Kroutil et al. [START_REF] Kroutil | Computer simulations of quartz (101)-water interface over a range of pH values[END_REF] for a 0.3 M NaCl solution. Note that the aforementioned values of the cation surface density correspond to one Na + ion per surface area of 20 Å × 20 Å that can hardly be considered as a continuous layer. Hydrogen bonds analysis. Table 3 reports the number of hydrogen bonds for molecules in the Z regions introduced above. Molecules in BIL form a smaller number of intermolecular H-bonds than the number of H-bonds in bulk neat water N B HB = 3.48. The lack of intermolecular contacts is in part compensated by the interaction with surface hydroxy groups. A water molecule of BIL at the CRB surface donates and accepts an equal number of hydrogen bonds to/from the surface with a slight preference for the interaction with the outer OH groups. A BIL molecule on the Q101 surface donates an almost the same number of H-bonds to both types of silanols, whereas it acts as an H-bond acceptor mainly for the outer OH groups (Table 3). This difference can be explained by a recessed position and specific orientation of the inner hydroxyls on the surface that allow the hydroxyl's H atom to readily interact with oxygen atoms of either neighbor outer silanols or of siloxane bridges. [START_REF] Dellostritto | Effect of ions on H-bond structure and dynamics at the quartz(101)-water interface[END_REF][START_REF] Skelton | Simulations of the quartz(10 11)/water interface: A comparison of classical force fields, ab initio molecular dynamics, and Xray reflectivity experiments[END_REF] The data in Table 3 show that molecules of BIL at the Q001 surface prefer to donate H-bond to the in-plane OH (IN) groups and to accept the H atom of the out-of-plane (OUT) hydroxyls. This is due to the fact that the H atoms of the IN silanols are involved in the intra-silanol H-bonds with neighbor OUT hydroxyls. It is interesting to note that the number of water-surface H-bonds does not correlate with the silanol density, but rather reflects the ability of the surface hydroxyls to form the H-bonds. Thus, the terminal silanols of the CRB surface are the most accessible to the solvent molecules that results in the formation of 1.35 water-surface H-bonds, whereas only 1.03 H-bonds are created between H 2 O molecules and silanols of the Q101 surface, although the silanol density is greater for the latter. Another factor influencing the water-surface interaction is the surface crystallinity. Despite the equal silanol density, the number of water-surface OH-bonds is twice as large for the crystalline CRB surface than for the amorphous one whose corrugated morphology hampers H-bonding between water and a part of surface OH groups. It is also worthy of noting small number of H-bonds, less than 0.1, between H 2 O molecules in BIL and the oxygen atoms of siloxane bridges on the surfaces. As one can expect, the number of such contacts decreases with the increase of OH because of lesser number of surface Si-O-Si bridges exposed to water. Table 3 shows that molecules in the DL and BLK regions have the number of hydrogen bonds equal to N B HB within the statistical uncertainty for all the surfaces and therefore, the local environment of H 2 O molecules typical of bulk liquid water is recovered at distances from the surfaces further than z BIL . These features are common to both silica/neat water and silica/1 M NaCl interfaces. The H-bonds analysis for the SiO 2 /1 M NaCl solution interfaces (Table 3) show that the presence of ions in the liquid phase reduces the number of intermolecular H-bonds in all Z regions. Note that the decrease of the number of H-bonds in the DL and BLK regions is due to the formation of solvation shells around both cations and anions, whereas only the solvation of Na + ions accounts for the change of number of intermolecular H-bonds in BIL. A number of intermolecular H-bonds N HB in an 1 : 1 electrolyte solution of concentration [C] can be estimated as

N HB ≈ N B HB - D C + D A 55.5 [C], (7)
where N B HB is the number of H-bonds in the neat bulk water (vide supra), and D C and D A stand for the number of intermolecular hydrogen bonds broken because of the formation of solvation shell around the cation and anion, respectively. Upper limits of these quantities can be taken as D C = 2N + and D A = N -with N + and N -being the number of water molecules in the first solvation shell of cation and anion, respectively. For the Na + and Cl -ions the MD simulations yield N + = 5.6 and N -= 6.6 and then, (7) gives N HB ≈ 3.16 for an 1 M solution. This lower limit estimate reasonably agrees with the number of intermolecular H-bonds in the DL and BLK regions of the systems.

In contrast to the intermolecular hydrogen bonds, the inspection of Table 3 shows that the number of water-surface H-bonds for the silica/electrolyte solution interfaces is very similar to that in the systems with the neat water. A little influence of ions on the number of water-silanol H-bonds is in line with a small concentration of Na + ions in BIL and suggests that the cations do not notably modify the watersurface interactions, at least at the given salt concentration.

Structure and dynamics of bonded interfacial layer. Two-dimensional probability density maps Q(r, u, Z ) (3) computed for BIL of the silica/water interfaces are displayed in Fig. 6. The maps show a high anisotropy of the orientations of DA vector and are clearly surface-specific. In the case of the CRB/water interfaces, most of intermolecular DA vectors in BIL are directed towards the liquid phase and form an angle greater than 120 • with the z-axis, Fig. 6. Intermolecular H-bonds in BIL on the Q101 surface have the Q(r, u, Z ) distribution in the same limits as that for the CRB surface, but most probable orientation corresponds to a DA vector parallel to the surface, as indicated by the spot at u ≈ 0. The orientation of intermolecular hydrogen bonds in BIL on the Q001 surface is characterized by a relatively narrow distribution that suggests a more ordered network of H-bonds compared to other surfaces. The position of probability density maximum shows that most of the intermolecular DA vectors in BIL on the Q001 surface form angles of ca. 120 • with the z-axis. In contract to the crystalline silica surfaces, the H-bond network in BIL on the amorphous surface has a featureless distribution map with the DA vector/z-axis angle in a wide interval from 40 • to 180 • , Fig. 6. The surface specificity of maps in Fig. 6 implies to that the surface morphology and the surface hydroxy groups play an important role in the formation of network of intermolecular hydrogen bonds by orienting water molecules a particular way. On the other hand, the analysis of Fig. 6 shows that the mean length of DA vector is independent on surface and therefore, the strength of intermolecular H-bonds is not affected by the water-surface interactions. Furthermore, a little or almost no influence of ions on the network of water-water H-bonds in BIL can be inferred from the comparison of Fig. 6a and Fig. 6b.

Figure 7 shows P(u 1 , u 2 , Z ) maps computed for BIL of water on the SiO 2 surfaces. The maps for the CRB surface shows that molecules in BIL tend to keep one OH bond directed to the surface (u 1 ≈ 0.5), whereas the orientation of the second bond varies in a wide range from 45 • to 180 • with respect to the z-axis. The most probable orientation of molecules in BIL corresponds to two OH bonds pointing towards the surface (u 1 , u2) ≈ (0.5, 0.5). Directions of OH bonds of BIL molecules on the Q101 surface are less definite, but a low probability density for u 1 < 0, u 2 < 0 indicates that the molecules keep from having both bonds directed into the liquid phase. The P(u 1 , u 2 , Z ) probability maps for BIL molecules of the Q001 surface have an apparent pattern showing an orientation of water molecules with one OH bond directed towards the surface and the seconds one pointing to the liquid phase and a narrow distribution suggests a high orientational ordering of OH bonds in the layer. Water in BIL on the amorphous surface has the most featureless P(u 1 , u 2 , Z ) probability maps with a pattern implying a paucity of molecular orientations with two OH bonds looking towards the liquid phase, similarly to the Q101/water interface. The comparison of the OH bond orientation maps for the silica/NW and silica/ES systems shows that the presence of ions in the liquid phase does not significantly affect the orientation of molecules in the layer. A small number of ions in BIL likely explains this result.

Figure 8 displays maps of water density in BIL of the crystalline SiO 2 /NW interfaces in the xy plane. In order to even out the amplitude of density fluctuations, the quantity used in the maps construction is ln(1 + ρ * ). As shown above, the presence of ions in the liquid phase does not impact the BIL structure and the maps corresponding to BIL of SiO 2 /ES interfaces can be found in ESI † . Figure 8a reveals the presence of two adsorption sites on the CRB surface, labelled A and B. Inspection of system con-figurations and of the map shows that molecules adsorbed in the site A interact with two OUT and one IN silanols, while molecules in the site B interact with two inner and one outer OH groups of the surface. Similar values of water density in the sites and an almost the same number of H-bonds with the two types of silanols (Table 3) imply a similar strength of water-surface interaction for these sites. Figure 9a displays a configuration of H 2 O molecules in BIL on the CRB surface. The molecules labelled BIL 1 and BIL 2 are in the sites A and B, respectively, and the orientation of BIL molecules accounts the corresponding patterns in Fig. 6 and Fig. 7.

The adsorption site A is also found on the Q101 surface which is characterized by chains of OH groups in the y direction similar to those on the CRB surface, cf. Fig. 8a and Fig. 8b. On the other hand, no site B is apparent in the xy density map for the Q101 surface. As it was mentioned above, the surface morphology likely impedes the interaction of H 2 O molecules with inner hydroxyls that results in the number of OH-bonds between water and the silanols twice as small as that with the outer OH groups (Table 3). The lack of water-surface H-bonds anchoring molecules to the surface is compensated by an increase of intermolecular bonding, and can explain a variety of OH bond orientations in BIL, Fig. 7. The xy density map for water in the first layer on the Q001 surface, Fig. 8c, shows a motif of adsorption sites that can readily be identified with molecules adsorbed between rows of geminal OH groups. Note that the sites form a hexagonal pattern with a site-to-site distance of ca. 4.9 Å. Combining the data in Table 3 and the maps in Fig. 6 and Fig. 7 for the system, one can conclude that an H 2 O molecule in BIL donates OH-bond to an IN silanol, accepts OH-bond from an OUT hydroxyl, and completes the H-bonds network by the interaction with neighbor molecules. Figure 9b displays a configuration of water molecules on the Q001 surface. Molecules labelled BIL 1 and BIL 2 are in the adsorption sites seen in Fig. 8c and the H-bonding and orientation of the molecules explain the maps shown in Fig. 6 and Fig. 7. The ρ * (z) profile of BIL on the surface features two subpeaks, Fig. 4a. Analysis of the positions of BIL molecules shows that the BIL 1 and BIL 2 molecules in Fig. 9b are at the origin of the subpeak in closest to the surface, whereas the molecule BIL 3 accounts for the second, more distant subpeak in the density profile.

The mobility of water molecules in BIL was characterized by computing the parallel-to-surface diffusion coefficient (D xy ) and the rotational relaxation time of molecular dipole (τ 4 reports the characteristics obtained for the crystalline SiO 2 /water interfaces and compares them to values typical of bulk liquid water. One sees that the translational and rotational mobility of H 2 O in BIL is strongly hindered. No obvious correlation of these quantities with the surface density of OH groups can be inferred from the data. On the other hand, one can notice that the most mobile H 2 O BIL molecules are on the Q101 surface, which has the smallest density of adsorption sites, compared with the CRB and Q001 surfaces (Fig. 8). A small dipole relaxation time computed for BIL water on the Q101 surface correlates with a wide distribution of molecules' orientation in Fig. 7, whereas the slowest rotational dynamics of H 2 O in BIL on the Q001 surface is in line with the maps in Fig. 6 and Fig. 7. The presence of ions in the liquid phase has a little influence on the diffusion coeffient and mainly affects the dipole relaxation time of water in BIL. The structure and dynamics of neat water on the (100), ( 101) and (001) surfaces of α-quartz was investigated by classical MD simulations by Notman and Walsh. [START_REF] Notman | Molecular dynamics studies of the interactions of water and amino acid analogues with quartz surfaces[END_REF] The authors observed the formation of at least two layers of ordered water on the surfaces with a more dense first layer and the calculations of molecular orientations revealed that molecules of the first layer point one H atom towards the surface. Notman and Walsh have concluded that positions of H 2 O molecules on the surface are dominated by the positions of hydroxy groups. Their calculations showed a more pronounced lateral ordering of water on the Q101 surface that was explained by the interaction of water with the oxygen atoms of Si-O-Si bridges more exposed to water on the (101) surface, compared to the other surfaces. [START_REF] Notman | Molecular dynamics studies of the interactions of water and amino acid analogues with quartz surfaces[END_REF] Periodic quantum-chemical calculations of water on the Q101 surface by Bandura et al. [START_REF] Bandura | Periodic density functional theory study of water adsorption on the α-quartz (101) surface[END_REF] showed that an isolated molecule favored to interact with an outer silanol and the adsorption was stabilized by three H-bonds of which the molecule donated two and accepted one. At a monolayer coverage, water formed a regular network with each molecule H-bonded to its three neighbors and to a surface hydroxyl. The authors found two types of H 2 O molecules interacting with the surface in different ways. Molecules of the first type pointed one of OH bonds down to an IN silanol and molecules of the second type lied nearly parallel to the surface and accepted an H-bond from an OUT hydroxyl. In line with the findings by Notman and Walsh, classical and DFT-based MD studies of the Q101/NW interface showed the formation of at least two layers of ordered water and revealed that both the outer and inner hydroxy groups play an important role in the water adsorption. [START_REF] Kroutil | Computer simulations of quartz (101)-water interface over a range of pH values[END_REF][START_REF] Skelton | Simulations of the quartz(10 11)/water interface: A comparison of classical force fields, ab initio molecular dynamics, and Xray reflectivity experiments[END_REF] Calculations of lateral water density in the first water layer on the Q101 surface have shown the presence of two adsorption sites. [START_REF] Kroutil | Computer simulations of quartz (101)-water interface over a range of pH values[END_REF] One of these sites closely matches the site A in Fig. 8b, whereas the second one, nearby the inner silanols, is not found in the present calculations (see ESI † ).

(2) µ ). † Table
The interaction of neat water with a hydroxylated (001) surface of α-quartz was investigated in a number of works [START_REF] Notman | Molecular dynamics studies of the interactions of water and amino acid analogues with quartz surfaces[END_REF][START_REF] Du | Molecular dynamics simulations of hydration, dissolution and nucleation processes at the α-quartz (0001) surface in liquid water[END_REF][START_REF] Yang | Water adsorption on hydroxylated α-quartz (0001) surfaces: From monomer to flat bilayer[END_REF][START_REF] Chen | Water thin film-silica interaction on α-quartz (0001) surfaces[END_REF][START_REF] Sulpizi | The silica-water interface: How the silanols determine the surface acidity and modulate the water properties[END_REF][START_REF] Gaigeot | Oxide/water interfaces: how the surface chemistry modifies interfacial water properties[END_REF][START_REF] Liu | Density functional theory study of water molecule adsorption on the α-quartz (001) surface with and without the presence of Na + , Mg 2+ , and Ca 2+[END_REF] . A classical MD study by Du and de Leeuw showed that that surface OH groups closely resemble water molecules and that the intermolecular structure in a region near the surface is similar to that of bulk liquid water. [START_REF] Du | Molecular dynamics simulations of hydration, dissolution and nucleation processes at the α-quartz (0001) surface in liquid water[END_REF] With the use of periodic DFT calculations, Yang and Wang [START_REF] Yang | Water adsorption on hydroxylated α-quartz (0001) surfaces: From monomer to flat bilayer[END_REF] obtained that an H 2 O molecule adsorbs in a bridging position between a pair of geminal surface OH groups; the same adsorption position with the molecule donating proton to an IN hydroxy group and accepting proton from an OUT silanol was found in a recent study by Liu et al. [START_REF] Liu | Density functional theory study of water molecule adsorption on the α-quartz (001) surface with and without the presence of Na + , Mg 2+ , and Ca 2+[END_REF] Contrary to the result by Du and de Leeuw, Yang and Wang obtained that increasing the surface coverage to a monolayer resulted in the formation of an ice-like structure. They noted that the ordering of interfacial water layer can be promoted by the compatibility of Q001 surface lattice constant with the a lattice parameter of I h ice. Quantum-chemical calculations by Chen et al. [START_REF] Chen | Water thin film-silica interaction on α-quartz (0001) surfaces[END_REF] also revealed that water in the first adsorbed layer forms a stable structure that remains unaffected by the adsorption of subsequent layers. When three and more layers were added on the surface, water structure consisted of the first layer and a bilayer made of the two next layers. Such a structural organization of the quartz/water interface was found to effectively shield the interaction between the surface and more distant water layers. In a DFTbased MD study, Gaigeot et al. [START_REF] Sulpizi | The silica-water interface: How the silanols determine the surface acidity and modulate the water properties[END_REF][START_REF] Gaigeot | Oxide/water interfaces: how the surface chemistry modifies interfacial water properties[END_REF] observed two types of water molecules in the first adsorbed layer. Molecules of the first type donated H-bond to IN silanols and the structural characteristics of the Hbonded complex were similar to those in liquid water. Molecules of the second type accepted H-bond from OUT hydroxyls and the complex had characteristics close to those in ice.

Influence of ions on the water behavior on neutral surfaces of α-quartz was investigated in ref. 22,23,25,30. DFT-based ab initio MD simulations by Kubicki and co-workers [START_REF] Dellostritto | Effect of ions on H-bond structure and dynamics at the quartz(101)-water interface[END_REF] showed that the Cl -ions adsorb on the Q101 surface in inner-sphere configuration, while the Na cations are less strongly bound to the surface. The authors found that when the cations were adsorbed on the surface in inner-sphere configuration they tended to increase ordering of interfacial water and to slow down the dynamics. In contrast, the presence both Na + and Cl -in the liquid phase increased the rotational freedom of H 2 O molecules and reduced the H-bond lifetime. With the use of DFT-MD simulations, Pfeiffer-Laplaud et Gaigeot [START_REF] Pfeiffer-Laplaud | Adsorption of singly charged ions at the hydroxylated (0001) α-quartz/water interface[END_REF] found that the Na + ions adsorb on the Q001 surface as inner-sphere complex, whereas the adsorption of Cl -occurs in a form between the inner-/outer-sphere configuration and the anions are 1.4 Å more distant from the surface than the cations. A discrepancy with the results of ref. 25 concerning the configuration of adsorbed ions is likely due to the different density of OH groups on the Q101 and Q001 surfaces. The authors of ref. 23 obtained that the solvation shells around the adsorbed cations and anions, where the silanol oxygens replaced those of water, were very similar to the solvation shells in bulk liquid water. The adsorption of Na cations was obtained to weaken the intra-silanol H-bonds and caused a dynamical reorientation of OH-groups. The investigations reported in ref. 22,30 mainly focused on the effect of surface charge on the structure and dynamics of the Q101/water interface.

Results of the present calculations generally agree with the outcome of the previous studies. Thus, Table 3 shows that a water molecule in BIL on the Q101 surface forms in average 2.3 H-bonds with neighbors and one H-bond with surface hydroxyls, and that both the OUT and IN silanols are of importance for the structural organization of BIL. This finding is in accord with results by Bandura et al. [START_REF] Bandura | Periodic density functional theory study of water adsorption on the α-quartz (101) surface[END_REF] . The position and orientation of H 2 O molecules in BIL on the Q001 surface nicely agrees with the data of ref. 65,69. In line with the outcome of previous studies, [START_REF] Yang | Water adsorption on hydroxylated α-quartz (0001) surfaces: From monomer to flat bilayer[END_REF][START_REF] Chen | Water thin film-silica interaction on α-quartz (0001) surfaces[END_REF] the data for this surface presented in Fig. 6, Fig. 7 and Fig. 8c suggest that the structure of the hydroxyl cover promotes the formation of a stable network of water-silanol and intermolecular hydrogen bonds. As pointed by Yang and Wang [START_REF] Yang | Water adsorption on hydroxylated α-quartz (0001) surfaces: From monomer to flat bilayer[END_REF] , the site-to-site distance commensurate with structural parameters of hexagonal ice can play an important role with this respect, cf. the hexagonal pattern of adsorption sites in Fig. 8c. The Q(r, u, Z ) maps for the Q001 surface in Fig. 6 qualitatively agree with such maps computed by Pezzotti and co-workers using ab initio MD. [START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF] A close inspection of the maps in Fig. 7 shows that the presence of ions slightly reduces the maximum value of the probability density that can point to weakening the H-bond network in BIL, in agreement with the results of ref. 23,25. Some discrepancies with results of previous studies can be explained by differences in computational methodologies. Thus, the CLAYFF force field 70 used in ref. 22 was shown to yield a more hydrophilic silica surface, compared to the INTERFACE potential model. [START_REF] Smirnov | Structure and sum-frequency generation spectra of water on uncharged Q 4 silica surfaces: a molecular dynamics study[END_REF] This can account for the presence of the second adsorption site found on the Q101 surface in that work, while such a site is indiscernible in Fig. 8b. It is also worth of noting that the thermal motion explicitely considered in MD simulations and the difference in system size, often impede the comparison of results of classical MD and static quantum-chemical calculations. On the other hand, differences in system size and/or length of MD trajectories hamper a detailed comparison with ab initio MD simulation results, which are also sensitive to the choice of exchange-correlation functional and to the treatment of dispersion interactions.

Structure of diffuse layer. The Q(r, u, Z ) maps for molecules of the diffuse layer of the silica/water interfaces are presented in Fig. 10. Compared to the maps for BIL in Fig. 6, the DL maps show less featured patterns. However, their inspection clearly points to that the isotropic orientation of DA vector typical of bulk liquid water is not yet completely recovered in the region, especially for the SiO 2 /NW systems. Greater Q(r, u, Z ) values at positive u point to a preferred orientation of the vector towards the surfaces. The tendency is particularly noticeable in the case of the (001) α-quartz surface for which the Q(r, u, Z ) map has a well defined density spot at u ≈ 0.5, Fig. 10. Additional calculations for the crystalline silica/water interfaces have shown that the anisotropy in the patterns in Fig. 10a mainly results from molecules of the second layer (Fig. 4), i.e. the layer next to BIL (ESI † ) Comparison of the maps computed for the neat water (Fig. 10a) and electrolyte solution (Fig. 10b) indicates that ions present in the liquid phase yield the distribution of DA vector orientations more isotropic, compared to that for the silica/neat water interfaces.

The P(u 1 , u 2 , Z ) maps computed for molecules in DL of the silica/water interfaces are shown in Fig. 11. Regardless the composition of aqueous phase, the maps reveal a significantly more uniform distribution of OH bonds orientation in DL compared to BIL, Fig. 7. For the SiO 2 /NW systems, a higher probability density at (u 1 , u 2 ) ≈ (0.5, 0.5) indicates that a part of molecules in DL tend to orient their OH bonds towards the surface. The trend is more apparent in the maps for DL on the Q101 and Q001 surfaces, and the preference in the OH bonds orientation towards the surface is noticeable even for DL on the amorphous silica surface, Fig. 11a. Once more, auxiliary calculations have revealed that this feature is primarily due to DL molecules in the region adjacent to BIL (ESI † ). The presence of ions in the liquid phase results in a more isotropic orientation of the OH bonds in DL for all systems, Fig. 11b, that is in line with the results for the H-bonds network, Fig. 10. Moreover, comparing Fig. 11a and Fig. 11b one can notice a slight tendency of OH bonds to be directed outwards the the surface in DL of the electrolyte solution.

Summarizing the above data, the analysis of OH bond orientations in the diffuse layer of the silica/neat water interfaces shows that molecules tend to direct their OH bonds towards the surface, Fig. 11a. The preference in the orientation is particularly apparent for molecules in the vicinity of BIL. This peculiarity can be explained by the formation of more stable hydrogen bonds between the DL and BIL molecules, compared to bonds between molecules in the diffuse layer. The stability of the BIL-DL bonding, which also accounts for the anisotropy of the DA vector orientation (Fig. 10a), likely arises from a low translational and rotational mobility of BIL molecules. [START_REF] Bouhadja | Dynamical properties of water and ions at the quartz (101)-water interface at a range of solution conditions: A classical molecular dynamics study[END_REF][START_REF] Notman | Molecular dynamics studies of the interactions of water and amino acid analogues with quartz surfaces[END_REF][START_REF] Argyris | Molecular structure and dynamics in thin water films at the silica and graphite surfaces[END_REF][START_REF] -V. Castrillón | Evolution from surfaceinfluenced to bulk-like dynamics in nanoscopically confined water[END_REF][START_REF] Gupta | Dynamics and vibrational spectroscopy of water at hydroxylated silica surfaces[END_REF] This effect inextricably related to the hydrophilic water-surface interactions and it is especially marked for molecules of DL on the Q001 surface because of a highly ordered and low-mobile BIL structure. In contrast to BIL, the characteristics of water in the diffuse layer are markedly influenced by the composition of liquid phase. The presence of ions in the solution results in a more isotropic DL, compared to the neat water and this finding is in line results of DFT-based MD study of Q001/1 M KCl interface [START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF] . The effect of ions can be due to the formation of labile solvation shells that disrupt intermolecular H-bonds and thus, the anisotropy of molecular orientations in DL. [START_REF] Dellostritto | Effect of ions on H-bond structure and dynamics at the quartz(101)-water interface[END_REF] Anticipating the following section, a screening of the surface electric field by ions can be another reason explaining the differences in DL structure of the silica/neat water and silica/electrolyte solution interfaces.

Nonlinear spectra.

Prior to discussing results of the nonlinear spectra calculations of the silica/water interfaces, it is worth to outline the interpretation of Im[χ ) ssp ] spectrum of water/vapor interface which is commonly considered as a reference system. Experimental Im[χ (2) ssp ] spectra of the interface are characterized by a positive peak with maximum at ca. 3700 cm -1 and a negative band centered at 3450 cm -1 . [START_REF] Ji | Characterization of vibrational resonances of watervapor interfaces by phase-sensitive sum-frequency spectroscopy[END_REF][START_REF] Nihonyanagi | Accurate determination of complex χ (2) spectrum of the air/water interface[END_REF] It was shown that Im[χ (2) ssp ] ∝ cos θ, where θ stands for the angle between the transition dipole moment of vibrational mode and the z-axis. [START_REF] Morita | A theoretical analysis of the sum frequency generation spectrum of the water surface. II. Time-dependent approach[END_REF][START_REF] Morita | A theoretical analysis of the sum frequency generation spectrum of the water surface[END_REF][START_REF] Morita | Theory of Sum Frequency Generation Spectroscopy[END_REF] For an H 2 O molecule at the interface with one OH bond pointing to the vacuum (free OH) and the second one to the liquid phase (bonded OH), directions of the transition dipole moment of two O-H stretching modes approximately correspond to the bond directions, because of the decoupling of the bond-stretching vibrations. Then, taking the positive direction of the z-axis outwards the surface, the transition dipole moment of the high-frequency mode due to the free OH bond has a positive cosine with the axis and such bonds are at the origin of the positive band in the Im[χ (2) ssp ] spectrum of water/vapor interface. Contrarily, the cosine of angle between the transition dipole moment of H-bonded OH and the z-axis is negative that accounts for the negative low-frequency band in the spectrum. The intensity of the bands is controlled by the transition dipole moment magnitude and by the angle θ that both depend on the intermolecular interactions. A very recent review by Tang et al. [START_REF] Tang | Molecular structure and modeling of water-air and ice-air interfaces monitored by sum-frequency generation[END_REF] summarizes the experimental and theoretical investigations of this and related systems.

Regarding the silica/water interface, the presence of silanols on the hydrophilic SiO 2 surface greatly complicates the matter because, in addition to the intermolecular H-bonding, H 2 O molecules can form hydrogen bonds of different strength with the surface. The structure of interfacial water layer is then controlled by the balance between the molecule-molecule and molecule-surface interactions depending on a number of factors such as the surface chemistry, morphology, distance from the surface, etc. Consequently, VSFG spectrum of the interface will consist of contributions arising from groups of molecules in different bonding configurations and differently oriented with respect to the surface.

Nonlinear spectrum of bonded interfacial layer. Figure 12 (upper row) presents the Im[χ

(2),R ⊥ ] spectra computed for BIL on the silica surfaces. The spectra have much in common and feature a negative band centered at ca. 3250 cm -1 and a positive band extending from 3400 to 3600 cm -1 . The shape of the negative band is almost identical for all the surfaces, whereas that of the positive band is manifestly surface-specific. Interestingly, the intensity of the negative band is very similar for all the crystalline surfaces. Making use of the interpretation of the Im[χ (2) ssp ] spectrum of the water/vapor interface, the negative band can be assigned to OH bonds directed towards the liquid phase and its universal shape and position can be explained by a sameness of intermolecular hydrogen bonds. The positive band can be attributed to OH bonds of BIL molecules pointing to the surface and the surface-dependent shape of the band suggests that the band bears information about water-surface interactions. The BIL spectra in Figure 12 show that increasing OH causes a decrease of band intensity from the high-frequency side and a downward shift of band maximum. Such a behavior can be related to the formation of more regular network of water-surface H-bonds while OH increases. The shape of the positive band can also be correlated to the surface roughness. Indeed, the band in the spectrum of the amorphous silica/water interface has lower intensity at the low-frequency side than at the high-frequency one, whereas the contrary is observed in the spectrum of BIL on the CRB surface, although the surfaces have the equal density of silanol groups. This difference of shapes can be attributed to a different number of water-silanol OH bonds. Table 3 shows that H 2 O molecules in BIL on the amorphous surface form two times less number of H-bonds with the surface than molecules of BIL on the surface of cristobalite and consequently, the smaller number of H 2 O molecules H-bonded to silanols of the amorphous surface accounts for the smaller intensity of the positive band from the low-frequency side. The spectrum of BIL on the amorphous silica surface deserves a particular comment. Figure 12 shows that the spectrum has a lower intensity than that for the crystalline surfaces. This can be due to two reasons. First, BIL layer on the amorphous surface is characterized by a lower orientational ordering compared to the crystalline surfaces (see Fig. 6 and Fig. 7) and one can then expect smaller χ (2) values for the amorphous silica/water interface. Second, because of the surface roughness, the first layer of interfacial water, as identified in z-profiles in Fig. 4, is formed by groups of molecules unequally distant from the surface. Depending on how far is the group from the surface, it will contribute to the spectrum with either BILlike or DL-like signature. Consequently, the spectrum of the first layer in the water density profile for the AS surface can be "contaminated" by a DL spectrum. It is also worth noting that the positive band in the BIL spectra in Fig. 12 has no feature from the high-frequency side at 3650-3680 cm -1 . A band at this position was observed in several works and assigned to surface hydroxyls [START_REF] Dalstein | The elusive silica/water interface: isolated silanols under water as revealed by vibrational sum frequency spectroscopy[END_REF] or to water weakly interacting with the oxygens of Si-O-Si bridges. [START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF][START_REF] Isaienko | Hydrophobicity of hydroxylated amorphous fused silica surfaces[END_REF] Calculations by Cyran et al. [START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF] have shown that such weakly adsorbed H 2 O molecules have indeed a vibrational mode with frequency in this spectral region. However, in order the band to be clearly visible in the spectrum the surface has to have large hydrophobic patches with a large number of bridging oxygens exposed to water. [START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF] For the systems studied in the present work, the H-bond analysis (Table 3) shows a very small number of water-siloxane bridge contacts and therefore, their spectral signature in the BIL spectra in Fig. 12 is submerged by the responses from the intermolecular and water-silanol interactions. Results of a recent modelling study of the VSFG spectra of water on a hydropobic silica surface [START_REF] Smirnov | Structure and sum-frequency generation spectra of water on uncharged Q 4 silica surfaces: a molecular dynamics study[END_REF] suggests that the band due the H 2 O molecules weakly interacting with the surface would appear in the spectra as a shoulder at ca. 3610 cm -1 .

The presence of ions has a little effect on the spectra of BIL, Fig. 12. This is in line with the results on a small influence of ions on the intermolecular structure in the layer. The most marked difference between the BIL spectra for the silica/neat water and silica/electrolyte solution interfaces is observed for the amorphous silica surface and it is likely due to the structural heterogeneity of BIL because of the surface roughness, as discussed above. One should mention that a qualitatively similar VSFG spectrum of BIL with a positive band above 3350 cm -1 and a negative band bellow this wavenumber, was obtained by Shen and co-workers in their experimental study of the lignoceric acid monolayer/water interface at pH = 2.5. [START_REF] Wen | Unveiling microscopic structures of charged water interfaces by surface-specific vibrational spectroscopy[END_REF] Nonlinear spectrum of diffuse layer. Figure 12 (middle row) displays Im[χ (2),R ⊥ ] spectra of the DL region of silica/water interfaces. The spectra of the SiO 2 /NW systems are characterized by an intense positive band from 3100 to 3400 cm -1 and a weak negative band from 3400 to 3650 cm -1 . The shape of the spectra implies the presence of sub-bands. However, as it is not clear whether these features are due to a specific structural organization of DL or stem from a limited sampling, the following discussion focuses on the main bands in the DL spectra.

The positive band changes its sign when the aqueous phase contains electrolyte ions. The shape of the DL spectra is very similar for the silica/water interfaces of each series and a minimum (maximum) of the band is at ca. 3250 cm -1 regardless the surface. A similar effect of ions on the sign of DL spectrum was obtained by Gaigeot et al. [START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF] in their modeling study of Q001/NW and Q001/KCl solution interfaces and the DFT-MD simulations of silica/water interfaces showed that DL spectra are sensitive to such factors as the surface charge and the formation of electric double layer. [START_REF] Pezzotti | What the Diffuse Layer (DL) Reveals in Non-Linear SFG Spectroscopy[END_REF][START_REF] Pezzotti | Deconvolution of BIL-SFG and DL-SFG spectroscopic signals reveals order/disorder of water at the elusive aqueous silica interface[END_REF] Although the present study deals with neutral silica surfaces, it was shown that a fully protonated (001) surface of α-quartz possesses an effective negative charge [START_REF] Pfeiffer-Laplaud | Adsorption of singly charged ions at the hydroxylated (0001) α-quartz/water interface[END_REF] and the same could therefore be supposed for other silica surfaces.

According to the scheme suggested by Shen et al. [START_REF] Wen | Unveiling microscopic structures of charged water interfaces by surface-specific vibrational spectroscopy[END_REF] to separate the BIL and DL contributions of charged solid/water interface, the pqr component of the nonlinear susceptibility of interfacial water slab χ (2) pqr,S can be written as χ

(2)

pqr,S = χ (2) pqr,BIL + χ (2) pqr,DL (8)
where χ

(2)

pqr,BIL and χ

(2)
pqr,DL are the second-order nonlinear susceptibilities of the BIL and DL regions, respectively. The second term in ( 8) is given by

χ (2) pqr,DL = χ (3) pqrz z BIL -∞ E z (ζ) dζ, (9) with χ (3)
pqr being the third-order nonlinear susceptibility (essentially the bulk quantity) and E z (ζ) standing for the z-component of electric field existing across the interface. The DL contribution χ (2) pqr,DL can be recast to

χ (2) pqr,DL = -χ (3) pqrz Φ(z BIL ), (10)
where Φ(z) is the surface electrostatic potential

Φ(z) = - z -∞ E z (ζ) dζ. (11)
One sees that spectrum of the DL regions depends on the magnitude and sign of the potential at the frontier between DL and BIL.

Figure 13 shows the electrostatic potential (11) computed for the silica/water interfaces. The electric field E z (ζ) in (11) was obtained as

E z (ζ) = E z (ζ) , (12)
where E z (ζ) is the z-component of electric field in a point with the z-coordinate ζ and E z (ζ) includes contributions due to both the solid and the liquid phase. The SiO 2 /NW interfaces are characterized by a negative potential that, according to (10), accounts for the positive band in the DL spectra in Fig. 12. The presence of ions in the liquid phase yields the potential positive and changes the sign of χ

If the χ (2)

pqr,DL spectrum would be due to the field-induced contribution only, the ratio χ (2) pqr,DL /Φ(z BIL ) would yield a spectrum identical for all the interfaces as the ratio gives the third-order nonlinear susceptibility χ (3) pqrz (10) which is the bulk quantity. Figure 14 displays -χ (2) DL /Φ(z BIL ) spectra of the silica/water interfaces. The spectra are in a qualitative agreement with the χ (3) sspz spectrum obtained by Shen et al. [START_REF] Wen | Unveiling microscopic structures of charged water interfaces by surface-specific vibrational spectroscopy[END_REF] using the spectra of charged surface/water interface at different pH and with

χ (3)
sspz spectra computed in ref. 29,34. However, a notable variation of the spectral intensity in Fig. 14, especially for the SiO 2 /NW interfaces, implies that there might be another contribution to the χ (2) susceptibility of the DL region, in addition to that by surface electric field.

This additional contribution to χ

DL can result from the fact that H 2 O molecules of diffuse layer "see" the bonded interfacial layer as a hydrophilic surface. One can then expect that even in the absence of surface field, DL will not be centrosymmetric and will therefore contribute to the χ (2) S spectrum of the interfacial water slab. Indeed, results of model simulations † show the existence of such an intrinsic DL spectrum which is characterized by a negative band in the region 3100 to 3700 cm -1 . The intensity of the band depends on the strength of surface-water interaction, and thus on the mobility of BIL molecules. The intrinsic spectrum can then explain the negative tail above 3400 cm -1 in the DL spectra in Fig. 12 and the intensity variation in the χ (2) pqr,DL /Φ(z BIL ) spectra in Fig. 14. The influence of ions on characteristics of interfacial water is commonly explained by the formation of Stern layer that shields more distant slab parts from the surface field. While the effect is expected to play an important role for charged surfaces and high electrolyte concentrations, such a marked influence of ions in the case of formally neutral surfaces is somewhat unexpected. Indeed, the calculations give a small concentration of Na cations in BIL (one Na + per 300-400 Å 2

⊥,DL /Φ(z BIL ) for the silica/neat water (a) and silica/1 M NaCl solution (b) interfaces. Filled shape shows the spectrum of χ viewed as the formation of continuous layer. On the other hand, the computed DL spectra unequivocally show that even such a small concentration has a drastic effect on the DL spectra. The formation of solvation shells around the ions in DL resulting in an increase of orientational disorder in the layer, can be another plausible reason, which in addition to the electrostatic screening, accounts for the influence of ions in the aqueous phase on the nonlinear spectra of the diffuse layer. The total χ (2) DL spectrum arises from a combination of all mechanisms and an unambiguous separation of different contributions is nontrivial task.

Nonlinear spectrum of total interfacial region. Figure 12 (lower row) presents the χ (2) spectra of the total interfacial region, i.e. both the BIL and DL zones. As these spectra are the sum of the corresponding contributions, one can readily see a cancellation effect. Due to the similar position and shape but the opposite sign of the low-frequency band in the BIL and DL spectra, the effect is especially notable for the crystalline/neat water interfaces in the spectral region below 3400 cm -1 . With this respect, the spectrum of the Q001/NW interface provides an ultimate example. Since the presence of ions in the liquid phase changes the sign of DL spectrum, the spectra of SiO 2 /ES interfaces are markedly different from those of the SiO 2 /NW systems. The data presented in Fig. 12 emphasize that the interpretation of the SFG spectra without the knowledge of separate contributions from different interfacial regions might be erroneous. Even for such a seemingly simple system as the water/vapor interface, Skinner and co-workers have concluded that the assignment of features in VSFG spectrum in analogy with the interpretation of IR and Raman spectra, can be dangerous because the χ (2) spectrum is a sum of signed contributions coming from different groups of molecules. [START_REF] Pieniazek | Interpretation of the water surface vibrational sumfrequency spectrum[END_REF] Roy and Hore came to the same conclusion after a computational study of nonlinear spectra of water on a model hydrophilic surface. [START_REF] Roy | Simulated structure and nonlinear vibrational spectra of water next to hydrophobic and hydrophilic solid surfaces[END_REF] Experimental VSFG investigations of silica/water interfaces have mostly reported the absolute spectra (compilation of experimental results can be found in ref. 18), while the phase-sensitive (heterodyne) VSFG experiments are much less frequent. [START_REF] Ostroverkhov | New information on water interfacial structure revealed by phase-sensitive surface spectroscopy[END_REF][START_REF] Myalitsin | Water structure at the buried silica/aqueous interface studied by heterodyne-detected vibrational sum-frequency generation[END_REF][START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF][START_REF] S.-H. Urashima | The topmost water structure at a charged silica/aqueous interface revealed by heterodyne-detected vibrational sum frequency generation spectroscopy[END_REF] Figure 15 compares spectra calculated for the amorphous silica/water interface to the heterodyne and homodyne VSFG spectra measured for the system at pH = 2.1. [START_REF] Myalitsin | Water structure at the buried silica/aqueous interface studied by heterodyne-detected vibrational sum-frequency generation[END_REF] One sees that the theoretical spectra of the SiO 2 /NW interface strikingly differ from the experimental counterparts. On the other hand, the spectra of the SiO 2 /1 M NaCl solution are in a good agreement with the measurements. Note that the aqueous phase in the experiment contained a buffer (phosphoric acid and phosphate salt) to keep the ionic strength at a constant value and therefore, the AS/ES interface better corresponds to the experimental conditions than the AS/NW system. Figure 15 provides a convincing evidence that modeling the nonlinear spectra of mineral/water interface has to explicitely take into account the chemical composition of the aqueous phase. Differences be- tween the shape of the theoretical and experimental spectra in Fig. 15, particularly in the low-frequency part, are likely due to limits of classical water model. ⊥ ] spectra of BIL of the silica/water interfaces including the spectral region of the surface OD groups. The inset zooms the zone from 2400 cm -1 to 2700 cm -1 .

Manifestation of surface hydroxy groups in SFG spectra of silica/water interface is a subject of ongoing debate. [START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF][START_REF] Dalstein | The elusive silica/water interface: isolated silanols under water as revealed by vibrational sum frequency spectroscopy[END_REF] As it was mentioned above, the present calculations used the deuterium mass for the H atoms of surface silanols in order to separate the vibrational dynamics of the OH groups and adsorbed water molecules. Figure 16 shows the χ

(2)
⊥ spectra of BIL of the SiO 2 /NW systems including the spectral region of the O-D stretching vibrations. One sees that the spectra feature a low intensity negative band at 2580 cm -1 arising from the polarization of H 2 O molecules by the deuterated silanols. The calculations show that the corresponding contribution induced by OH groups would appear at ca. 3650 cm -1 , where it will interfere with an intrinsic band due to the O-H stretching vibration of surface hydroxyls with a sign depending on the angle between the OH bonds and the normal to the surface. Dalstein and co-workers [START_REF] Dalstein | The elusive silica/water interface: isolated silanols under water as revealed by vibrational sum frequency spectroscopy[END_REF] found that the silanol band has the negative sign that means that both the induced and inherent contributions will add up. It is noteworthy that the position of the OH band can also overlap with the band due to molecules interacting with hydrophobic parts of silica surface. [START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF][START_REF] Isaienko | Hydrophobicity of hydroxylated amorphous fused silica surfaces[END_REF] Phase-sensitive VSFG spectra of silica/water interface have a positive band in this spectral region which, as it was discussed above, has been attributed to H 2 O molecules weakly interacting with the surface. [START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF] Additional experimental and theoretical investigations are necessary to disentangle the silanol and water contributions in this region.

Rationalization of nonlinear spectrum of silica/water interface. With the above results in hands, one can now attempt to rationalize the nonlinear VSFG spectra of the silica/water interface at low pH values. The contributions to the spectrum identified in the work are depicted in Fig. 17. The spectrum is a sum of BIL and DL contributions, and of a contribution arising from the surface silanol groups. The BIL spectrum is characterized by a positive and negative bands above and below 3400 cm -1 , respectively. The position of the maximum of the positive band reflects the density of surface OH groups and the band shape says about the surface crystallinity. A positive peak from the high-frequency side of the BIL band indicates the presence of hydrophobic patches on the surface. [START_REF] Cyran | Molecular hydrophobicity at a macroscopically hydrophilic surface[END_REF][START_REF] Isaienko | Hydrophobicity of hydroxylated amorphous fused silica surfaces[END_REF] A negative signal of silica surface OH groups can attenuate the signal of the weakly interacting H 2 O molecules in this spectral region. [START_REF] Dalstein | The elusive silica/water interface: isolated silanols under water as revealed by vibrational sum frequency spectroscopy[END_REF] The negative band below 3400 cm -1 in the spectrum is to a large extent insensitive to characteristics of the underlying surface as it arises from the intermolecular interactions of BIL water with molecules of diffuse layer. The BIL spectrum is mainly insensitive to the chemical composition of the liquid phase. Spectrum of diffuse layer is not less intriguing. [START_REF] Pezzotti | What the Diffuse Layer (DL) Reveals in Non-Linear SFG Spectroscopy[END_REF] Thus, the spectrum has an intrinsic part due to the fact that the bonded interfacial layer plays a role of hydrophilic surface for DL molecules. The intrinsic DL spectrum has the negative sign and its intensity and shape are related to the mobility and ordering of molecules in BIL. Then, the electrostatic potential across the interface make the DL layer spectroscopically active via the third-order bulk contribution. Sign of the χ (3) contribution is directly related to the sign of the potential at the DL/BIL frontier and the DL spectrum was also found to provide information on the electric double layer formation. [START_REF] Pezzotti | What the Diffuse Layer (DL) Reveals in Non-Linear SFG Spectroscopy[END_REF] Furthermore, the spectral response of diffuse layer is sensitive to the presence of solvated ions that disrupt intermolecular H-bonds and make the layer structure more isotropic. All these mechanisms contribute to the spectral response of DL and account for the dependence of the SGF spectrum of silica/water interface on the surface charge and composition of the liquid phase. These considerations agree with the conclusions by Darlington at al. [START_REF] Darlington | Separating the pH-dependent behavior of water in the Stern and diffuse layers with varying salt concentration[END_REF] who found that the highfrequency part of VSFG spectrum probes BIL region (Stern layer) of silica/water interface, whereas the spectral part below 3400 cm -1 reflects the spectral response of diffuse layer of the interface.

Conclusions

The structural organization and vibrational sum-frequency generation spectra of silica/water interfaces were investigated by classical molecular dynamics simulations. The solid surfaces were represented by crystalline surfaces of α-cristobalite and α-quartz, and by the surface of amorphous silica structure with the surface density of hydroxy groups varied from 4.7 to 9.6 OH groups per nm 2 . No deprotonation of the surface silanols was considered and the simulations corresponded to the point of zero charge conditions, i.e. pH ∼ 2 -4 of aqueous phase. The liquid phase was modeled with either a neat water or 1 M NaCl solution. While the silica/neat water interface represents an idealized case study, the silica/NaCl solution interface mimics a real experimental system, where the adjustment of liquid phase pH necessitates the addition of acid or base to the solution.

Calculation of profiles of water structural characteristics along the normal to the surface has allowed to identify bonded interfacial layer (BIL), diffuse layer (DL) and the bulk region of the liquid phase. According to previous works, [START_REF] Pezzotti | Structural definition of the BIL and DL: a new universal methodology to rationalize non-linear χ (2) SFG signals at charged interfaces, including χ (3) contributions[END_REF][START_REF] Wen | Unveiling microscopic structures of charged water interfaces by surface-specific vibrational spectroscopy[END_REF] the bonded interfacial layer is related to a layer of water in the immediate proximity to the surface. The diffuse layer is defined as a region extending to 7 Å beyond BIL and a more distant part of aqueous phase corresponds to the bulk region. The presence of ions in the liquid phase marginally affects the spatial extent of the regions.

Joint probability density maps characterizing the intermolecular H-bonding network and the orientation of O-H bonds in the interfacial regions have shown that the structure of BIL depends on the surface morphology and density of OH groups. The accessibility of the silanols to H 2 O molecules plays a crucial role in the water-surface interactions. The most structured BIL layer was found on the (001) surface of α-quartz, where the high surface density of silanol groups and their specific orientation favor the formation of highly ordered water layer. For all interfaces, both the water-silanol H-bond network and the structure of BIL are largely insensitive to the presence of ions in the aqueous phase.

The structural descriptors computed for DL indicate that structure of the layer is to a large extent surface-independent. For the silica/neat water interfaces, DL water molecules were obtained to preferentially orient their OH bonds towards the surface and this is directly related to an effective negative charge of formally neutral silica surfaces. The presence of ions in the liquid phase strongly affects the intermolecular hydrogen bonding in the layer. The analysis of H-bonds in DL shows a lower number of the bonds in the presence of ions that indicates breaking the H-bonds network upon the formation of solvation shells. Moreover, small size Na + ions approach the surface closer than large size Cl anions. Despite a low surface density of the cations (0.25 Na + /nm 2 at most), it appears that this quantity suffices to alternate the electrostatic potential and to influence the water orientation. Therefore, the difference in the structural organization of DL in the silica/neat water and silica/electrolyte solution interfaces can be ascribed to the chaotropic effect due to the ions and the influence of ion-induced electric field.

The calculation of spectra of the second-order nonlinear susceptibility, Im[χ (2),R ⊥ (ω)], reveal that the spectra of BIL feature low-frequency negative and high-frequency positive bands. The position of maximum of the positive band reflects the strength of water-surface interactions, while the shape of the band is indicative of surface crystallinity. The position and shape of the low-frequency negative band were obtained to be common to all interfaces as the band is due to O-H oscillators directed towards the liquid phase. The BIL spectra are almost insensitive to the presence of solvated ions, with the exception of spectrum of the amorphous silica interface. The surface roughness leading to the absence of net separation between the BIL and DL regions accounts for the peculiarity.

The DL spectra have a shape very similar for all systems and the spectra are dominated by the contribution from the third-order susceptibility. Consequently, the sign of the spectrum is directly related to the sign of electrostatic potential at the frontier between the diffuse and bonded interfacial layers. [START_REF] Pezzotti | What the Diffuse Layer (DL) Reveals in Non-Linear SFG Spectroscopy[END_REF][START_REF] Wen | Unveiling microscopic structures of charged water interfaces by surface-specific vibrational spectroscopy[END_REF] In addition to the third-order contribution, the χ (2 DL spectrum has an intrinsic part due to the fact that BIL acts as a hydrophilic surface for H 2 O molecules of DL. The intrinsic spectrum has the negative sign and its characteristics depend on the mobility of molecules in BIL.

The calculated intensity and Im[χ

(2),R
⊥ ] spectra of silica/neat water interfaces differ from the conventional and phase-resolved VSFG spectra of silica/water interfaces measured at pH values corresponding to the point of zero charge conditions. Contrarily, the spectra of amorphous silica/1 M NaCl solution interface are in a good agreement with the phase-sensitive and homodyne SFG spectra of wa-ter on the surface of fused silica. Since in the experiments, the liquid phase always contains solvated ions, even at low pH values, this results suggests that modeling the nonlinear spectra of silica/water interface has to take into account a real composition of liquid phase by explicit including solutes in the simulations. The ability of classical MD technique to model interfacial systems of thousands of atoms for long simulation trajectories at modest computational expenses, opens a way of obtaining statisically reliable both linear (IR, Raman) and nonlinear (VSFG) theoretical spectra of the systems, and of establishing the relationship between the spectra and structure of the interfaces. The comparison of the results in Fig. S3 -Fig. S6 shows an excellent agreement of both the structural characteristics and nonlinear spectra computed in the two setups. Small differences between the results visible in Fig. S3 and Fig. S6 can be ascribed to a limited statistical averaging. S6 Intrinsic χ (2) spectrum of diffuse layer

In order to find out whether DL of a solid/water interface has a nonzero χ

(2)

⊥ susceptibility in the absence of surface electric field, the following numerical experiment was conducted. 1152 SPCFw water molecules were enclosed between two surfaces of graphite slab of ten graphene layers. Dimensions of the simulation cell were 29.532 Å, 29.8382 Å and 73.4731 Å along the x-, yand z-axis, respectively. The interactions between the molecules and the carbon atoms were described with a Lennard-Jones (12-6) potential with parameters taken from ref. 5. In a first series of simulations, hereafter labelled as "Series 1", all interaction parameters had their original values. In the second series referred to as "Series 2", the D 0 parameter (depth of potential) of the C-Ow atoms pair was multiplied by four for H 2 O molecules of the first interfacial layer (BIL), while all other parameters were left unchanged. The stronger interaction between the BIL molecules and the surface limited their mobility and the "sticky" H 2 O layer mimicked BIL of a silica/water interface.

Simulations of "Series 1" were carried out as described in Section 3.3 of main article. The computational protocol of "Series 2" simulations was slightly different. First, the final configuration of the "Series 1" runs was used to identify BIL water molecules and the C-Ow interaction parameter were changed for these molecules as discussed above. Then, 50 MD simulations of 160 ps length starting from the configuration were performed to produce results for "Series 2". The initial velocities in each of these runs were set independently and the first 120 ps (240000 time-steps) were used as an equlibration stage. The scheme allowed a better sampling the phase space while preserving "sticky" BIL at time t = 0 in each run. The simulations were performed for T = 293 K. ⊥ ] spectra computed for different regions of the interface in the two series of simulations. The DL spectrum is shown in Fig. S11b and it is characterized by a negative band in a wide wavenumber range from 3100 to 3700 cm -1 . One can see that the spectrum has a notable intensity already for the "normal" BIL (Series 1). When BIL molecules become less mobile (Series 2), a low-frequency part of DL spectrum grows in intensity and its minimum shifts to lower wavenumbers. This leads to a downward shift of the negative band in the spectrum of the total interfacial region, Fig. S11c. A gray shaded contour in Fig. S11b shows Im[χ (2),R ⊥ ] spectrum calculated for the bulk region. As water structure in this region is isotropic, the intensity of the calculated spectrum can be used as an estimation of intensity error due to finite size effects, insufficient sampling, etc. Figure S11b clearly shows that the intensity of DL spectrum in both the series of simulations is largely superior this error estimate.

Figure 1 :

 1 Figure 1: A side view of hydroxylated crystalline silica surfaces in the xz plane: a) CRB surface; b) Q101 surface; c) Q001 surface. Colors legend: red -oxygens, peach -silicons, white -hydrogens. The panel a) also shows the atom notation common to all surfaces.

Figure 2 :

 2 Figure 2: System layout used in modeling silica/water interfaces.

Figure 3 :

 3 Figure 3: Definition of coordinates for computing the joint probability density distributions Q(r, u, Z ) (3) and P(u 1 , u 2 , Z ) (4) .

Figure 4 :

 4 Figure 4: a) -z-profile of relative density ρ * , b) -z-profile of S 2 orientational order parameter for dipole andc) -z-profile of S 2 parameter for HH vector of water molecules on the silica/water interfaces; silica surface is specified by its acronym in the upper left angle of graph. The vertical dashed line at z = 0 denotes the position of bottommost Oh atoms, except the AS surface where z = 0 stands for a mean z-coordinate of the silanol oxygens; the negative z values correspond to the liquid phase. The vertical dotted lines delimit Z regions labelled BIL, DL and BLK in the panels (a), see text for details.

Figure 5 :

 5 Figure 5: z-profiles of distribution ions in the liquid phase of the Q101/ES (a) and AS/ES (b) systems. The ρ * (z) density profiles of water molecules are shown by blue lines for comparison.

Figure 6 :Figure 7 :

 67 Figure 6: Two-dimensional maps of joint probability distribution Q(r, u, Z ) for water molecules in BIL of silica/water interfaces: a) -silica/NW, b) -silica/ES; the silica surfaces are specified by its acronym in the upper left angle of map. Each distribution was normalized to its maximum value.

Figure 8 :

 8 Figure 8: Density maps of water in the xy plane in BIL of crystalline silica/neat water interfaces: a) -CRB/NW, b) -Q101/NW, c) -Q001/NW; quantity used in the maps construction is ln(1 + ρ * ). Dashed circles labelled A and B in the panels (a) and (b) indicate two adsorption sites discussed in the text.

Figure 9 :

 9 Figure 9: Configurations of water molecules adsorbed on the CRB (a) and Q001 (b) surfaces. The quadrilateral represents the surface; blue and light pink (panel b) spheres are the oxygen atoms of BIL and DL (panel b) water molecules, respectively; black dashed lines show H-bonds. See Fig. 1 for other color legends.

Figure 10 :Figure 11 :

 1011 Figure 10: Two-dimensional maps of joint probability distribution Q(r, u, Z ) for water molecules in DL of silica/water interfaces: a) -silica/neat water, b) -silica/1 M NaCl solution. Each distribution was normalized by its maximum value.

Figure 12 :

 12 Figure 12: Spectra of Im[χ (2),R ⊥ (ω)] of the silica/water interfaces, upper row -BIL region, middle row -DL region, bottom row -total interfacial region (BIL + DL); each graph is labelled by the acronyms of the silica surface and Z region in the upper left angle.

Figure 13 :

 13 Figure 13: Φ(z) potentials(11) in the liquid phase of the silica/water interfaces. The origin of the z axis is set to z BIL , two vertical dashed lines represent limits between the BIL/DL (z BIL ) and DL/BLK (z DL ) regions.

Figure 14 :

 14 Figure 14: Spectra of -χ

( 3 )

 3 sspz susceptibility of bulk water, adapted from Fig.2bin ref.41. Note the difference in magnitude of the calculated spectra in the panels a) and b).

Figure 15 :

 15 Figure 15: Comparison of the Im[χ (2),R ⊥ ] (a) and |χ (2),R ⊥ (ω)| 2 (b) spectra computed for the AS/neat water and AS/1 M NaCl solution interfaces with the experimental spectra of fused silica/water interface. The experimental spectra were adapted from ref. 16.

Figure 16 :

 16 Figure 16: Im[χ (2),R

Figure 17 :

 17 Figure 17: Schematic rationalization of the Im[χ (2) ssp ] nonlinear spectrum of silica/water interface at low pH values, indicating different contributions and factors influencing the spectral features.

Figure S1 :Figure S2 :

 S1S2 Figure S1: Distributions of mean polarizabilities (a) and dipoles of molecules in bulk liquid water for a system of 500 SPCFw molecules in DIM and r-DIM. Vertical red dashed lines indicate the quantity values of an isolated water molecule.

Figure S3 :

 S3 Figure S3: z-profiles of the relative density (a) and of the orientational order parameter S 2 for dipole (b) and HH vector (c) computed in two setups of simulation cell for the CRB/NW interface, see main article for details.

Figure S4 :Figure S5 :Figure S6 :

 S4S5S6 Figure S4: 2D maps of probability distribution Q(r, u, Z ) (3) for the CRB/NW interface: a) -slit pore layout, b)slab layout of simulation cell; Z regions are indicated with acronyms in the upper right corners.

Figure

  FigureS7and FigureS8present the two-dimensional maps of joint probability distributions Q(r, u, Z ) and P(u 1 , u 2 , Z ), respectively, computed for a second water layer (L2) on the crystalline silica/neat water interface. The L2 layer is easily recognizable in the z-profiles of ρ * density in Fig.4of the main article. Comparison of Figure10aand Figure11a(main article) with FigureS7and FigureS8, respectively, clearly shows that the anisotropic patten in the maps for DL mainly stems from the layer next to BIL of the interfaces.

Figure S7 :

 S7 Figure S7: Two-dimensional maps of joint probability distribution Q(r, u, Z ) (3) computed for the second water layer on crystalline silica/neat water interfaces: a) -CRB, b) -Q101, c) -Q001.

Figure S8 :

 S8 Figure S8: Two-dimensional maps of joint probability distribution P(u 1 , u 2 , Z ) (4) computed for the second water layer of crystalline silica/neat water interfaces: a) -CRB, b) -Q101, c) -Q001.

Figure

  Figure S9 compares the Im[χ (2),R ⊥ ] spectra of BIL and DL with the L2 spectrum for the interfaces.From comparison of Fig.S9band Fig.S9c, one sees that the L2 spectrum, to a large extent, accounts for the spectrum of the diffuse layer. Yet, a part of DL more distant from the surface contributes to the spectrum of the diffuse layer.

Figure S9 :

 S9 Figure S9: Im[χ (2),R ⊥ ] computed for the BIL (a), L2 (b) and DL (c) regions of crystalline silica/neat water interfaces.

Figure

  FigureS10displays maps of water density in the xy-plane in BIL of crystalline silica/1 M NaCl solution interfaces; the maps can be compared with those presented in Fig.8of main article. The comparison shows that the presence of ions in solution does not affect the adsorption sites seen in Fig.8. It should, however, be noted that compared to Fig.8b, the density map in Fig.S10breveals adsorption sites in proximity of inner OH groups (see a red circle in Fig.8b), similar to the second type of adsorption sites found in ref.4 for water on the Q101 surface. It seems that in the case of the surface, the ions reduce mobility of some water molecules and stabilize their H-bonding to the IN silanols.

Figure S10 :

 S10 Figure S10: Density maps of water in the xy plane in BIL of crystalline silica/1 M NaCl solution interfaces: a) -CRB/ES, b) -Q101/ES, c) -Q001/ES; quantity used in the maps construction is ln(1 + ρ * ). Dashed circles labelled A and B in the panel (a) indicate two adsorption sites discussed in the text of main article. Dashed circle A in the panel b) shows an adsorption site akin to the site in the pannel a), whereas a red dashed circle marks a second type of adsorption sites, see text.

Figure S11 :

 S11 Figure S11: Im[χ (2),R ⊥ ] spectra computed for the different interfacial regions of the graphite/water system in the two series of MD simulations as explained in the text, a) -BIL, b) -DL, c) -total interfacial region. Shaded contour in the panel b) shows the spectrum of the bulk region.

Figure

  Figure S11 presents Im[χ (2),R

Table 1 :

 1 Characteristics of silica/water systems.

		CRB surface	Q101 surface	Q001 surface	AS surface
	L x × L y a (Å)	25.026 × 24.885 27.516 × 29.480 25.531 × 24.567 40.315 × 41.432
	OH (nm -2 )	4.7	5.6	9.6	4.7
	No. of atoms in SiO 2 layer b No. H 2 O molecules c	1170 (60) 500	1656 (96) 625	1530 (120) 500	2928 (168) 1754
		482	603	482	1692
	No. of Na + /Cl -ions	9/9	11/11	9/9	31/31
	∆z d (Å)	24.91	24.06	25.90	32.68

a Simulation cell dimensions in the surface plane. b Number of OH groups is given in parentheses.

c

The first and second lines correspond to number of water molecules in the silica/NW and silica/ES systems, respectively. d Interlayer spacing. Defined as difference of z-coordinates of the topmost and bottommost Oh atoms of opposite silica surfaces for the crystalline structures, and as difference of mean z-coordinates of Oh atoms of opposite surfaces for the AS model.

Table 2 :

 2 Parameters of interatomic potentials.

	Nonbonded potentials				Eqn. (1)
		Ow a	Hw	Si	Ob	Oh	Hh	Na +	Cl -
	q (|e -|)	-0.82 0.41 1.10 -0.55 -0.675 0.40	1.0	-1.0
	R 0 (Å)	3.5532	-	4.15	3.47	3.47	1.085 2.424	5.422
	D 0 (kcal/mol)	0.1554	-	0.093 0.122 0.122 0.015 0.3526 0.0128
	O-H potential		E = E 0 exp(-2γ(r -r 0 )) -2 exp(-γ(r -r 0 ))
			Ow-Hw			Oh-Hh
	E 0 (kcal/mol)		115.658			115.658
	γ (Å -1 )		2.1866			2.2431
	r 0 (Å)		1.012			1.0
	H(Si)-O-H potential				E = 1 2 K α (α -α 0 ) 2
			Hw-Ow-Hw			Si-Oh-Hh
	K α (mdyn•Å)		0.727			0.4744
	α 0 (deg.)		113.24			126.00

Table 3 :

 3 Number of hydrogen bonds per water molecule in different regions a of silica/water interface.

	DA pair b	CRB surface Q101 surface Q001 surface AS surface
			BIL
	Ow-Ow c	2.13	2.35	2.40	2.74
		2.02	2.15	2.20	2.51
	Ow-Oh	0.32/0.28 d	0.34/0.26	0.04/0.60	0.30
		0.31/0.27	0.32/0.25	0.04/0.54	0.29
	Oh-Ow	0.30/0.28 e	0.35/0.06	0.41/0.03	0.22
		0.29/0.27	0.36/0.06	0.36/0.03	0.23
	Ow-Ob	0.07	0.02	0.00	0.07
		0.06	0.02	0.00	0.08
	Total in BIL	3.38	3.38	3.48	3.33
		3.22	3.16	3.16	3.11
			DL
	Ow-Ow	3.47	3.47	3.48	3.48
		3.22	3.25	3.22	3.23
			BLK
	Ow-Ow	3.47	3.47	3.48	3.48
		3.23	3.27	3.27	3.27

a See text and Fig.

4

for the region definition. b H-bond donor -acceptor pair. c For each DA pair, the first and second lines correspond to SiO 2 /NW and SiO 2 /ES interfaces, respectively. d Entries correspond to the Ow-Oh(OUT)/Ow-Oh(IN) DA pairs, see Fig. 1. e Entries correspond to the Oh(OUT)-Ow/Oh(IN)-Ow DA pairs, see Fig. 1.

Table 4 :

 4 Mobility characteristics of water in BIL of the crystalline silica interfaces.

		CRB surface Q101 surface Q001 surface Bulk water
	D xy , Å 2 /ps a	0.028	0.051	0.026	0.172
		0.026	0.035	0.030
	τ µ , ps (2)	11.39	9.07	22.78	2.44
		13.83	20.98	34.98

a For each surface, the first and second lines correspond to SiO 2 /NW and SiO 2 /ES interfaces, respectively.

  at most) that can hardly be

	Amplitude (arb. unit)	3000 0.025 0.000 0.025 0.050 0.075 0.100 0.125 a)	3200	3400	3600 CRB/neat water 3800 Wavenumber (cm 1 ) 3000 0.000 0.005 0.010 Q101/neat water Q001/neat water AS/neat water 3 ( ) liquid water 0.015 0.020 b)	3200	3400	3600 CRB/1M NaCl 3800 Q101/1M NaCl Q001/1M NaCl AS/1M NaCl 3 ( ) liquid water

Acknowledgments

The author is grateful to Prof. Daniel Bougeard for critical reading the manuscript.

Conflicts of interest

There are no conflicts to declare.

Electronic Supplementary Information

Structure and sum-frequency generation spectra of water on neutral hydroxylated silica surfaces Konstantin S. Smirnov

Univ. Lille, CNRS, UMR 8516 -LASIRe -Laboratoire Avancé de Spectroscopie pour les Interactions la Réactivité et l'Environnement, F-59000 Lille, France

S1 Computation of system dipole and polarizability

The calculation of spectrum of the resonant part of nonlinear susceptibility χ (2),R by eq. ( 6) requires the knowledge of the dipole and polarizability of system. For a slab of water molecules these are obtained by summing up contributions of individual molecules. In the present work, the molecular characteristics were computed in a restricted dipole interaction model (r-DIM). Dipole interaction model (DIM) considers each molecule as a polarizable point dipole. 1 For a system of N such interacting dipoles, a vector of dipoles m (length 3N) can be written as

where m 0 is vector non-interacting dipoles, a stands for matrix of dipole's polarizabilities, the T matrix is built of dipole-dipole interaction tensors, and E 0 is a vector of external electric field. The solution of (S1) is given by 1,2 m = m 0 + a eff (E 0 + Tm 0 , (S2)

where

From the computational viewpoint, the inversion of 3N × 3N matrix in (S3) is clearly the limiting step of model efficiency. To overcome the bottleneck, r-DIM takes the matrix of effective polarizabilities as a eff = a and the computation of molecular dipoles is then performed as

where the quantities in m 0 and a are obtained with models described in ref. 3. Figure S1 displays distributions of mean polarizability and dipole values of molecules in bulk liquid water computed with DIM and r-DIM. The data were obtained for a system of 500 SPCFw molecules in an MD simulation for a temperature of 293 K. One sees in Fig. S1a that the local field correction has a relatively small influence on the polarizabilities. Figure S1b shows that despite an underestimation of dipole values by ca. 10 %, the r-DIM model correctly captures the increase of dipoles due to the external field,.

Figure S2 compares χ

(2),R

⊥ spectra for the water/vapor interface obtained with the two models. One sees that taking a eff = a in (S2) has a little impact on both the shape and intensity of the spectra that justifies the simplification.

S5 Mobility characteristics of interfacial water

The mobility of water molecules was characterized by the parallel-to-surface diffusion coefficient (D xy ) and the rotational relaxation time of molecular dipole (τ (2) µ ). The quantities were obtained from the analysis of the mean-squared displacement of molecules in the xy-plane (MSD xy (t)) and of the dipole rotational time correlation function C µ 2 (t), respectively,

where P 2 () stands for the 2nd-order Legendre polynomial and u is the unit vector along the molecular dipole µ.

where Γ() stands for the gamma function.