C Humbert

O Pluchery

E Lacaze

A Tadjeddine

B Busson

Multiscale description of molecular adsorption on gold nanoparticles by nonlinear optical spectroscopy

Nonlinear optical Sum and Difference-Frequency spectroscopies are used to probe and model the surface of thiophenol-functionalised gold nanoparticles grafted on a Si(100) substrate through two different silanization procedures. By scanning the [980-1100 cm -1 ] infrared spectral range with the CLIO Free Electron Laser, ring deformation vibrations of adsorbed thiophenol are investigated. Quantitative data analysis addresses three levels of organization: microscopic, nanoscopic and molecular. Grafting with paminophenyl-trimethoxysilane shows an increase of around 40% in surface density of nanoparticles (Ns) as compared to 3-aminopropyl-triethoxysilane. The relative amplitudes of the resonant and nonresonant contributions to the SFG and DFG spectra are discussed in terms of Ns, Fresnel reflectivity factors and local amplification of the nonlinear signals by coupling to the surface plasmon of the particles. They are shown to quantitatively scale with Ns, as measured by atomic force microscopy. Vibration mode assignment is performed through a critical analysis of literature data on IR and Raman spectroscopies coupled to DFT calculations, for which a methodology specific to molecules adsorbed on gold atoms is discussed.

Introduction

Nonlinear optical infrared-visible Sum-Frequency Generation (SFG) spectroscopy is attracting a growing interest for the study of gold nanomaterials as demonstrated in various and recent works. 1- 3 Since we recently showed that molecular SFG signal can be strongly enhanced when molecules are grafted on gold nanoparticles (AuNps) in presence of plasmonic exaltation, 3 a study by SFG of adsorbed molecules becomes possible even at very low coverage. SFG, an intrinsic vibrational spectroscopy of interfaces, offers promising perspectives on such nanostructured materials. On one hand, such investigations performed by SFG spectroscopy are most of the time limited to the characterisation of methyl and methylene moieties of organics layers composing the probed interface due to the difficulty to maintain widely tuneable tabletop infrared (IR) sources (Optical Parametric Oscillators (OPO), Optical Parametric Amplifiers) outside the 2.5-4 µm spectral range. Some exceptions relate to the studies of solid/air and liquid/air organic interfaces; 4-7 electrochemical or catalytic analyses of cyanide, carbon monoxide and C60 fullerene monolayers 8 up to 9 µm; amide I modes of proteins at various interfaces 9 and more recently to the case of surface phonons of quartz at 12 µm. 10 No work was performed so far on nanomaterials functionalized by organic components in the mid-IR spectral range, which however contains the fingerprint of organic molecules. On the other hand, it should be mentioned that a SFG study was performed around 10 µm on the coupling of surface states of Ag nanoclusters embedded in a Si3N4 matrix by tuning the visible excitation wavelength to the surface plasmon resonance (SPR) of the particles. 11 This exemplifies the way often used to overcome the IR range limitation by coupling the SFG spectrometer to a Free Electron Laser (FEL). 12,13 The latter offers a wide IR domain ranging from 5 to 150 µm, with higher power than any tuneable IR source based on a tabletop laser at this date. 14 In this paper, the coupling of nonlinear optics to the surface optical properties of nanomaterials is investigated with a unique Sum and Difference-Frequency Generation (SFG/DFG) setup coupled to the CLIO FEL. We present SFG and DFG measurements performed on two kinds of samples made of gold 55 spherical nanoparticles grafted on silicon with either paminophenyl-trimethoxysilane (APHS) or 3-aminopropyltriethoxysilane (APTES) to monitor the surface coverage (Scheme 1). AuNps are functionalized with thiophenol (C6H5SH) molecules and we show that reliable information on both the 60 molecules and the particles can be obtained even for a volume density as low as 4%. To reach that purpose, we perform a detailed three scale description (Figure 1) of the various contributions to the SFG/DFG spectra, which allows to evidence specific chemical and physical properties of the nanostructured samples at all scales: 65 microscopic, where their optical properties are described according to the 3-layer model, with the help of atomic force microscopy (AFM); nanoscopic, where the individual nanoparticles are taken into account to enlighten the differences induced by the grafting molecules and the enhancement effects due to a coupling to the 70 SPR; molecular, where the vibrational response of individual thiophenol molecules is analysed by DFT to understand the nature and activities of the three observed vibration modes, addressing controversy in this rich spectral range.

75 Scheme 1 Molecules of p-aminophenyl-trimethoxysilane (APHS, left) and 3-aminopropyl-triethoxysilane (APTES, right) used to graft the gold nanoparticles on the silicon substrate.

Results and discussion

Vibrational fingerprint in nonlinear optical spectroscopy SFG and DFG experimental data of thiophenol adsorbed on AuNPs are presented in Figure 2 for the two grafting silane molecules, i.e. APHS (2A and 2B) and APTES (2C and 2D), respectively. Three observations can be made from these measurements. Firstly, all SFG/DFG spectra display three 10 vibration modes (labelled 1, 2 and 4 with reference to the discussion below) related to the phenyl ring of thiophenol, located at 10002, 10232 and 10712 cm -1 , interfering with a constant background generated by the AuNPs. Secondly, the interference induces different shapes in SFG and DFG spectra. In the SFG case, the interference pattern produces a peak for the 1000 cm -1 mode, a dip for the 1071 cm -1 mode and a Fano profile for the intermediate 1023 cm -1 mode. On the contrary, for the DFG case, all modes appear as peaks. This is related to the change in sign of the damping constant (Гn) of the Lorentzian oscillators between SFG and DFG, detailed below (see equations (1-2) as well as Ref. 15). The third interesting feature relates to the silane used to graft the AuNPs on the silicon substrate. We observe that SFG/DFG intensities of the AuNps and the 3 vibration modes contributions are stronger in the case of APHS silanisation. This is an indication of a higher coverage of the silicon substrate with APHS, as will be proved further. 1) and ( 2). SFG and DFG scales are the same to ease the comparison between the spectra.

To extract quantitative information from these measurements, (2) as developed in the next part. The experimental SFG/DFG nonlinear susceptibility of an interface is generally modelled by: 3

χ eff,SFG (2) = C SFG e iϕ SFG + ∑ a n (ω IR -ω n +iΓ n ) n χ eff,DFG (2) = C DFG e iϕ DFG + ∑ a n (ω IR -ω n -iΓ n ) n (2)
where the first complex term in each sum is a non-resonant (NR) 55 contribution (i.e. NR with the IR frequency) related to the substrate, more specifically in our case to the AuNPs. The Lorentzian oscillators describe the three vibration modes of the thiophenol molecules, with n and n their frequencies and damping constants, respectively. Their amplitudes are linear 60 combinations of products of the components of the IR transition dipole moment and Raman polarizability tensor. 16,17 Equations ( 2) are used to fit the four experimental spectra of Figure 2. Prior to the fits, spectra were corrected from stray light background, account was taken of the actual visible and IR powers (equation ( 1)) and spectra were normalized to a common reference sample (ZnS crystal). There are several non-equivalent sets of parameters which may produce a correct fit, even if they do not correspond to a global minimum of the fitting algorithm. In order to fully explore the sets of possible fit parameters, we used the procedure already applied to thiophenol in [18] by fitting the data with complex lorentzian amplitudes, getting the eight alternate sets of parameters and analysing them to reach the best fits with real amplitudes. In order to choose the correct fit, we have therefore imposed constraints on the parameters. First, we have looked for coherent signs for the amplitudes of the vibrations. The only possibility showed opposite signs for mode 1 on one side and modes 2 and 4 on the other side. This is coherent with the SFG literature on these modes. 14,18 A plausible explanation is that the IR transition moment of mode 1 has an opposite sign upon adsorption on AuNPs with respect to the two others as shown in Table 1. Fitting the curves in this condition leads to a NR phase (SFG and DFG) of the AuNps contribution analogous to Refs. 3,18 on silver and AuNPs. In order to rationalize the comparison between the curves, we fixed the n to their mean values, resp. 3.5, 3.5 and 5.5 cm -1 . We note that the larger 4 was already observed in Refs. 12,18. Finally, the most stringent constraint lies on the amplitudes. Considering that any external effect on the amplitudes of the modes (Fresnel factor, local field, coupling to the SPR) must affect the three modes in the same way, the ratios of the three amplitudes must remain constant, as displayed in Table 1. We remark that the mandatory theoretical coherence between vibrational SFG and DFG spectra recorded on a same sample imposed severe constraints on our fits, leading to a unique solution for the parameters which in fact differs from the global minimum found by the fitting program. The usefulness to record DFG in parallel to SFG is evidenced in this case as it allows to calibrate the fits of SFG data. 15,18 The following conclusions may be drawn from the fit parameters. Firstly, the values for the NR background of AuNps remain coherent with the FG response of bulk gold. 15 This effect explains well the interference patterns observed in the FG experimental spectra and relates to the excitation of interband transitions in the AuNps, resonant with the FG energies. 3 Actually, the NR amplitude (CSFG, CDFG) is higher for APHS than for APTES, with ratios 1.19 for SFG and 2.41 for DFG, respectively. In addition, the NR phases () are coherent between APHS and APTES. Mean values lie around 80° for DFG and 140° for SFG, the latter being a bit higher than usual for analogous materials. It means that the nonlinear signal from AuNps is not far from phasequadrature with respect to the vibration modes of the thiophenol molecules, and follows the trends observed for molecules adsorbed on bulk gold. 15 This hints at analogous contributions from AuNps, implying in turn analogous local environments between both samples. This conclusion is confirmed by the observation that the overall positions of the three peaks remain remarkably constant within experimental uncertainties, even when molecular adlayer is changed from APHS to APTES. Secondly, it is possible to compare the absolute resonant amplitudes from the values of fit parameters a1. For a given sample, the ratio DFG/SFG gives an estimate of the local amplification due to a resonance phenomenon with the surface plasmon of the particles at DFG and SFG frequencies, respectively. This ratio is 1.72 for APHS and 1.60 for APTES. The resonant and NR amplitudes ratios (APHS/APTES) will be discussed in more details further. 4 (cm -1 ) 5.5 5.5 5.5 5.5

Optical and morphological factors at the nanoscale

As mentioned above, in order to extract the chemical and physical 65 information from our FG spectroscopic measurements, we need to detail the role of the various contributions to the FG intensity defined in equations ( 1) and ( 2). The χ eff (2) third rank complex tensor can be expressed in ppp polarisation combination as:

χ eff,ppp (2)
= -F x (ω FG )F x (ω vis )F z (ω IR )cosθ FG cosθ vis sinθ IR χ xxz (2) -F x (ω FG )F z (ω vis )F x (ω IR )cosθ FG sinθ vis cosθ IR χ xzx (2) +F z (ω FG )F x (ω vis )F x (ω IR )sinθ FG cosθ vis cosθ IR χ zxx (2) +F z (ω FG )F z (ω vis )F z (ω IR )sinθ FG sinθ vis sinθ IR χ zzz (2)
(3) 70 which will be written in short notation including the angular projections:

χ eff,ppp (2)
= -F xxz χ xxz (2) -F xzx χ xzx (2) + F zxx χ zxx (2) + F zzz χ zzz (2) (4)
The Fresnel factors Fi (i=x,y,z) describe the surface reflectivity, and we calculate them in a three-layer model: Air/AuNps/Silicon, 75 with respective refractive indices n1, n* and n2, as depicted in Figure 1(A). In fact, the silane monolayer and the thiophenol submonolayer do not affect significantly the sample reflectivity in the probed spectral range with respect to the Si substrate or AuNps layer, which contributions largely dominate in our samples. For an 80 isotropic surface, the explicit complex contribution of the Fresnel factors can be calculated as 16

F x (ω) = 2n 1 (ω)cosθ 2 n 1 (ω)cosθ 2 +n 2 (ω)cosθ 1 F y (ω) = 2n 1 (ω)cosθ 1 n 1 (ω)cosθ 1 +n 2 (ω)cosθ 2 F z (ω) = 2n 2 (ω)cosθ 1 n 1 (ω)cosθ 2 +n 2 (ω)cosθ 1 ( n 1 (ω) n * (ω) ) 2 (5)
where 1 and 2 refer to the angles of incidence of the -frequency beam in media 1 and 2, knowing that n 1 (ω)sinθ 1 = n 2 (ω)sinθ 2 .

85

Fi should be calculated for each optical frequency contributing to the FG process, i.e. for =IR, vis, SFG and DFG. Index of air is taken as unity and, for the n2 refractive index of silicon, we used a double exponential fit of the literature data. 19 The difficulty lies in the derivation of an expression for the refractive index of the AuNp layer (n*) in the desired spectral ranges. To address this point, we evaluate the dielectric constant *=(n*) 2 by considering a layer of nanoparticles embedded in a host matrix (air in our case, with h=1). We choose here the Maxwell-Garnett formalism 20 for the effective dielectric constant:

ε * = ε h 1+2rv Au 1-rv Au , r = ε Au -ε h ε Au +2ε h (6)
where vAu is the volume ratio of the space occupied by the AuNps with respect to the total volume of the layer. We consider an analytical expression of the gold dielectric constant Au in the visible range 21 and its value at ~1000 cm -1 for the IR range. 19 In order to determine the AuNP density (hence the volume ratio), we performed atomic force microscopy (AFM) measurements for the two systems, i.e. nanoparticles attached to silicon through APHS and APTES as displayed in Figure 3. In a first step, at large scale (10µm), we always obtained images similar to the ones shown on Figure 3a for APHS and 3b for APTES. These AFM images allow to establish the proportions of AuNP monolayer on the substrate. In a second step, we zoomed on the monolayers (see inserts in Figure 3a and3b). The local densities were thus obtained as an average on six different areas. Combining this value to the monolayer proportion finally allowed to establish the average AuNP surface density (Ns), leading to a value of 3.12x10 10 cm -2 for APHS and 2.28x10 10 cm -2 for APTES. For particles with a 172 nm mean diameter, this leads to volume densities within the intermediate layer vAu = 4.72% and 3.45% for APHS and APTES grafting, respectively.

The resulting Fresnel contributions to the SFG and DFG signals (i.e. related to χ eff,ppp

) appear in Table 2 for AuNps grafted on silicon through APTES or APHS, with the dominant components underlined. For the sake of comparison, the same quantities are calculated for the Air/Thiophenol/Au(111) 3-layer model interface within the same spectral range. In this second case, the organic layer cannot be dismissed because the nonlinear process occurs inside this layer. We used n1=1, n 2 = √ε Au from Refs. 19,21 for Au(111) and n*=1.2 as calculated from the Lorentz-Lorentz model of an organic monolayer. 16 It is worth noting that the |Fzzz| contribution largely dominates in the SFG/DFG intensity of the AuNps with respect to the other components. This is mainly due to the high refractive index of silicon at all wavelengths, which favours the Fz component. On the contrary, for Au(111), the main effect stands in the high refractive index of gold in the IR range, which only favours the Fz(ωIR) component. Another interesting feature is that the dominant Fresnel factors of the AuNps are of the same order of magnitude as the Au(111) ones, as underlined in Table 2. This result shows that a direct comparison between 50 intensities of SFG/DFG spectra can be performed for thiophenol grafted on Au(111) and AuNps. This validates the comparison performed in Ref. 3 which evidenced an enhancement as high as 25 due to plasmonic exaltation.

From the comparison between the APTES and APHS samples 55 in Table 2, the differences follow the values of the surface densities. The Fresnel contribution decreases when vAu increases.

It is in fact very sensitive to the surface densities of particles through the (n*)² factor in Equation ( 5). We also observe that, for a visible wavelength set at 523.5 nm, the Fresnel contribution to 60 the SFG signal is slightly stronger than the DFG one as a consequence of the optical properties of gold in the visible range. As a consequence of the above analysis, the tensor decomposition of the ppp contribution defined in equation ( 4) simplifies into:

χ eff,ppp (2)
≈ F zzz χ zzz

(2)

(7)

This proves that only components of the molecular response perpendicular to the surface, i.e. IR transition moments and Raman 70 polarizabilities, contribute to the signal, even for spherical AuNps for which one could have expected equal contributions of all directions. The presence of the silicon surface beneath the AuNps leads to this selection rules, which is an important result of the present calculations.

75

Consistency between microscale and nanoscale

In order to validate the FG data analysis at the microscale and nanoscale used so far, we have to check whether the whole set of data is consistent. We focus on the comparison between APTES 80 and APHS samples, which in principle only differ by the surface density of adsorbates. The ratios of the fit parameters for APHS and APTES grafting is displayed in Table 3, together with the χ zzz (2) ratios deduced for the resonant and NR contributions after correction from Fresnel effects (equation ( 7)). has not yet been fully addressed by experimental methods, but it seems reasonable to consider that the visible laser excites the SPR, inducing an enhancement equal for both SFG and DFG signals. In addition, the SFG and DFG beams may also couple to the SPR, generating an additional factor to the amplification. This FG factors depends on the SPR properties at the FG frequencies, accounting for the differences in amplitude between SFG and DFG responses. As already evidenced in a previous study, 3 the SPR of isolated particles (around 520 nm) should favour the SFG response. However, gold nanoparticles tend to aggregate upon adsorption on the substrate (Figure 3), leading to the appearance of a second red-shifted plasmon resonance, which eventually boosts the DFG response. This has already been shown for similar deposition procedure of AuNPs on silicon substrates. 3 This explains why, in the present study, the local amplification favours DFG rather than SFG. We estimate the DFG/SFG amplification ratio at 1.89 for APHS and 1.73 for APTES when corrected from Fresnel effects, i.e. calculated on χ zzz

(2) ratios.

As a summary of the analysis developed so far for the FG 40 intensity coupled to AFM measurements and to calculations of the Fresnel contributions, we may draw the four following conclusions: 1) the four sets of SFG and DFG data are consistent, involving the contributions of three identical vibration modes in the direction normal to the silicon surface; 2) the amplitude ratios of the three modes remain constant whatever the nonlinear vibrational spectroscopy and the silane, hinting at analogous molecular geometries on the particles and local environments; 3) the differences between absolute values of the amplitudes (resonant and non resonant), after correction of Fresnel reflectivity After the investigation of the FG intensity at the microscale and nanoscale, the last step in such an analysis, when done on a flat surface, would be to relate the microscopic χ (2) to the molecular response, provided that at least two polarization contributions could be measured (usually ppp and ssp). This is routinely done 65 under the conditions that the surface is isotropic and all molecules share the same tilt (even if dual angle distributions appear necessary for complex cases 22 ). An additional approximation on the structure of the molecule (i.e. cylindrical symmetry around the direction of the IR transition moment of a given vibration) is 70 commonly added. 16 As a consequence, a multiscale description relates molecular IR and Raman activities of the vibration modes to the microscopic SFG response, giving access to unknown parameters, commonly the molecular tilt angle. A complete simulation of an SFG spectrum, including several vibration modes, 75 must take into account the relative IR and Raman activities and their full decomposition in the molecular frame. It requires usually a molecular DFT calculation. 6,7 The situation in the case of AuNPs is more complicated than the flat surface archetype, as further effects blur the molecular 80 response building up the microscopic FG signal and its analysis. Orientation averaging represents a challenge indeed when the molecules decorate gold nanoparticules. Molecular orientation on one particle is not constant but rather distributed along a sphere. In addition, the density of the molecular packing at the surface of the 85 particle is unknown, as is therefore the degree of order at this scale. As for the particles themselves, their distributions of shapes and sizes complicate the analysis. In addition, the distribution functions of the molecular tilt and twist angles are unknown. It would seem rather hazardous at this stage to draw hypotheses on these 90 parameters. Finally, for SFG and DFG spectroscopies, Fresnel factors involving a polarization of light parallel to the surface (i.e. ssp contribution) are very unfavourable owing to the high refractive index of silicon, as calculated in Table 2. As a consequence, we consider that it is not possible, without additional 95 data on the structure and organization of the molecules at the interface, to use FG data for a detailed orientation analysis. Nevertheless, we demonstrate that DFT analysis of the geometry and vibration modes of the thiophenol molecule, both free and bonded to gold, leads to important information when linked to the 100 SFG and DFG data.

Several publications have focussed on the analysis of thiophenol vibration modes by IR, 23,24 Raman 24-27 and SFG 5,12 spectroscopies. A comparison evidences that there are still some discrepancies in their interpretations. In the early literature, thiophenol was 105 described in the convenient C2v geometry as an extension of the analysis of monosubstituted benzenes to make the analytical calculations easier. 24 This means that the H atom was considered negligible and merged into a SH pseudo-atom. For the free thiophenol molecule, experimental data show that there are four 110 active modes in the 1000-1100cm -1 range, which qualitative IR and Raman activities are shown in Table 4. The positions and activities of the modes are taken from Refs. 24 and 25, even if these values change slightly within the literature and with the experimental tool. In the same way, thiophenol linked to a metal atom, either silveror gold-benzenethiolate, has been investigated by IR, Raman, SERS and SFG. A summary of the main bands of interest is also shown in Table 4, from Ref. 25. Only three modes survive upon adsorption on metal, but their origin and symmetries are debated. Within the C2v approximation, in-plane modes may be symmetric/a1 or antisymmetric/b2 with respect to the plane

DFT analysis of free and adsorbed thiophenol molecule

To clarify the vibrational properties of thiophenol, we have calculated by DFT its molecular geometry and vibration modes in the 1000-1100cm -1 region. We first made use of the most usual functional (B3LYP) and basis sets (6-311++G(d,p)) for such molecular systems to account for our data. In order to validate the calculations, we first tried to reproduce the IR and Raman data before analysing the SFG intensities. As will be illustrated below, we faced two problems to account for the literature results. Firstly, Table 1 and the discussion above show that adsorption on gold plays a great role on the energies of the vibrations in this region, which implies that it has to be taken into account in the calculations. In order to simulate adsorption onto the gold 40 nanocrystal on the geometric and vibrational parameters, calculations have been performed both on free thiophenol and on a model thiophenol molecule linked to one gold atom (goldbenzenethiolate, GBT) as illustrated on Scheme 2. Such a procedure was used for the simulation of electronic properties of 45 CO on platinum 28 and CNˉ on transition metals, 29 although more accurate predictions require the use of metal clusters or periodic slabs and heavier calculation methods. 30 It appeared that B3LYP functional produces inconsistent mode assignment with GBT (see below and Supplementary Information). In order to produce 50 consistent results, functional PBE0 was successfully used on both molecules. However, such basic DFT results were not satisfactory enough to reproduce experimental data, even with a correction factor. Therefore, a second-order perturbation anharmonicity step was added to frequency calculations, 31 which resulted in a very 55 nice accuracy on the predicted frequencies (we note that this method does not provide anharmonic corrections for the IR and Raman activities. Details of the DFT calculations of the present paper are given in the Supplementary Information). As illustrated on Scheme 2, the ground states geometries differ 65 appreciably. Thiophenol has a planar configuration whereas GBT converges to a perpendicular geometry, where the S-Au bond rotates by 90° to reach the symmetry plane perpendicular to the phenyl ring. As for the latter, it is highly probable that the actual rotation angle around the C-S bond in a rather densely packed zone 70 lies between the two extreme values as a consequence of molecular packing and low free rotation barriers. For this system, additional steric and molecular packing effects at the surface of the gold nanoparticle are not taken into account here and further restrict the free rotation. Both geometries possess Cs symmetry.

75

All calculations confirm the presence of four vibration modes between 1000 and 1150cm -1 (Table 5). The harmonic approximation for the vibration energies is, as expected, too high (see Supplementary Information). We could check that a uniform scaling factor, even for as few as these four modes, does not 80 satisfactorily account for the experimental frequencies. On the contrary, a mere anharmonic calculation limited to the four modes under study gives very satisfactory results without any scaling factor for B3LYP, and with a uniform additional scaling factor 0.992 for PBE0 (Table 5). For GBT, it appears that B3LYP is 85 unable to handle the metal atom in a correct way within our calculation scheme, leading to inconsistent results for the assignment and position of the modes. PBE0 proves more robust in that case, giving very satisfactory results after anharmonic correction and with the same 0.992 scaling factor as for thiophenol.

90

From Tables 4 and5, we may draw the following conclusions: i) modes 1, 2 and 4 are SFG-active for GBT; ii) mode 3 is SFGinactive as a consequence of its vanishing Raman activity in GBT; iii) upon adsorption on gold, mode 4 experiences a strong red shift in position (experimental value 18 cm -1 , calculated value 24 cm -1 ), 95 which could lead to its confusion with mode 3. The explanation of Ref. 27 is correct, as this mode encompasses a strong contribution of the C-S stretch, which is directly perturbed by the replacement of the terminal H by the Au atom; iv) DFT predicts very well the positions of the modes and their evolution upon adsorption as long as anharmonicity is considered; v) the coarse prediction of qualitative IR and Raman activities is correct when compared to experimental values; vi) DFT lacks accuracy on the fine prediction of actual values for IR and Raman activities, especially for GBT.

Table 5. Calculated (anharmonic) vibration modes of thiophenol and GBT.

Frequencies for B3LYP are not scaled. A 0.992 scaling factor is applied for PBE0, compared to experimental ones in Refs. 24 and 25. Considering Table 5, the effects of adsorption on the IR and Raman activities inferred from DFT can be summarized as follows.

Upon adsorption, the IR activity of mode 1 is enhanced, whereas it decreases for modes 2 and 4 and remains constant for mode 3. On the contrary, Raman activity of modes 2 and 4 rises upon adsorption, whereas it vanishes for mode 3. However, a direct comparison remains hindered by the geometry change between thiophenol and GBT. In addition, as seen before, we have to consider at least a partial free rotation around the C-S bond. In order to separate the effects of the bonding to gold and of the rotation around the C-S bond, we also have performed similar calculations on the two complementary conformers of higher ground state energy: thiophenol in a perpendicular configuration and GBT in a planar geometry. The results show that around 2/3 of the energy shift of mode 4 are due to the gold bond and 1/3 to the conformation. Taking into account the effect of the atom exchange alone leads to a ~16 cm -1 shift, very close to the experimental value. The ~5 cm -1 slight shift of mode 3 is mainly due to conformation. As for the IR and Raman intensities, the effects are less clear. It appears that the Raman enhancement of modes 2 and 4 and the damping of mode 3 are mainly due to the gold atom, even if conformation also favours a decrease of the Raman activity for mode 3. The geometric factor is the cause of the decrease in the IR activity of mode 4, whereas both effects participate to the increase for mode 1 and mode 2. In all cases, it appears obvious that the adsorption process must be explicitly taken into account when modelling vibrational properties of molecules adsorbed on metal.

DFT proves very useful for the assignment of the modes and the comparison of their linear optical activities (i.e. IR and Raman activities). The three vibration modes appearing on Figure 2 are therefore identified with modes 1, 2 and 4 (Figure 4 below), respectively, their position being in full agreement with the calculated ones (Tables 1 and5). Considering the plane of symmetry of the GBT molecule, they are all totally symmetric (inactive mode 3 is antisymmetric). Their IR transition dipole moments have therefore two components, one dominant along the C1-C4 axis and the other perpendicular to the 60 ring plane. The Raman tensor has four independent nonvanishing terms, building up a Raman polarizability ellipsoid with three nonvanishing principal axes. Focussing on GBT, the transition dipole of mode 1 presents a transition dipole moment parallel to the phenyl plane, oriented in the opposite direction with respect to 65 those of modes 2 and 4. As for Raman polarisabilities, mode 4 is opposite to modes 1 and 2. As a consequence, DFT calculations on GBT cannot reproduce correctly the signs of the SFG activities of the three modes (Table 1). The same analysis performed on the free thiophenol molecule accounts for the correct sign alternation (Figure 4). However, all modes also have a minor contribution perpendicular to the ring, and the Raman polarizability contributions include a total of four non vanishing terms, which makes uncertain the conclusions on the final signs of the resonant hyperpolarizabilities. In addition, as stated above, the various orientations of the molecules at the surface of the particles require caution when correlating the individual molecular properties to the averaged FG response, in particular when one considers that the 10 actual geometry of the molecule on the surface is not necessarily exactly the perpendicular one.

As far as SFG intensities are concerned, we note that their qualitative behaviour cannot be correctly reproduced by DFT. In fact, considering an estimate of SFG activity from the product of IR and Raman activities as a first order approximation, it appears that GBT molecule overestimates the contribution of mode 4, whereas in the free molecule mode 1 is largely underestimated (which is coherent with the dramatic changes induced by adsorption). In fact, nor the free molecule neither the gold-bonded 20 one leads to a molecular SFG activity of mode 1 bigger than that of mode 2 in the molecular frame. Of course, orientation averaging will lead to partial cancellation of some hyperpolarizability contributions to the microscopic FG activity as a consequence of symmetry, so that the latter cannot be inferred directly from the product of IR and Raman activities. Nevertheless, considering that 1) DFT already encounters troubles to account for the IR and Raman data from literature; 2) a calculation of the effects of orientation averaging seems out of reach as a consequence of the lack of knowledge on the actual orientation distributions of the molecules on the nanoparticles; 3) discrepancies between experimental data and DFT simulations (e.g. mode 1 vs mode 4) are too strong to be due only to orientation effects; 4) a careful analysis of individual components of the IR dipole moments and Raman polarizabilities in the molecular frame shows that they are 35 roughly identically weighted for the three modes and that, even by linear combination of components, it is not possible to recover the activity ratios (for example, all Raman polarizability components of mode 2 are basically twice those of mode 1); we consider that, even with a complete description of the structure at the nanoscale, 40 the present DFT results may not quantitatively describe the FG spectroscopy results. The reasons lie in the amount of approximations done during the calculations: choice of the exchange-correlation functionals (B3LYP and PBE0 might not be the most appropriate for the Au atom), modelling of the Au surface by only one atom, no frequency dispersion. The fact that we did not take into account the effects of local electric fields (either from the laser beams or from the SPR) may also partially account for this result. Surprisingly, we note that good estimates of the SFG intensity ratios can be obtained by combining the IR activity of the GBT molecule to the Raman activity of the free molecule. In particular, it becomes possible to recover the leading SFG intensity of mode 1. However, at this stage, no theoretical background justifies such a methodology.

Literature confirms that an accurate calculation of the 1000-55 1100cm -1 region of phenyl ring vibrations by DFT is in general not a straightforward task. This zone seems rather impervious to standard DFT analysis. Adsorption on metal further induces strong modifications of the vibrational parameters (energy position, IR and Raman activities), at least for vibration modes impacted by the 60 S-metal link. Such discrepancies between DFT calculations and experimental SFG data have already been pointed out in this zone in the past. 7 The authors even extended at the time their comment to the whole fingerprint mid-infrared region [1000-2000cm -1 ]. In a general way, the few published attempts to reproduce SFG results 65 by DFT give contrasted results: an overall agreement may sometimes be found 6 although accuracy remains out of reach. Again, this is probably due to the inadequacy of the exchangecorrelation parts of the functionals used in the literature for such analyses (mainly B3LYP). However, contrary to the mid-infrared 70 situation, it seems possible to account rather nicely for the high energy range [2000-3000cm -1 ]. 7,32 To summarize, we have shown that the free molecule alone could not account for all the experimental properties of the system. Only a comparison between the calculations performed on 75 thiophenol and GBT could elucidate most of the properties. More precisely, mode assignment and energies were successfully reached on GBT through the addition of anharmonic corrections. The free molecule leads to better relative signs of the SFG activities, whereas none accounted satisfactorily for the SFG 80 relative amplitudes. Considering the rather high amount of work necessary to account for the essential vibrational properties of adsorbed thiophenol, coupled to the difficulty encountered to quantitatively predict SFG intensities from the single molecule, caution seems recommended when using DFT results as an input 85 for SFG predictions, especially when conclusions are drawn on molecular orientation of adsorbates on a surface.

Experimental

All chemicals were purchased from Aldrich, except other mention. AuNPs are synthesized from 1.7 mg of gold salts (HAuCl4.3H2O 90 99.999 %) dissolved in 20 ml Millipore water following Turkevich method. The solution was then heated up until boiling and stirred vigorously. At that point, 0.8mL of 8.5x10 -4 M trisodium citrate (Na3C6H5O7) was quickly added at one time while heating and stirring for 30 minutes. The resulting solution turned successively 95 from light yellow to deep gray and finally dark red after a few minutes. The corresponding UV-visible absorbance spectrum displayed the expected absorption peak at 524 nm usual for AuNPs in water. 3 Wafers of ultrasonically cleaned n-doped silicon (1x1 cm 2 , Siltronix) were silanized in an absolute methanol solution. 100 After being grafted on such a wafer in methanol, AuNps were functionalized with a 10 -3 M thiophenol (C6H5SH) solution dissolved in dichloromethane (CH2Cl2). The use of two different silanes, differing by their chain length and nature (p-aminophenyltrimethoxysilane (APHS), H2N(C6H4)Si(OCH3)3 or 3-105 aminopropyl-triethoxysilane (APTES), H2N(CH2)3Si(OC2H5)3), aimed at monitoring and comparing the different surface coverages by the AuNps.

SFG and DFG measurements were performed in a ppppolarisation scheme for the SFG/DFG, visible and IR beams, 110 respectively. The sketch of the experiments is depicted in Figure 1. It is clear that the only change between SFG and DFG experimental setups resided in the detection angle and the measured frequency of each nonlinear signal, both related to the SFG/DFG rules of momentum and energy conservations (kSFG/DFG = kVis  kIR and ωSFG/DFG = ωVis  ωIR), as detailed in Ref. 3. The SFG/DFG setup coupled to the CLIO Free Electron Laser (FEL) European Facility is described in details elsewhere. 14 Briefly, its dedicated pulsed and intense IR laser beam (~1ps, ~10µJ/pulse reduced to 2µJ/pulse for SFG/DFG experiments, ~5 cm -1 resolution) was tuned around 10 µm in order to probe specific vibration modes of ring deformation of thiophenol. The visible laser beam at 523.5 nm (~5 ps, 4µJ/pulse, ~6 cm -1 resolution) originated from the frequency doubling of a pulsed and amplified Nd:YLF laser source at 1047 nm through a nonlinear BBO crystal.

The laser sources had the same repetition rate for the macropulses (25 Hz) and micropulses (62.5 MHz) fixed by the CLIO FEL clock, in order to synchronize temporally the IR and visible beams at the same point of the probed samples. In these conditions, SFG and DFG coherent nonlinear processes were always generated simultaneously, provided that spatial and temporal overlap of the input beams was achieved at the interface. The angles of incidence of the visible and IR beams were 55° and 65° with respect to the surface normal, respectively. The visible beam fixed at 523.5 nm was near the maximum of the SPR of the 17±2 nm diameter 20 AuNPs. SFG/DFG signals were recorded by photomultipliers after spatial and spectral filtering through a monochromator. Data were normalized by the SFG/DFG signal of a ZnS reference crystal to take account of all laser fluctuations and absorptions on the laser paths.

The samples surface density was analyzed with AFM (Digital Instrument, DI3100) in tapping mode. The silicon tips used in our measurements have 130kHz working frequency and their curvature radius at the apex is around 10nm. This latter does not allow to obtain sufficient lateral resolution to image correctly isolated nanoparticles but it allows to identify each gold nanoparticle and give their mean diameter (172nm). Finally AFM allows to determine the number of nanoparticles grafted on the APHS and APTES surface in order to deduce the AuNps surface density used in the calculations of the Fresnel factors in the SFG

Conclusions

We have presented IR-vis SFG/DFG measurements on thiophenolfunctionalized AuNps grafted on a silicon substrate through two different silanes. This is the first experimental evidence, to the best 40 of the authors' knowledge, of a vibrational nonlinear activity in the fingerprint spectral range observed on the surface of spherical metallic nanoparticles (17±2 nm diameter). This has been made possible thanks to the sensitivity of the vibrational FG technique, based on a setup coupled to the CLIO FEL. However, the levels of M. Ortega for their scientific and technical support on the SFG optical setup and the CLIO Free Electron Laser; P. Archirel for his help in DFT calculations. The research leading to these results has received funding from the European Community's Seventh Framework Programme (FP7/2007-2013) under grant agreement 100 n.°226716.

Notes and references

a Laboratoire de Chimie Physique, Université Paris-Sud, CNRS, Bâtiment 201 Porte 2, 91405 Orsay, France. Fax: 33169153328; Tel: 3316915 3290; E-mail: christophe.humbert@u-psud

Fig. 1

 1 Fig.1 Sketch of the SFG and DFG processes on the surface of AuNps at different scales: (A) microscopic; (B) nanoscopic; (C) molecular.

Fig. 2

 2 Fig.2 Left: SFG (A) and DFG (B) spectra of thiophenol molecules (insert) adsorbed at the AuNps/APHS/Si(100) interface. Right: SFG (C) and DFG 30

I

  35we develop hereafter an analysis of the various parameters contributing to the SFG/DFG intensity.DiscussionSFG/DFG spectroscopy at the microscale The intensity I(FG), where FG stands for either SFG (ωSFG = ωvis+ 40 ωIR) or DFG (ωDFG = ωvis -ωIR), of interfaces is modelled by: 16 IVis are the intensities of the incident IR and visible laser beams, respectively; ωIR, ωvis and ωFG the IR, visible and sum/difference frequencies of the three laser beams involved in the 45 nonlinear process; IR, vis and FG the angles of incidence of the incoming and outgoing beams; n1() the refractive index of the upper medium at frequency ; χ eff (2) the effective nonlinear second order susceptibility of the probed interface, which includes Fresnel contributions and microscopic local-field corrections to the 50 molecular nonlinear second order susceptibilities χ ijk

Fig. 3

 3 Fig.3 AFM images of AuNps grafted on Si(100) through a) APHS and b) APTES silanization at different scales (10x10µm 2 and 500x500nm 2 ).

  50effects, are due to the surface density of nanoparticles, which differs as a function of the silane used to graft them on the silicon; 4) an additional enhancement factor, related to SPR specific of each sample, accounts for the difference between absolute SFG and DFG amplitudes. The great sensitivity of the technique to the 55 surface density finally illustrates the potential role of such spectroscopies in the quantitative comparison of analogous composite materials.Thiophenol IR and Raman properties at the molecular scale 60

  10 perpendicular to the phenyl ring and containing the S atom. Refs. 24 and 25 on one side, Ref. 26 on the other side disagree on the attribution of modes 3 and 4. For the latter, mode 3 is symmetric and 4 antisymmetric, vice versa for the former. The origin of the SERS and silver-phenylthiolate Raman activity at 1075cm -1 is also debated, as Ref. 26 attributes it to an increase of the activity of peak 3 whereas Refs. 25 and 27 talk about a strong shift of the position of the Raman-active mode 4. Finally, Ref. 25 reports a fourth mode at 1065cm -1 , only IR-active. Such uncertainties on the actual active modes have rather strong impact on the interpretation of SFG data, 20 either for absolute Raman and IR activities or for geometry and orientation purposes.

60Scheme 2

 2 Scheme 2 Molecules of thiophenol (left, called planar in the text) and goldbenzenethiolate (right, two different views of the same geometry, called perpendicular in the text) in their lowest energy conformations.

Fig. 4

 4 Fig.4 Top: Illustration of the vibrational activity of mode 1 (998.0cm -1 ), 45

55

 55

Table 1

 1 Fit parameters for the SFG and DFG data of Figure 2, based on Equations (2).

			SFG	DFG	SFG	DFG
			(APHS)	(APHS)	(APTES)	(APTES)
	AuNp	C  (°)	0.147 137.57	0.0630 77.7342	0.124 143.345	0.0261 81.161
	Mode 1	1 (cm -1 ) a1 (a.u.) 1 (cm -1 )	1000.2 -0.339 3.5	1001.3 -0.584 3.5	999.5 -0.289 3.5	999.9 -0.461 3.5
	Mode 2	2 (cm -1 ) a2/a1 2 (cm -1 )	1023.6 0.91 3.5	1023.3 0.91 3.5	1023.8 0.91 3.5	1022.4 0.91 3.5
	Mode 4	4 (cm -1 ) a4/a1	1070.4 0.50	1072.9 0.50	1070.8 0.50	1070.1 0.50

60

Table 2 .

 2 Fresnel factors amplitudes (ppp polarization combination) for incoming wavelengths 523.5 nm (visible) and 10 µm (IR).

	SFG	APTES	APHS	Au(111)
	|Fxxz|	0.106	0.102	0.516
	|Fxzx|	0.108	0.102	1.655 10 -2
	|Fzxx|	0.118	0.113	1.616 10 -2
	|Fzzz|	1.015	0.885	0.645
	DFG	APTES	APHS	Au(111)
	|Fxxz|	0.111	0.107	0.537
	|Fxzx|	0.113	0.107	1.723 10 -2
	|Fzxx|	0.109	0.103	1.979 10 -2
	|Fzzz|	0.941	0.805	0.790
	65

Table 3 .

 3 Summary of the parameters deduced for APHS and APTES samples. The Table shows the APHS/APTES ratios.

SFG (

res) SFG (NR) DFG (res) DFG (NR)

  that it is the case indeed for SFG (resonant and NR) and for the DFG resonant signals within 8% error, which remains compatible with error amplitudes in the experimental measurements and fit parameter estimations. Only the NR DFG value deviates from the three others. This is due to the fact that the 15 NR backgrounds in DFG suffer from their low absolute values which hamper their quantitative analysis. As for the comparison between DFG and SFG for both samples, it includes a different local field amplification factor related to frequency dependence of the plasmon resonances. The coupling between FG and SPR

		|𝛘 𝐞𝐟𝐟 (𝟐) |	1.17	1.19	1.27	2.41
		|𝛘 𝐳𝐳𝐳 (𝟐) |	1.34	1.32	1.48	2.82
		Ns		1.37
		The quantity χ zzz (2) is in principle proportional to the surface
		density of FG sources on the sample, 16 in our case thiophenol
		molecules for the resonant part and AuNps for the NR part.
		Considering that the surface density of molecules scales as the
		surface density of particles, the χ zzz (2) ratio should be directly
		proportional to the ratio of Ns quantities, provided that local field
	10	amplification is constant for a given spectroscopy. We can check
		in Table 3

20

processes

Table 4 .

 4 Experimental vibrations of free and metal-bonded thiophenol, with qualitative IR and Raman activities (from Refs. 24 and 25). w: weak,

	Mode E (cm -1 )	IR activity	Raman activity	E (cm -1 )	Raman activity
	1	1002	w	vs	1000	vs
	2	1026	s	s	1025	s
	3	1072	m	0	(1065)	0
	4	1093	s	m	1075	s

m: medium, s: strong, vs: very strong Free molecule SERS (Ag/Au) and AgSPh

  .fr 105 b Institut des Nanosciences de Paris, Université Pierre et Marie Curie, 4 place Jussieu, 75005 Paris, France. † Electronic Supplementary Information (ESI) available: [Details of all DFT calculations of the present paper]. See DOI: 10.1039/b000000x/ 110 1 T. Kawai, D. J. Neivandt and P. B. Davies, J. Am. Chem. Soc., 2000, 122, 12301; C. Humbert, B. Busson, J.-P. Abid, C. Six, H. H. Girault and A.Tadjeddine, Electrochimica Acta, 2005, 50, 3101; C. Weeraman, A. K. Yatawara, A. N. Bordenyuk, A. V. Benderskii, J. O. Pluchery, C. Humbert, E. Lacaze and B. Busson, Phys. Chem. Chem. Phys., 2009, 11, 7729. B. Busson and A. Tadjeddine, J. Phys. Chem. C, 2008, 112, 11813. A.A. Mani, Z.D. Schultz, Y. Caudano, B. Champagne, C. Humbert, L. Dreesen, A.A. Gewirth, J.O. White, P.A. Thiry and A. Peremans, J.

Acknowledgements

The authors acknowledge C. Six, A. Gayral, J.-P. Berthet and J.-95

FG intensities measured above make us confident that analogous measurements should be possible even using tabletop laser sources in the [1000-2000cm -1 ] range. For an extension towards lower energies, the FEL remains mandatory on such materials.

A careful analysis of experimental results shows that FG Provided that an absolute calibration is achieved and the method extended towards more commonly accessible frequency regions (e.g. 2500-3500cm -1 ), it becomes an interesting tool for the 60 analysis of such composite materials, complementary to microscopy (AFM, STM) and linear optical techniques (IR, Raman). This opens the door to more routine analysis of the surface chemistry of nanoparticles using nonlinear optical spectroscopies. In addition, derived methods which lower the 65 detection threshold, like heterodyne SFG or doubly resonant SFG, 33 should prove helpful for that purpose.

In the investigated samples, we have shown that grafting gold nanoparticles on silicon with APHS is more efficient than with APTES, with an increase of the surface density by around 40%, 70 confirming results of Ref. 2. Apart from this parameter, which impacts the absolute FG intensities, both grafting silanes lead to very similar structures at the molecular and nanoscopic scale, since the SFG activities of the molecular vibrations and gold nanoparticles share the same properties.

75

We intend to continue this systematic work of understanding the optical and chemical properties of complex materials based on metallic nanoparticles. There are still further steps to undertake, concerning an accurate description of the interplay between the surface plasmon of the particles and the molecular vibrations 80 mediated by nonlinear optical spectroscopies (SFG/DFG). Improving theoretical and computational methodology is mandatory to bridge the gap between molecular properties deduced from DFT calculations and a description at the nanoscale of the orientation distribution of molecules at the surface of the particles.

85

In the future, applications including more complex molecular systems adsorbed on AuNps and related, for instance, to heterogeneous catalysis and biosensors could be accurately described at the nanoscale with SFG/DFG spectroscopy. In fact, we believe that the experimental and theoretical procedures 90 developed in this work could be used to systematically monitor the surface chemical properties of such nanomaterials, whatever their sizes, shapes or organization on a given substrate.


Evelyne Martin

Pier Luca Palla

Fabrizio Cleri

Assil Bouzid

Guido Ori

Sébastien Le Roux

Mauro Boero

Carlo Massobrio

S Gabardi

M Bernasconi

On the occurrence of size effects in the calculation of thermal conductivity by first-principles molecular dynamics: the case of glassy GeTe 4

Keywords: Thermal Conductivity, Glass, First-Principles Molecular Dynamics

The thermal conductivity of a glass can be obtained by first-principles molecular dynamics provided we exploit a methodology that has been termed the approach-to-equilibrium molecular dynamics (AEMD) [1] [2, Chap. 8]. In the present work, we investigate the occurrence of size effects by comparing the thermal conductivity of two g-GeTe 4 atomic models of different sizes. This issue is far from being trivial since, in principle, size effects are not expected to occur in disordered systems beyond a few interatomic distances. For this reason, it is important to search unambiguous pieces of evidence substantiating this point. The first system of length L = 18 Å contains 185 atoms. By duplicating it along one direction we form the second system, that contains 370 atoms and features a double length L = 36 Å and an identical cross section. The thermal conductivity increases by a factor 3 from L = 18 to 36 Å, thereby approaching the experimental value. Our investigation exemplifies the crucial role of the system size to take full advantage of the AEMD methodology and bring the calculated values in better agreement with experiments.

Introduction

Molecular dynamics (MD) simulations are widely used to obtain the thermal conductivity of bulk systems, nanostructures and interfaces [START_REF] Schelling | Comparison 220 of Atomic-Level Simulation Methods for Computing Thermal 221 Conductivity[END_REF][START_REF] Lampin | Thermal Con-215 ductivity from Approach-to-Equilibrium Molecular Dynamics[END_REF][START_REF] Volz | Molecular Dynamics Simulation of Ther-223 mal Conductivity of Silicon Nanowires[END_REF][START_REF] Zaoui | Fourier-Like Conduc-226 tion and Finite One-Dimensional Thermal Conductivity in Long 227 Silicon Nanowires by Approach-to-Equilibrium Molecular Dy-228 namics[END_REF][START_REF] Schelling | Kapitza Conduc-230 tance and Phonon Scattering at Grain Boundaries by Simulation[END_REF][START_REF] Lampin | Thermal 233 Boundary Resistance at Silicon-Silica Interfaces by Molecular 234 Dynamics Simulations[END_REF]. In the literature, the forces between atoms are described by classical force fields with few exceptions. This is due to the need of extended MD trajectories to determine the thermal conductivity and obtain a fully converged autocorrelation function in the Green-Kubo approach [START_REF] Zwanzig | Time-Correlation Functions and Transport Coef-236 ficients in Statistical Mechanics[END_REF], or to reach the steady state in the direct method [START_REF] Schelling | Comparison 220 of Atomic-Level Simulation Methods for Computing Thermal 221 Conductivity[END_REF]. However, it is highly desirable to go beyond the use of classical force fields, in particular in the case of glasses, since for disordered structures the use of quantitative models has found to be crucial [START_REF] Massobrio | Structural Properties of Amor-239 phous GeSe 2[END_REF]10,11]. This is specifically relevant for the class of chalcogenide glasses, for which empirical potentials often fail to quantitatively describe their structure and chemical interactions [12]. Recently [13], we have demonstrated that it is possible to determine the thermal conductivity of a glass, g-GeTe 4 , Email address: evelyne.lampin@univ-lille.fr (Evelyne Martin) [START_REF] Zhang | Thermal Conductivity and Specific Heat of Bulk Amorphous Chalcogenides Ge 20 Te 80-x Se x (x=0[END_REF]. However, it 39 is known that the maximum heat carrier free paths can 40 be considerably larger than the average value extracted 41 from the thermal conductivity using the kinetic theory.

It is therefore desirable to ascertain whether or not the thermal conductivity of g-GeTe 4 features a length dependence. This is exactly the objective of the present work. As a byproduct, one would like to understand whether the underestimate thermal conductivity found in [13] with a system of 185 atoms can be at least partially overcome by resorting to a longer system in the direction of the heat flux.

The article is organized as follows. In a first section, the two atomic models of g-GeTe 4 and the MD methodology are presented. In a second section, the MD trajectories are analyzed, and the thermal conductivity of the two systems are extracted and discussed. Then, we critically analyze our results before drawing a final set of conclusions.

Method

The atomic model of g-GeTe The atomic structure of the two systems is identical, as evidenced by the comparison of the total pair correlation functions g(r) (Fig. 2).

The details of the CPMD calculations are as follows. The Becke, Lee, Yang and Parr (BLYP) [START_REF] Becke | Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior[END_REF][START_REF] Lee | Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density[END_REF] exchange-correlation functional is used, with a norm-conserving pseudo-potentials [START_REF] Troullier | Efficient Pseudopotentials for Plane-Wave Calculations[END_REF] describing the valence-core interaction. Long-range dispersion forces K is reached faster for the smaller system (Fig. 3).

∆T (K) Time (ps) L=18 Å L=36 Å L=18 Å L=36 Å L=18 Å L=36 Å L=18 Å L=36 Å L=18 Å L=36 Å L=18 Å L=36 Å L=18 Å L=36 Å MD τ=4.0ps τ=7.0ps
Figure 4: Time evolution of the temperature difference in phase 2 for the two systems and adjustment by a decaying exponential (the decay time is mentioned in label). [START_REF] Lampin | Thermal Con-215 ductivity from Approach-to-Equilibrium Molecular Dynamics[END_REF]:

108 κ = L 2 4π 2 C V .ρ τ (1

  from a density functional (DFT) based first-principles 20 MD (FPMD). The calculations have been performed by 21 applying a methodology that we recently developed, the 22 approach-to-equilibrium MD (AEMD) [1] [2, Chap. 8]. 23 In this methodology, a transient regime is created and 24 exploited to calculate the thermal conductivity in a way 25 that proved to be faster and more effective than conven-26 tional approaches. 27 Past experience of application of the AEMD method-28 ology provided evidence of the sensitivity to the system 29 size of thermal conductivities obtained by AEMD. In 30 particular, it has been shown that in crystalline Si (cSi), 31 the thermal conductivity changes with the length of the 32 system even up to 1.2 µm. This length dependence is 33 intimately related [14] to the existence of long phonon 34 free paths in cSi (larger than 10 µm) [15] [16]. In a 35 disordered material like g-GeTe 4 , the free paths of heat 36 carriers are expected to be very short. At low temper-37 ature (T ≈ 130 K), a mean free path (MFP) of 4-5 Å 38 has been measured by Zhang et al.

  4 has been obtained by Car-Parrinello molecular dynamics (CPMD) [18] using the CPMD code 1 . This model yields structure factors and pair correlation functions in very good agreement with measurements [10]. In our first implementation of AEMD to the case of g-GeTe 4 [13] the cubic simulation box contained N = 185 atoms with a side L =18 Å.Here the cubic box is doubled in the direction of thermal transfer to form a rectangular box containing N = 370 atoms, with a length L = 36 Å and a cross section S = 18 × 18 Å 2 (Fig.1).

Figure 1 :

 1 Figure 1: The atomic models containing respectively 185 and 370 atoms.

Figure 2 :

 2 Figure2: Pair correlation function g(r) for the system of length L = 18Å containing 185 atoms and the system of length L = 36Å containing 370 atoms. The number of first neighbors as obtained from the integration of the first peak is equal to 3.8 in both cases.

101Figure 3 :

 3 Figure 3: Time evolution of the hot block (red line) and cold block (blue line) during the two AEMD phases. Upper graph: L = 36 Å. Lower graph: L = 18 Å.

Figure 5 :

 5 Figure 5: Thermal conductivity as a function of the length and for 3 values of the initial difference in temperature between hot and cold blocks (∆T 0 = 200, 100 or 50 K) for atomic models containing respectively 185 and 370 atoms.

134[ 5 ]

 5 There are several examples of length dependence 135 of the thermal conductivity obtained using the AEMD 136 methodology. In bulk silicon and germanium[START_REF] Lampin | Thermal Con-215 ductivity from Approach-to-Equilibrium Molecular Dynamics[END_REF] [14], 137 the thermal conductivity changes on the entire range of 138 lengths studied, from tens of nanometers up to 1.2 µm. 139 In α-quartz [14], the thermal conductivity varies up to ≈ 140 100 nm before reaching a plateau. In silicon nanowires 141 , the thermal conductivity always reaches a plateau 142 only after a length that increases with the nanowire di-143 ameter. In silicon membranes [27], the same behavior 144 is observed, with the plateau obtained at shorter length 145 in presence of holes. All these cases were studied using 146 classical force fields. 147 We have shown [14] that this behavior is due to the 148 truncation at length L of the distribution of phonon 149 mean free paths. In crystals, phonons encounter scatter-150 ing events that make their lifetime finite, or the distance 151 they cover finite. The average mean free paths Λ ph is ν being the sound velocity. Nevertheless, phonon mean 155 free paths can be distributed up to a maximum Λ max con-156 siderably larger than Λ ph . For example in silicon, Λ ph = 157 40 nm and Λ max > 10µm . 158 In amorphous materials and glasses, heat carriers are 159 not expected to cover long distances due to the disor-160 der of the atomic structure that could be seen schemati-161 cally as a distribution of defects beyond the first neigh-162 bor shell. Indeed, Zhang et al. [17] have estimated Λ ph 163 to 4-5 Å from their measurements of the thermal con-164 ductivity of g-GeTe 4 . In amorphous silicon studied by 165 classical force fields [14], the thermal conductivity does 166 not exhibit any length dependence above 5 nm. 167 The variation of thermal conductivity obtained in the 168 present work, ranging from L = 18 to 36 Å could be 169 due to the fact that L is still smaller than the maxi-170 mum phonon mean free path Λ max . This conjecture can 171 only be validated by looking for the saturation in the 172 calculated values of the thermal conductivity. This re-173 quires the calculation of the thermal conductivity on a 174 third, even more extended system. It remains true that 175 the present calculations have the pleasing effect of ap-176 proaching the AEMD value to the experimental coun-177terpart. Also, it has to be observed that the size scaling of the CPMD method is prohibitive and by no means 179 linear. For this reason, it might be of interest to resort to 180 alternative schemes (as the so-called second generation 181 CPMD[START_REF] Kühne | Second Generation Car-Parrinello Molecular Dynamics[END_REF]) in order to speed up the achievement of the 182 planned calculation for a third, larger system size. 183 5. Conclusion 184 The thermal conductivity of a glass, g-GeTe 4 has 185 been calculated from first-principles molecular dynam-186 ics. In particular, the dependence on the system length 187 has been the focus of the work, since the AEMD 188 methodology has shown such a dependence in many 189 other materials. Although the disordered nature of g-190 GeTe 4 and an average phonon mean free paths mea-191 sured to 4-5 Å in the literature were in favor of an ab-192 sence of length dependence, we have obtained an in-193 crease by a factor 3 of the thermal conductivity when 194 the cell length along the heat path changes from 18 to 195 36 Å. This amounts to a calculated value much closer 196 to experimental data. Therefore, it appears that such a 197 length dependence has some physical origins that can 198 be traced back to the existence of maximum mean free 199 paths larger than the system size. Further investigations 200 for even larger sizes will be instrumental to demon-201 strate that the combination of FPMD with the AEMD 202 approach is able to provide fully size-converged, reli-203 able values for the thermal conductivity, at least in the 204 case of a prototypical disordered system.

  205

See http://www.cpmd.org, copyright

2000-2017 jointly by IBM Corporation and by Max Planck Institute, Stuttgart.

Acknowledgements

estimated from the kinetic theory following: 213