Vibrational dynamics of CO on Pd(111) in and out of thermal equilibrium

Vibrational dynamics of CO on Pd(111) in and out of thermal equilibrium Ra´ul Bomb´ın,1, 2, 3, ∗A. S. Muzas,2Dino Novko,4J. I˜naki Juaristi,1, 2, 5, †and Maite Alducin2, 5, ‡ 1Departamento de Pol´ımeros y Materiales Avanzados: F´ısica, Qu´ımica y Tecnolog´ıa, Facultad de Qu´ımicas (UPV/EHU), Apartado 1072, 20080 Donostia-San Sebasti´an, Spain 2Centro de F´ısica de Materiales CFM/MPC (CSIC-UPV/EHU), Paseo Manuel de Lardizabal 5, 20018 Donostia-San Sebasti´an, Spain 3Departament de F´ısica, Universitat Polit`ecnica de Catalunya, Campus Nord B4-B5, E-08034, Barcelona, Spain 4Centre for Advanced Laser Techniques, Institute of Physics, Bijeniˇcka 46, 10000 Zagreb, Croatia 5Donostia International Physics Center (DIPC), Paseo Manuel de Lardizabal 4, 20018 Donostia-San Sebasti´an, Spain (Dated: July 18, 2023) Using many-body perturbation theory and density functional perturbation theory, we study the vibrational spectra of the internal stretch (IS) mode of CO on Pd(111) for the bridge and hollow adsorption structures that are experimentally identified at 0.5 ML coverage. Our theoretical treatment allows us to determine the temperature dependence of the IS vibrational spectra under thermal conditions as well as the time evolution of the non-equilibrium transient spectra induced by femtosecond laser pulses. Under thermal conditions (i.e., for equal electronic Teand phononic Tl temperatures), the calculated lifetimes at 10-150 K are mostly due to nonadiabatic couplings (NC), i.e., first-order electronic excitations. As temperature increases, also the contribution of the secondorder electron mediated phonon-phonon couplings (EMPPC) progressively increases from 25% at low temperatures to 50% at 300 K. Our calculations for the laser-induced non-equilibrium conditions comprise experimental absorbed fluences of 6-130 J/m2. For fluences for which Te>2000 K, the transient vibrational spectra are characterized by two different regimes that follow the distinct time-evolution of Teand Tland are respectively dominated by NC and EMPPC processes. At lower fluences, the initial fast regime becomes progressively negligible as Tedecreases and only the steady second regime remains visible. Qualitatively, all these spectral properties are common to the both adsorption structures studied here. I. INTRODUCTION The internal stretch (IS) mode of polar molecules adsorbed on metal surfaces, being decoupled in energy from the rest of vibrational modes of the system, can be monitored during a reaction in time-resolved vibrational spectroscopy experiments [1]. In these kind of setups, the adsorbate dynamics is initiated by a pump laser pulse and the frequency shift and the linewidth of the IS mode are tracked with femtosecond time resolution. It has been suggested that these changes can provide insight into the initial steps of the induced reaction [2–5], including the transient charge flow that the pump laser may induced in organic adsorbates [6,7]. Commonly, after the pump pulse heats the surface, the IS mode exhibits a rapid initial redshift followed by a slower (but also fast) recovery (e.g., CO/Ru(0001) [3], NO/Ir(111)[8], CO/Ir(111) [9], CO/Pt(111) [10–13], CO/Cu(100) [5,14], CO/Cu(110) [15]). However, for the case of CO molecules coordinated to ruthenium tetraphenylporphyrin on a Cu(110) surface a blueshift has been reported [16]. Closely related, thermal heating experiments represent another example in which an increase of the linewidth and a frequency redshift in the IS mode is expected as ∗[email protected] †[email protected] ‡[email protected] temperature rises. The latter has been observed for CO/Cu(100) [14,17–19], CO/Pt(111) [10,20,21], CO/Ru(0001) [22], CO/Ru(0101) [23], as well as for CO/Ir(111) and NO/Ir(111) [9]. However, Persson et al. reported a thermal blueshift for CO/Ni(111) [24] when varying the temperature from 25 to 300 K. Under similar temperature changes, a blueshift has been also observed in a high coverage phase, with a rather more complex structure, of CO/Pd(100) [17]. All these findings suggest that an accurate microscopic description of the different mechanisms that participate in the vibrational dynamics of polar molecules on the different metallic surfaces is needed in order to give a correct interpretation to experimental measurements. In our previous work we studied the frequency shift that is induced on the IS mode of the 0.5 ML of CO on a Pd(111) surface under pump-probe conditions [25]. There, we showed that due to the sharp shape of the electron density of states (DOS) around the Fermi level, the electron-phonon (e-ph) coupling is effectively screened [25], giving place to an anomalous transient blueshift in the initial few hundreds of femtoseconds upon arrival of the pump laser pulse. The time evolution of the observed blueshift resembles that of Ref. [16], what may hint at a similar mechanism behind their origin. Here we extend that study to the whole range of thermal and pump-probe non-equilibrium experimental conditions, characterizing not only the frequency-shift but also the changes that are induced in the linewidth. As there has been a long debate about the adsorption structure of the 0.5 ML of CO on Pd(111), we study arXiv:2304.10845v4 [cond-mat.mtrl-sci] 17 Jul 2023 2 both the fcc-hcp and bridge structures that have been experimentally observed [26,27]. All the calculations presented here are performed with the methodology of Refs. [19,28], which is based on many-body perturbation theory (MBPT) and density functional perturbation theory (DFPT). This robust theoretical framework has been successfully applied to explain the transient subpicosecond vibrational spectral changes measured for CO/Cu(100) [19,28]. The paper is organized as follows. In Sec. II we describe the method that we use to evaluate nonadiabatic effects in the IS vibrational mode of the CO adlayer. In the MBPT framework this is done in terms of the phonon self-energy. In Sec. III we characterize the electronic and phononic properties of the system under study. Further details for the Pd(111) surface and CO/Pd(111) system that are not essential for the discussion presented in the main text are provided in Appendices Aand B. Our results for the frequency shift and linewidth of the CO IS mode under thermal and non-equilibrium conditions are shown in Secs. IV and V, respectively. A comparison to other theoretical and experimental studies on the dynamics of the IS mode of polar molecules in other transition metal surfaces is presented in Sec.VI. Finally, the main conclusions are summarized in Sec. VII. II. METHODS Electronic structures, dynamical matrices, and phonon perturbation potentials are calculated with density functional theory (DFT) and DFPT using the quantum espresso first principles package [29,30]. The exchangecorrelation interaction is described with the Bayesian error estimation functional with van der Waals correlation (BEEF-vdW) [31], which is specifically designed for treating surface science problems. We use the ultrasoft pseudopotentials by A. dal Corso [32] with a kinetic energy cutoff of 1360 eV and a Gaussian smearing of 0.27 eV for the electronic state occupancies. The electronic states and charge densities are evaluated by sampling the Brillouin zone with a Γ-centered 12×12×1 Monkhorst-Pack (MP) mesh [33]. The CO/Pd(111) system is modeled with a c(2×√3)rect surface unit cell that includes four Pd-layers with 12.9 ˚ A of vacuum. The desired 0.5 ML of CO coverage is achieved by including two CO molecules in the supercell. All the atoms in the supercell are relaxed until the Hellmann-Feynman forces on each atom are below 2×10−5Ry/aB(∼2×10−4eV/˚ A). See Appendices Aand Bfor more details about these structures. The dynamical matrices and perturbation potentials are calculated on a Γ-centered 6×6×1q-grid. We use the Electron-phonon Wannier (epw) code to evaluate the e-ph matrix elements [34,35]. From the electronic states and charge density computed in the Γcentered k-grid we construct a set of 108 maximally localized Wannier functions (MLWFs), using the selected columns of density matrix procedure [36]. The e-ph matrix elements are calculated in the kand qgrids that have been used to sample the electronic and phononic Brillouin zones. Finally, taking advantage of the properties of the MLWFs [37–40], the e-ph matrix elements are interpolated on finer kand q-grids to obtain the desired accuracy in the properties that we study. The vibrational spectra of CO on Pd(111) are studied within MBPT by calculating the phonon self-energy due to e-ph coupling πλ(q, ωq,λ) (λ,q, and ωq,λ are the index, momentum, and energy of the phonon mode, respectively). The corresponding phonon linewidth is determined by taking the imaginary part of this quantity as, γq,λ =−2Imπλ(q, ωq,λ), while the renormalization of the phonon frequency is obtained from the real part as, ω2 q,λ =ω2 A+2ωARe [πλ(q, ωq,λ)−πλ(q,0)], with ωAbeing the adiabatic phonon frequency obtained in the DFPT calculation [41,42]. As nonadiabatic corrections to the phonon energy are usually small compared to the adiabatic energy (i.e., ωq,λ −ωA≪ωA), the renormalized phonon frequency can be approximated by ωq,λ ≈ωA+Re [πλ(q, ωq,λ)−πλ(q,0)]. Here we are interested in the vibrational spectra probed by infrared light. Thus, only the long wavelength part of the phonon selfenergy, denoted hereafter as πλ(ω0,λ)≡πλ(q≈0, ω), will be taken into account. Following previous works [19,28], we consider that the phonon self-energy consists of firstand secondorder terms in the e-ph coupling, πλ(ω0,λ)=π[1] λ(ω0,λ)+ π[2] λ(ω0,λ). Later it will be clear that π[1] λand π[2] λcorrespond to dominant interband and intraband contributions, respectively. The first order accounts for nonadiabatic coupling (NC) to electron-hole pairs, while the second order corresponds to the so-called electron mediated phonon-phonon coupling (EMPPC). This theory has been able to reproduce and explain the thermal [19] and ultrafast [28] vibrational spectra of CO/Cu(100) measured in Refs. [14] and [5], respectively. Nonadiabatic coupling.— The expression for the firstorder NC term that exclusively accounts for the electronhole pair (de)excitations reads [19,41,42] π[1] λ(ω0,λ;Te) = X µµ′kσgµµ′ λ(k,0) 2f(ϵµk)−f(ϵµ′k) ω0,λ +ϵµk−ϵµ′k+iη , (1) where µ,k, and ϵµkare the electron band index, momentum, and energy, respectively; gµµ′ λ(k,q) are the e-ph matrix elements; and the summation over σaccounts for the spin degree of freedom. We have introduced the Fermi-Dirac distribution function f(ϵµ,k)= 1/(eβ[ϵµk−µ(Te)]+1), where β=1/(kBTe), kBis the Boltzmann constant, Teis the electronic temperature, and µ(Te) is the chemical potential that equals the Fermi level at Te=0, i.e., εF=µ(0). The chemical potential, is calculated self-consistently by solving numerically the equality, Ne=RDOS(ϵ)f(ϵ, Te, µ(Te))dϵ, that assures conservation of the number of electrons Ne. Strictly speaking ηis an infinitesimal parameter and one obtains the exact 3 first-order phonon linewidth in the limit η→0. However, ηis usually fixed to a small value representing the broadening of the electronic states due to scattering processes [41,43]. In our case we fix it to a small, physically motivated value of 30 meV [44,45]. Note that the terms µ=µ′vanish so that π[1] λonly contains interband contributions. As a final remark, note that the expression for the linewidth derived from Eq. (1) is equivalent to that obtained from the Fermi golden rule that has been used by other authors [46–54]. Electron mediated phonon-phonon coupling.— Early studies on the vibrational relaxation of polar molecules on metal surfaces already accounted for two-phonon processes that involved anharmonic coupling between the adsorbate vibrational modes [48]. However, the fact that the IS mode is decoupled in energy from the other vibrational modes makes this mechanism quite inefficient. Nonetheless, as shown in Ref. [19], the presence of electron-hole pairs in the metal surface may compensate this mismatch in energy. The IS mode can excite electron-hole pairs that experience a second scattering process with a low energy phonon mode. In the MBPT language the EMPPC is accounted by evaluating the second order intraband phonon-self-energy that reads [19,25,28] π[2] λ(ω0,λ;Te;Tl)=−X µµ′kσ,λ′k′|gµµ λ(k,0)|2 1−gµ′µ′ λ(k′,0) gµµ λ(k,0) !gµµ′ λ′(k,q′) 2 ×X s,s′=±1 f(ϵµ,k)−f(ϵµ′,k′−s′sωq′,λ′) ϵµ,k−(ϵµ′,k′−s′sωq′λ′) snb(sωq′λ′)+f(s′ϵµ′,k′) ω0,λ ω0,λ +iη+s′(ϵµ,k−ϵµ′,k′)+sωq′λ′,(2) where q′=k′−kand nb(ωq,λ)=1/(eβωqλ−1) is the BoseEinstein distribution, with β=1/(kBTl) and Tlthe lattice temperature. Equation (2) couples the studied (q≈0, λ) mode with other modes (q′,λ′) via electron-hole pairs. In the results presented in this work, the vertex correction [second term in the square bracket of Eq. (2)] is neglected to reduce the computational cost. Nonetheless, it has been suggested that the contribution of the vertex correction is small [55,56]. In Eqs. (1) and (2) the electronic temperature Teenters the Fermi-Dirac distribution, while the lattice temperature Tlis included via the Bose-Einstein distribution that only appears in Eq. (2). This allows us to study the frequency shift and linewidth of the IS mode of CO adsorbed on the Pd(111) surface under different temperature conditions. In section IV we consider thermalized conditions, that is, the Te=Tlcase. Instead, in section V we consider the ultrafast non-equilibrium conditions that are created when the system is excited with a pump femtosecond laser pulse. The theoretical treatment of the dynamical evolution of the system under the effect of an ultrafast pump pulse usually relies on a model for the degrees of freedom that are excited. Among these approaches one of the simplest and most popular techniques is the two-temperature model (TTM) [57,58]. It assumes that the energy of the pump pulse is initially absorbed by the electronic bath on the surface, creating hot electrons. Those hot electrons then thermalize transferring energy to other electrons or to the lattice via electron-electron and e-ph interactions, respectively. This approximation needs as input some information about the system, namely, the electronic and lattice heat capacities, the electronic thermal conductivity, the e-ph coupling constant and their dependence with temperature. However, the lack of experimental data for these quantities for the high electronic temperatures that are generated in pump-probe experiments forces to rely on some approximations. In previous works these parameters have been included in the model under some approximations, such as making use of low temperature analytical expressions and extrapolating the available experimental data to the high temperature regime. Alternatively, these quantities can be evaluated fully from first principles [28,58]. The latter has been done recently by Li and Ji [59] for different transition metals including Pd in a wide rage of temperatures. In this work, we make use of those results to construct a TTM. In Appendix Cwe provide further details about the TTM that we employ. III. THE SYSTEM The adsorption structure of CO on Pd(111) has stood as a puzzle for a long time. For the 0.5 ML coverage, it has been determined that two c(4×2)-2CO arrangements coexist at temperatures in the range of hundreds of Kelvin [26,27]: the hollow configuration CO(hollow)/Pd(111), in which the CO molecules adsorb on fcc and hcp hollow sites, and the bridge configuration CO(bridge)/Pd(111) with CO at bridge sites. To give insight on how much the adsorption structure of the CO adlayer may affect the vibrational dynamics, we study the two situations independently. Further details about the atomic and electronic structures for these two configurations are provided in Appendix B. The electronic band structures show only small differences between the two 4 configurations. On the contrary, clear differences arise in the vibrational spectra. For this reason, we focus on the vibrational spectra and the coupling between electronic and phononic degrees of freedom in the main text and leave the discussion about the electronic band structures for the appendices. Figures 1(a) and (b) show the phonon DOS for the hollow and bridge configurations, respectively. The solid blue area represents the total DOS, while the orange area corresponds to the projection onto the adsorbate normal modes. In each projected DOS, we label the energy regions in which the IS and the other vibrational modes appear, namely, external stretch (ES), frustrated rotations (FR), and frustrated translations (FT). The label “mixed” indicates that different adsorbate modes overlap in that energy range. The IS mode, located at ω≈1800−1900 cm−1, is decoupled in energy from the rest of vibrational modes (notice that the x-axis is shortened in the range ω∈[450,1600] cm−1to facilitate visualization). The Pd modes appear in the range ω∈[0,230] cm−1 and overlap with the CO FT modes. At intermediate energies, the vibrational spectra is dominated by FR and ES modes. It is precisely in this range where clear differences between the hollow and bridge configurations emerge. In the former case, the peak at ω=351 cm−1corresponds to FR, while the double-peak structure around ω=300 cm−1is mainly a mixture of FR and ES modes. In contrast, for the bridge configuration three peaks are clearly identified in the region ω∈[250,415] cm−1that correspond, from the highest to lowest energy, to FR, ES, and a mixture of FT, FR, and ES modes (labeled as mixed) in the figure. The Eliashberg spectral function α2F(ω), which is shown in Fig. 1(c) for the hollow (red) and bridge (green) structures together, provides information on the effect of the e-ph coupling in the phonon spectra. The Eliashberg function is here calculated as α2F(ω)= 1 πN(ϵF)X qλ Im hπ[1] λ(q, ωqλ)i ωqλ δ(ω−ωqλ),(3) where N(ϵF) is the DOS at the Fermi level, δ(x) is the Dirac delta function, and π[1] λ(q, ωqλ) = X µµ′kσgµµ′ λ(k,q) 2f(ϵµk)−f(ϵµ′k+q) ωqλ+ϵµk−ϵµ′k+q+iη (4) is the q-dependent bare phonon self-energy for the ωqλ mode. Notice that Eq. (1) is directly obtained by setting q=0 in the above expression. Importantly, in order to obtain α2F(ω) one needs to calculate Imπ[1] λin doubledelta approximation [41,42,60]. Following common practice, the results in Fig. 1(c) correspond to evaluating α2F(ω) at T=0. Comparing for each structure the shape of α2F(ω) to that of F(ω) in the region where the Pd modes appear (ω∈[0,230] cm−1) and in the CO modes region, one concludes that the coupling to the latter modes is weaker. This is inferred by 0.00 0.20 0.40 F ( ) Pd + FT Mixed FR IS (a) 0.00 0.20 0.40 F ( ) Pd + FT Mixed ES FR IS (b) 0 150 300 450 (cm 1) 0.00 0.06 0.12 0.18 2 F ( ) 1600 1750 1900 (c) FIG. 1. Phonon DOS for CO/Pd(111) in (a) hollow and (b) bridge configurations. Total DOS in blue and DOS projected onto the CO modes in orange. (c) Eliashberg function for CO/Pd(111) in the hollow (red) and bridge (green) configurations. A Gaussian smearing of 0.5 meV (∼4 cm−1) is employed. noting that the relative height of the CO modes to the Pd modes in α2F(ω) is smaller than in F(ω). The coupling as defined with α2F(ω) is particularly small for the CO IS mode in which we are interested in this work, however, for our purposes this comes from the high frequency of the IS mode and not the small value of Fermi-surfaceaveraged phonon linewidth Imπ[1] λ. IV. THERMAL CONDITIONS In this section we analyze the temperature dependence of the IS mode spectra under thermal conditions (Te=Tl=T). Figure 2(a) shows the IS frequency shift ∆ω(T)=ω(T)−ω(0) in the range T∈[0,300] K. As expected for moderate temperatures (kBT≪ω0λ) [19], the first order NC gives a negligible contribution (∆ω[1] < 0.1 cm−1, not shown). Thus, the monotonically increasing frequency redshift that we observe in the figure as Tincreases is driven by the EMPPC processes and, specifically, by the temperature dependence of the BoseEinstein distribution [see Eq. (2)]. The redshift is slightly larger for the bridge structure, being the difference in 5 ∆ω(T) between both structures 0.60 cm−1at T=300 K. The total IS linewidth γ(solid lines) and its dependence on Tare shown in Fig. 2(b) for the two structures together with their corresponding first-order NC contributions γ[1] (dashed lines). The results for both structures are qualitatively very similar, the difference being the slightly larger IS phonon linewidth of the bridge structure. In each case, the NC term clearly follows the mentioned constant behavior that is expected for moderate temperatures. This term dominates the linewidth at very low temperatures (T <100 K) but there is also a substantial EMPPC contribution of about 25% that cannot be neglected. As temperature increases, the EMPPC term γ[2](T), which is responsible for the dependence of γon T, increases with Tand becomes equally important. Similarly to CO/Cu(100) [19], the dependence of γon T goes as T3at very low temperatures and changes to linear in the range T∈[160,300]. As shown in Table I, the main contribution to γ[2](T) comes from the coupling to the surface vibrational modes, followed by the coupling to the FR, ES, and FT. Furthermore, we also find that the IS dephasing (i.e., coupling of the IS phonon at q=0 with all other IS phonons at q=0) is almost negligible. The relative importance of each of these contributions to γ[2](T) is in accordance with the discussed properties of the Eliashberg spectral function. In other words, the main contributions to the EMPPC processes come from those modes that, as described by α2F(ω), couple more efficiently to electron-hole pairs. In contrast, for CO/Cu(100) [19] it was found that the coupling to the surface modes is weaker than to ES and FR modes. Altogether, the results in Fig. 2(b) show that the IS linewidths of the hollow and bridge structures vary, respectively, within the range γ=3.55−5.1 cm−1and γ=3.60−5.4 cm−1as Tincreases from 10 to 300 K. Intriguingly, very similar linewidth increase with temperature was measured for CO/Pt(111) [20]. The corresponding lifetimes (τ=ℏ/γ) at 10−150 K for the hollow (1.70−1.04 ps) and bridge (1.49−0.98 ps) structures are comparable to the ones measured in Pt(100) TABLE I. First-order NC and second-order EMPCC contributions to the IS mode linewidth γ(T) for the hollow and bridge configurations under two thermal conditions, T=100 K and 300 K. The mode-resolved contributions to EMPPC (IS dephasing, FT, ES, FR, and Pd phonons) are given in the bottom part of the table. All linewidths are in cm−1. hollow bridge T=100 K T=300 K T=100 K T=300 K γ[1] 2.55 2.53 2.83 2.80 γ[2] 1.25 2.6 1.16 2.60 IS 0.025 0.030 0.021 0.021 FT 0.067 0.130 0.038 0.082 ES 0.23 0.38 0.17 0.27 FR 0.35 0.54 0.27 0.40 Pd 0.58 1.52 0.67 1.83 7 6 5 4 3 2 1 0 (cm 1) (a) Hollow Bridge 0 40 80 120 160 200 240 280 320 360 T (K) 0 1 2 3 4 5 6 (cm 1) (b) FIG. 2. Temperature dependence of the CO IS vibrational mode for CO(hollow)/Pd(111) (blue squares) and CO(bridge)/Pd(111) (orange squares). (a) Total frequency shift ∆ω(T). (b) Total linewidth γ(T). For each configuration the dashed-line connected circles show the corresponding NC contribution γ[1]. (1.5±0.5 ps [61]), Pt(111) (2.2 ps [62]) and Cu(100) (2±1 ps [63]). The exception is the Au(111) surface where the CO is weakly physisorbed and the measured lifetime is larger (49±3 ps [64,65]). Note also that our values of the first-order NC linewidth γ[1] are close to the theoretically obtained ones in Ref. [52] for various surfaces. Unfortunately, up to our knowledge there is no experimental data available for the CO/Pd(111) system. V. NON-EQUILIBRIUM CONDITIONS In this section we assume that the CO/Pd(111) system, which is initially in thermal equilibrium at T=100 K, is pumped with a 450 nm laser pulse with 100 fs duration that hits the surface at t=0.1 ps. The considered absorbed fluences are in the range [6,130] J/m2 that comprises the experimental values used in existing femtosecond pump-probe vibrational spectroscopy experiments. As explained in Sec.II, the non-equilibrium conditions created by the pump pulse are described with the TTM, from which we calculate the time-dependent electronic and phononic temperatures entering the firstand second-order phonon self-energy expressions [Eqs.(1) and (2)]. Following experiments, we calculate the measured transient frequency shift due to e-ph coupling as ∆ω(t)≈Re[πλ(ω;t)]−Re[πλ(ω; 0)]. Obviously, the contributions to ∆ω(t) from the firstand second-order terms 6 0.0 1.0 2.0 3.0 T (103 K) (a) 6.5 4.5 2.5 0.5 1.5 (cm 1) (b) 01234 t (ps) 0.0 0.4 0.8 1.2 1.6 2.0 (cm 1) (c) FIG. 3. Transient changes in the CO IS vibrational spectra of CO(hollow)/Pd(111) induced by a 450 nm laser pump pulse with 100 fs duration, absorbed fluence 40 J/m2, and peak intensity hitting the system at t=0.1 ps. The initial temperature is 100 K. (a) Electron Te(t) (blue) and lattice Tl(t) (orange) temperatures, (b) frequency shift, and (c) linewidth change as a function of time. In (b) and (c), the black line is the sum of the first order NC term (blue dashed) plus the second order EMPPC term (orange dashed). Red and green dotted lines are the CO modes and surface modes contributions to EMPPC, respectively. will be ∆ω[1],[2](t)≈Re[π[1],[2] λ(ω;t)]−Re[π[1],[2] λ(ω; 0)]. In the same manner, the transient linewidth change and corresponding firstand second-order contributions are obtained as ∆γλ=−2Imπλ(ω;t)+2Imπλ(ω; 0) and ∆γ[1],[2] λ=−2Imπ[1],[2] λ(ω;t)+2Imπ[1],[2] λ(ω; 0), respectively. We start considering the case in which the pump laser deposits an absorbed fluence of 40 J/m2in the CO(hollow)/Pd(111) system. Figure 3(a) shows the time evolution of Te(t) and Tl(t). The time interval of the figure extents a few picoseconds after the instants at which Teand Tlare maxima (≈0.15 ps and ≈2 ps, respectively). Figure 3(b) shows the transient frequency shift ∆ω(t) that is induced in the CO IS mode (black line). It exhibits an initial blueshift that reaches a maximum value of 1.5 cm−1that coincides in time with the maximum electronic temperature at t≈150 fs. After the first few hundreds of femtoseconds the blueshift vanishes and gives place to a redshift that reaches its maximum at t≈1.7 ps (∆ω=−6.5 cm−1). The same qualitative behavior of the IS frequency shift was obtained for the bridge configuration [25]. As shown in that case, the analysis of the ∆ω[1](t) (blue dashed curve) and ∆ω[2](t) (orange dashed curve) contributions reveals that the initial blueshift is created by the first-order NC processes, while the subsequent redshift is due to the EMPPC processes. It is clear from Fig. 3(a) that ∆ω[1](t) and ∆ω[2](t) follow the time evolution of Teand Tl, respectively. As demonstrated in Ref. [25], the nature of the predicted IS blueshift in CO/Pd(111) is purely electronic and related to the abrupt change in the Pd(111) electron DOS around the Fermi level. In this respect, it contrasts with the blueshift observed under thermal conditions that is explained in terms of couplings to other CO modes [15,24]. The detailed analysis that supports the electronic nature of the blueshift in the hollow structure follows the arguments of our previous work for the bridge structure [25] and are not repeated here. Finally, Fig. 3(c) shows that the transient linewidth change ∆γ(t) at F=40 J/m2starts to exhibit the two distinct regimes dictated by the different time evolution of Te(t) and Tl(t). The NC term (blue dashed line) is responsible for the initial fast increase, while the slower but dominant change in the linewidth is due to the strong Tldependence of the EMPCC term (orange dashed line). Quantitatively, the largest contribution of the EMPCC processes ocurring at t>1.5 ps (∼2.0 cm−1) is more than twice the maximum NC contribution of about 0.8 cm−1achieved at t∼100 fs. The effect that the absorbed fluence has on the transient changes of the IS spectra for the hollow and bridge structures is shown in Fig. 4for the following fluences: F= 6, 19, 40, 80 and 130 J/m2. The corresponding Te(t) and Tl(t) curves are represented in panels (a) and (b), respectively. The results for the transient frequency shift of the IS mode in the hollow configuration are shown in Fig. 4(c). At the lowest fluence (F=6 J/m2), the maximum electronic temperature of about 1000 K seems insufficient to induce a non-negligible NC contribution. In this case, only a small redshift in the frequency, due to EMPP processes, is observed as the lattice temperature increases. For intermediate fluences (F= 19 J/m2and the above discussed 40 J/m2), a blueshift appears during the first 300 fs that is progressively followed by a steady redshift of larger magnitude than the initial blueshift. Notably, for larger fluences (F= 80 and 130 J/m2) the transient blueshift is suppressed in a faster manner than for the intermediate fluences. This behavior is a consequence of the nonmonotonic dependence of ∆ω[1] on Te that will be discussed below. Figure 4(d) shows the results for the transient frequency shift in the bridge structure. The time evolution of ∆ωfor the three lowest fluences (F= 6, 19, and 40 J/m2) is qualitatively the same in the two structures, although the initial blueshift, which is due to the NC term, is quantitatively more pronounced in the bridge than in the hollow configu- 7 1 4 7 Te (103 K) (a) 1 4 7 Tl (102 K) (b) 15 12 9 6 3 0 3 (cm 1) (c) 01234 t (ps) 0 1 2 3 4 5 6 (cm 1) (e) 15 12 9 6 3 0 3 (cm 1) (d) 01234 t (ps) 0 1 2 3 4 5 6 (cm 1) (f) FIG. 4. Transient changes in the IS vibrational spectra for CO(hollow)/Pd(111) (left) and CO(bridge)/Pd(111) (right) induced by a 450 nm laser pump pulse with 100 fs duration and peak intensity hitting the system at t=0.1 ps. The initial temperature is 100 K. (a) Induced electron temperature Te(t), (b) lattice temperature Tl(t), (c), (d) frequency shift, and (e), (f) linewidth change as a function of time. Purple, blue, green, orange and red lines correspond to absorbed fluences of 6, 19, 40, 80, and 130 J/m2, respectively. ration. More remarkable are the differences appearing at the largest fluences, which are also due to the different Te-dependence of the NC term between the two structures (see below). For F=80 J/m2, the blueshift persists for longer times in the bridge than in the hollow configuration. For F=130 J/m2the blueshift is only partially suppressed in the bridge configuration, while in the hollow configuration not only the blueshift is smaller but a redshift peak appears in the time interval in which Te>5700 K. The transient changes in the linewidth for the hollow and bridge configurations are shown in Figs. 4(e) and (f), respectively. Both the time evolution and the fluence dependence of the linewidth change are qualitatively the same for both structures. For the two lowest fluences, only the temperature dependence in EMPP contributes to the linewidth changes (not shown) and, as a result, ∆γ(t) follows the evolution of the lattice temperature. The contribution of the NC processes to ∆γ(t) starts to manifest at F= 40J/m2as an incipient peak emerging at t∼100 fs that becomes more pronounced as fluence increases. In particular, for the largest fluences of 80 and 130 J/m2, the well-defined peak structure that we observe for t<0.5 ps demonstrates that the linewidth changes in this time interval are dominated by the NC contribution. As a final remark, note that the quantitative differences between the two structures appear during the first 1-1.5 ps and are mainly due to the differences in the first order term, which is larger for the hollow structure. To give further insight into why the electronic and lattice temperatures determine the dynamics of the IS mode in the femtosecond (t<0.5 ps) and picosecond (t>1 ps) regimes, respectively, we explicitly study the effect that Teand Tlhave on the NC and EMPPC terms. Figure 5shows the contributions to ∆ωand the ∆γdue to the NC (dashed lines) and the EMPPC (solid lines) processes as a function of Tefor the hollow configuration. In the case of the EMPPC terms (∆ω[2] and ∆γ[2]) different lattice temperatures are considered. The ranges of temperatures shown in the figure, Te∈[100,8000] K and Tl∈[100,900] K, cover the temperatures induced by the laser pulses studied in this section. As shown in Fig. 5(a), the NC contribution ∆ω[1] produces a frequency blueshift for Te>900 K that starts to decrease its magnitude for Te>3000 K and gives place to a redshift for Te>6000 K. 8 5 0 5 10 15 20 (cm 1) (a) 01234567 Te (103 K) 0 2 4 6 8 10 12 14 16 18 (cm 1) (b) NC EMPPC - Tl =100K EMPPC - Tl =200K EMPPC - Tl =300K EMPPC - Tl =500K EMPPC - Tl =700K EMPPC - Tl =900K FIG. 5. The NC and EMPPC contribution to the change in frequency ∆ωand linewidth ∆γof the IS mode evaluated for the CO(hollow)/Pd(111) structure as a function of Tefor different Tl. Variations are given with respect to the thermal situation Tl=Te=100 K. Values of the EMPPC that are relevant for the regime t<(>)1 ps shaded in blue (red). This is the aforementioned non-monotonic behavior that causes the rapid quenching of the initial blueshift at F=80 and 130 J/m2in the hollow configuration. [See in Fig. 4(a) that at those fluences Tegoes well above 3000 K.] The dependence of ∆ω[1] on Teis qualitatively the same for the bridge configuration but the maximum blueshift occurs at 4000 K and it turns into a redshift at about 7000 K [25]. These qualitative differences explain that the blueshift is expected to be more pronounced and persistent in the bridge configuration [see Figs. 4(c) and (d)]. For fixed Tl, also the EMPPC term blueshifts as Te increases. This slight blueshift is however considerably smaller than the redshift of about −20 cm−1that occurs when Tlincreases from 100 to 900 K at fixed Te. As a result, the Tl-induced redshift dominates the EMPPC term for any realistic pair of temperatures (Te, Tl) describing the ultrafast regime, as explicitly shown for instance by the results of Fig. 3(b). Regarding the temperature dependence of the NC and EMPCC contributions to the linewidth change, Fig. 5(b) shows that the NC term ∆γ[1] is basically zero for Te< 1000 K. Next, it increases monotonically with Te, taking a value of ≈4 cm−1at 7000 K, which is the maximum Tereached at the highest fluence of 130 J/m2here considered [see Fig. 4(a)]. The almost constant dependence of ∆γ[2] on Te, which is expected for kBTe≪ω0λ[19], extends here up to ∼4000 K. Thus, for Te<4000 K, we still observe the nearly linear dependence on Tlthat is expected for kBTe≪ω0λ(see Sec. IV). Actually, we observe this linear dependence for the whole lattice temperature range shown in the figure, even if the dependence of the EMPPC term on Teis visibly enhanced for Te>4000 K. That enhancement is due mainly to the increasing importance that the IS dephasing process acquires at such high Te. In this respect, it is worth noticing that the peaks observed in ∆γat t=0.1 ps for the largest fluences of 80 and 130 J/m2[see Fig. 4(e)] are not only caused by the NC term, as it is the case for F=40 J/m2[see Fig. 3(c)], but also by the EMPPC term via the IS dephasing. VI. COMPARISON TO OTHER CO/METAL SYSTEMS The thermal [19] and non-equilibrium [28] IS vibrational spectra of CO adsorbed on the Cu(100) surface have been studied within the same theoretical framework that we use here. But before comparing the results of the IS vibrational spectra between the two systems, it is worth analyzing the intrinsic differences between the Cu and Pd surfaces that are likely to affect the adsorbate vibrational spectra. First, the laser-induced Te(Tl) for equal absorbed fluences is lower (higher) in Pd(111) than in Cu(100) because of the larger e-ph coupling constant in bulk Pd (λ= 0.4) than in bulk Cu (λ= 0.14). Second, the electron DOS around the Fermi level are very different for the two surfaces. Whereas in Cu(100) the sp-band crosses the Fermi level and the DOS above and below εF is rather constant, in Pd(111) the Fermi level lies at the edge of the d-band continuum and, as a result, the DOS is much larger below than above εF(see Fig. 8). As already discussed, the latter has important implications for the laser-induced IS dynamics due to the significant change of the chemical potential with temperature [25]. Under thermal conditions, the IS linewidth at T=10 K for CO/Cu(100) (γ≈2.9 cm−1[19]) is slightly smaller than for CO/Pd(111) (γ≈3.6 cm−1). Interestingly, the first-order contribution γ[1] is larger in CO/Pd(111) probably due to the larger DOS around εF. In contrast, the second-order contribution γ[2] is smaller due to a weaker electron-mediated coupling of the IS mode to the Pd(111) phonons. A possible explanation is the largest mass of the Pd atoms with respect to the Cu mass (notice that gµ,µ′ λ∼1/√Mλwith Mλthe mass associated with the phonon λ). Finally, the change of the IS linewidth with temperature in Cu(100) and Pd(111), 9 namely, ∆γ(T)≈0.9 cm−1and 0.8 cm−1, respectively, for T∈[10,150] K, are very similar. This temperature dependence is dominated in both systems by the coupling with the surface modes. We note once again that the measured linewidth increase of about 3 cm−1between 50 K and 300 K for CO/Pt(111) [20] is close to our results for CO/Pd(111). Regarding the properties of the transient vibrational spectra induced by femtosecond laser pulses, there are significant differences between CO/Cu(100) and our results for CO/Pd(111). One important difference, which is directly related to the small Tlthat the laser-excited electrons induce in Cu(100) as compared to Pd(111), refers to the time evolution of the transient spectra. The largest frequency and linewidth changes in CO/Cu(100) occur during the initial Te-driven subpicosecond regime in which Te≫Tl. In the subsequent Tl-driven picosecond regime (i.e., when Te≃Tl), the temperatures are typically much smaller (e.g., about 400 K for F=170 J/m2) and so the frequency and linewidth changes. In contrast, the spectral changes in CO/Pd(111) are more pronounced during the steady picosecond regime that in this system is characterized by large lattice temperatures (e.g., ∼1000 K for F=130 J/m2). Only for the largest fluences (with Te>4000 K), the transient changes in the fast and steady regimes are comparable. Furthermore, there are additional qualitative differences between the two systems that it is worth remarking. In CO/Cu(100) [28], the pronounced IS frequency redshift that occurs during the first picosecond is ruled by the first-order NC processes. The transient linewidth change however is dominated by the EMPPC processes. In particular, the largest change is due to the IS dephasing that takes place within the first picosecond. In CO/Pd(111), an initial frequency blueshift driven by the NC processes is expected to occur during the first 0.3 ps for absorbed fluences large enough to raise Teabove 1000 K. The blueshift progressively turns into a large and steady redshift, which is caused by the EMPPC processes leading the late picosecond regime. The latter dominates the linewidth change. It is only for the largest fluences (Te>4000 K) that the subpicosecond regime controlled by the NC processes is also visible. It is also meaningful to compare our results to timeresolved vibrational spectra measured in other surfaces with pump-probe vibrational spectroscopy. In Ref. [12], the experimental linewidth change of the IS mode of CO on Pt(111) exhibits the same behavior that we obtain here for CO/Pd(111). Specifically, the Tl-driven late steady regime appears at all the experimental fluences (40, 80, and 130 J/m2), while the fast Te-driven subpicosecond regime becomes increasingly important at F≳80 J/m2. Regarding the frequency change, the late steady frequency redshift and its dependence on Fare also quite similar to our CO/Pd(111) results. However, a redshift instead of a blueshift seems to occur in the subpicosecond regime. Nonetheless, there is an unusual nonmonotonic dependence of ∆ωon the absorbed fluence as found in later experiments with improved subpicosecond resolution [13]. We speculate that the differences between CO/Pt(111) and CO/Pd(111) may be caused by the differences in their electron DOS. Even if qualitatively the DOS of the two systems are similar, the decay of the d-band edge above εFis sharper in Pd(111) than in Pt(111). Thus, we may expect a somehow softer dependence of the chemical potential on temperature that would affect the frequency sift behavior. The fact that the transient IS mode spectra follows the dynamics of the lattice temperature in the picosecond regime is common to other experiments in transition metal surfaces. Lane et al. experimentally studied the ultrafast vibrational dynamics of CO and NO on the Ir(111) surface [8,9]. The authors found that the IS frequency shift of NO (CO) follows the time evolution of the lattice temperature when the absorbed fluence is below 8 J/m2(20 J/m2). Similarly, in the case of CO/Ru(1000) it has been shown that the IS dynamics follows the lattice temperature at least up to a fluence of 55 J/m2[3]. Finally, it is worth remarking the blueshift observed for CO molecules coordinated to ruthenium tetraphenylporphyrin adsorbed on Cu(110) with pump-probe vibrational spectroscopy [16]. The values of the blueshift at different fluences are very similar to the results presented here for CO/Pd(111). For sufficiently small fluences, only a redshift that evolves in time with the lattice temperature is observed. However, for larger values of the fluence, a blueshift in the subpicosecond regime emerges. VII. CONCLUSIONS In summary, we have studied the vibrational relaxation of the internal stretch mode of CO adsorbed on the Pd(111) surface using density functional theory in combination with many-body perturbation theory. In particular, by evaluating the phonon self-energy up to second order in the electron-phonon (e-ph) coupling matrix elements, we have characterized the two mechanisms that participate in the vibrational relaxation of the CO internal stretch mode: the first-order interband nonadiabatic coupling (NC) and the second-order intraband electron mediated phonon-phonon (EMPPC) coupling. Under thermal conditions both mechanisms contribute significantly to the linewidth γ(T=0)≈3.6 cm−1, with the interband contribution γ[1](T=0)≈2.7 cm−1being the largest. The redshift and change in the linewidth that are induced as temperature increases are explained mainly in terms of the coupling of the IS mode to the low energy phonon modes of the system. We have also studied the subpicosecond and picosecond transient changes that are induced in the IS phonon frequency and linewidth when the CO/Pd(111) system is excited with a femtosecond laser pulse (non-equilibrium conditions). For large electronic temperatures (with Te≫Tl), the particular electronic structure of Pd may give rise to an effective screening of the e-ph interaction, inducing an uncon- 16 4116 (2006). [77] D. Karlsson, R. van Leeuwen, Y. Pavlyukh, E. Perfetto, and G. Stefanucci, Fast green’s function method for ultrafast electron-boson dynamics, Phys. Rev. Lett. 127, 036402 (2021). [78] V. P. Zhukov, F. Aryasetiawan, E. V. Chulkov, and P. M. Echenique, Lifetimes of quasiparticle excitations in 4d transition metals: Scattering theory and LMTO-RPAGW approaches, Phys. Rev. B 65, 115116 (2002). [79] N. W. Ashcroft and N. D. Mermin, Solid State Physics (Holt-Saunders, 1988). [80] Y. V. Petrov, N. A. Inogamov, and K. P. Migdal, Thermal conductivity and the electron-ion heat transfer coefficient in condensed media with a strongly excited electron subsystem, JETP Letters 97, 20 (2013). [81] P. Szymanski, A. Harris, and N. Camillone, Temperature-dependent electron-mediated coupling in subpicosecond photoinduced desorption, Surf. Sci. 601, 3335 (2007). [82] A. P. Miiller and B. N. Brockhouse, Crystal dynamics and electronic specific heats of palladium and copper, Can. J. Phys. 49, 704 (1971). [83] C. Kittel, Introduction to solid state physics, 6th ed. (Wiley & Sons, New York, 1986). [84] P. Johnson and R. Christy, Optical constants of transition metals: Ti, V, Cr, Mn, Fe, Co, Ni, and Pd, Phys. Rev. B 9, 5056 (1974). [85] S.-Y. Hong, P. Xu, N. R. Camillone, M. G. White, and N. Camillone, Adlayer structure dependent ultrafast desorption dynamics in carbon monoxide adsorbed on pd (111), J. Chem. Phys. 145, 014704 (2016).