Energy dissipation at metal surfaces

Abstract Conversion of energy at the gas–solid interface lies at the heart of many industrial applications such as heterogeneous catalysis. Dissipation of parts of this energy into the substrate bulk drives the thermalization of surface species, but also constitutes a potentially unwanted loss channel. At present, little is known about the underlying microscopic dissipation mechanisms and their (relative) efficiency. At metal surfaces, prominent such mechanisms are the generation of substrate phonons and the electronically non-adiabatic excitation of electron–hole pairs. In recent years, dedicated surface science experiments at defined single-crystal surfaces and predictive-quality first-principles simulations have increasingly been used to analyze these dissipation mechanisms in prototypical surface dynamical processes such as gas-phase scattering and adsorption, diffusion, vibration, and surface reactions. In this topical review we provide an overview of modeling approaches to incorporate dissipation into corresponding dynamical simulations starting from coarse-grained effective theories to increasingly sophisticated methods. We illustrate these at the level of individual elementary processes through applications found in the literature, while specifically highlighting the persisting difficulty of gauging their performance based on experimentally accessible observables. Graphical Abstract


Introduction: concepts & open questions
In chemical reactions at solid surfaces, different forms of energy are converted into one another. Chemical energy is released or consumed by breaking and making individual chemical bonds of adsorbed species or intermediates formed in the course of the reaction. Along the way, parts of this chemical energy are converted (transiently) into translational or vibrational energy of these surface species. Energy exchange can also occur with substrate degrees of freedom (DOFs), which (ultimately) will lead to a thermalization of adsorbates prevailing for a sufficiently long time at the surface. Net energy flow out of the adsorbate/surface fringe is thereby denoted as dissipation. Apart from (vibrationally excited) desorbing species, energy dissipation proceeds prominently into the substrate bulk. Largely, this occurs through the excitation of lattice vibrations, so-called substrate phonons [1,2]. At metal surfaces, however, the non-adiabatic excitation of electron-hole (eh) pairs provides another competing energy dissipation channel. The continuous distribution of electronic states across the Fermi level allows in principle to even excite lowest-energy eh-pairs -an argument that has often been used to suspect a high relevance of this channel at metal surfaces [3][4][5][6][7].
However, the real role and relative importance of both dissipation channels and how this varies across systems are presently largely unclear [3,8]. In fact, on this microscopic level we still know very little about these mechanisms of energy dissipation in general. This is rather intriguing, considering that chemical reactions at metal surfaces drive important applications and technologies like heterogeneous catalysis or surface growth. In the established microkinetic theories used in these fields [9][10][11] one, for instance, generally assumes that any reaction heat released in an exothermic reaction step is quasi-instantaneously dissipated away to ensure an immediate thermalization of the surface species. This motivates a Markovian view of an overall surface chemical reaction as a sequence of independent elementary steps such as adsorption, diffusion, reaction, and desorption as depicted in Figure 1. Even though unanimously applied today [12], the general validity of this Markovian picture is less clear from the microscopic perspective of the phononic and eh-pair dissipation mechanisms. In a typical exothermic surface reaction like dissociative oxygen adsorption at a transition metal surface heats of the order of a few electron volts (eV) are released. Considering that, e.g. the energy scale of phonons is meV, one may start to wonder how instantaneous the phononic energy uptake and concomitant thermalization of surface species really is.
Additionally fueled by sustainability considerations concerning waste heat recovery or general heat management in heterogeneous catalysis, this has motivated fundamental research to arrive at a better microscopic understanding of energy dissipation at metal surfaces. Experimentally, the key approach are dedicated surface science experiments that study individual elementary processes An impinging molecule may, for instance, be directly reflected from the metal substrate, nevertheless exchanging energy and momentum with the surface. It may, however, as well dissociate on the surface -potentially through a vibrational precursor state -with the excess kinetic energy leading to hyperthermal motion of the fragments, so-called 'hot diffusion'. In all of these processes, energy exchange with the metal occurring either through the excitation of lattice vibrations (phonons) or non-adiabatic electron-hole (eh)-pairs may significantly influence the resulting dynamics.
at well-defined single-crystal surfaces in ultra-high vacuum [13,14]. The question of energy dissipation is thereby a sub-topic in the broader context of longresearched gas-surface dynamics (GSD). A central problem here is that key observables accessible in GSD experiments are often ensemble averages (vide infra). At present it is largely unclear in how much these observables are actually sensitive to details of the microscopic energy dissipation mechanisms. In fact, an unambiguous interpretation of the measured data is typically impossible without detailed modeling efforts. The dramatic increase in computer and algorithmic power has in this respect led to a strong surge of corresponding, in particular first-principles based theoretical work in recent years. This work has its own limitations though. It centrally still struggles with the necessity to simultaneously provide a reliable account of the surface electronic structure (and ensuing energetics) and adequately follow or sample the on-going surface dynamics. At the electronic structure level, the difficult task is to procure a description of localized orbitals of surface species on the one hand and the highly delocalized metallic band structure on the other [15,16]. Considering dynamical simulations, issues arise in turn from the extensive ensemble averaging required to compute key experimental observables, and from excessive system sizes when aiming to explicitly resolve the energy transfer to the multitude of electronic and phononic DOFs of an extended metal (Me) substrate.
These challenges can presently only be met by numerically efficient, effective theories. They typically rely on density-functional theory (DFT) with semi-local exchange-correlation functionals [16][17][18][19][20] and still often have to resort to bathtype treatments of the substrate DOFs. Open questions these theories try to address at the level of an individual elementary process (cf. Figure 1) include the relevance of each dissipation channel under specific conditions, whether these channels interact and influence each other or simply provide additive contributions, and how this picture will change from one system to another. Since all theories potentially able to consistently address these questions in practice are by nature approximate, validation by comparison to experiment is obviously vital. As such, another important aspect is also how experimentally accessible quantities in different types of GSD experiments at single crystal surfaces can be calculated. This not only from the perspective of reproducing these quantities accurately, but also to understand to which extent the details of the underlying energy dissipation mechanisms actually matter for the target properties of interest. In this topical review we survey the state-of-the-art of corresponding theories and their application in practice. For this we focus in subsequent sections on work done to elucidate the role of energy dissipation in adsorption processes (Section 2), in vibrational motion (Section 3) and in surface diffusion (Section 4). Emphasis is placed on concepts and in how much the impressive amount of work carried out in particular over the last years allows already to derive some general insights and ruling principles. An extensive list of references guides the interested reader to the methodological and technical details. The scope of this work is instead to introduce this lively field, its accomplishments and challenges to a broader audience.

Inelastic scattering & adsorption
Energy dissipation in adsorption processes is suitably studied in molecular beam experiments that expose an initially clean single-crystal surface to a beam of gasphase species of defined initial kinetic energy and impinging from a defined initial angle. Detailed information about energy dissipation can then, e.g. be derived from an analysis of the translational and quantum state resolved vibrational and rotational energy distributions of inelastically scattered molecules [14,21]. For the ensemble of gas-phase molecules in a molecular beam, another central experimentally accessible kinetic quantity is the (initial) sticking coefficient, i.e. the fraction of molecules that has lost sufficient energy to remain adsorbed at the surface [22]. The latter is usually a function of several variables such as beam incidence energy and angle, substrate temperature and specifically prepared rotational and vibrational quantum state of the impinging molecules [23]. Sticking coefficients are thus beam-ensemble averages over the binary measure of whether an individual molecule adsorbs on a surface or not. Sogained information thus naturally convolutes effects of energy dissipation from several successive elementary processes that occur, e.g. during a dissociative adsorption process, such as the actual dissociation event, vibrational motion in a potential trapped precursor state, and subsequent diffusion of the reaction products (cf. Figure 1). As such, sticking coefficients are a prime example of GSD observables for which it is clear from the start that they will at best provide only very indirect insight into the detailed energy dissipation mechanisms. A substantial amount of first-principles modeling has nevertheless been devoted to compute these quantities. On the one hand, this is because high-accuracy references from experiment are available for a wide range of systems. On the other hand, sticking coefficients are fundamental kinetic quantities describing adsorption processes in microkinetic models. As such, it is most intriguing and relevant to understand which aspects of microscopic energy dissipation propagate through to this more coarse-grained, technological level.

GSD within the frozen surface approximation
Addressing sticking coefficients from a modeling perspective poses a significant challenge to theory. While impinging on the surface, the adsorbing species sample a wide range of configurations far away from the (ultimate) equilibrated adsorption geometry. This requires to accurately describe the potential energy surface (PES) representing the adsorbate-surface interaction over corresponding wide regions. Even for simplest diatomic adsorbates just accounting for the molecular DOFs (the position of the molecule's center of the mass above the surface, its bond length and angular orientation) leads already to a six-dimensional PES that needs to be computed and represented. Explicitly treating at least all these molecular DOFs has thereby been established as a complete necessity by numerous studies which have specifically outlined the dangers of a reduced dimensionality treatment and shown how (intuition-based) simplifications over the adsorbate coordinates can yield dramatically wrong results with respect to extracted dynamical properties [2,[24][25][26][27][28][29][30][31]. A proper description of the beamensemble requires on top of this extensive statistical sampling of typically several tens of thousands of impinging molecular trajectories starting from varying initial conditions. A key concern thus foremost lies in procuring a numerically efficient, yet sufficiently accurate description of the aforementioned adsorbatesurface interaction potential.
The twofold nature of this challenge has promoted a likewise two-step, 'divideand-conquer'-like strategy [25,[32][33][34]: First, a continuous PES function is constructed by interpolating or fitting to a large set of ab initio data points that sample the adsorbate DOFs over a rigid, or so-called 'frozen' surface (FS) (depicted in blue in Figure 2). The evaluation of energies and forces on this continuous PES function then comes at a practically negligible computational cost which allows for extensive molecular dynamics (MD) simulations in a second step. While computationally tractable, at this FS level the inherent neglect of substrate mobility obviously precludes any (phononic) energy exchange with the lattice. Focusing for the time being also exclusively on the adiabatic Born-Oppenheimer PES, there is thus no energy dissipation mechanism explicitly considered in the model. Sticking must correspondingly be concluded by some ad hoc criterion, typically assuming a molecule to stick after its trajectory has prevailed for a sufficiently long time close to the surface or upon exceeding a threshold interatomic separation for the case of dissociative chemisorption.
Generally one would expect this crude approximation to work the better, the weaker the adsorbate couples to the lattice vibrations. Due to the large mass mismatch, this is best fulfilled for lightest H 2 molecular adsorbates at transition metal surfaces, which for many reasons has been a prototypical GSD system [35] that has received significant attention over the years anyway [24][25][26][27][28][35][36][37]. It is thus not altogether surprising that already FS simulations based on DFT-derived, six-dimensional PESs have been found to largely reproduce experimental H 2 sticking curves [8,[37][38][39]. An escalating activity is nevertheless being directed toward heavier and more strongly chemisorbing adsorbates in recent years. For these systems, in particular energy exchange with the lattice is expected to play an increasingly important role and requires that some (at least approximate) account of substrate mobility is included in the theory [40,41].

Effective models for phononic energy dissipation
Opening the phononic dissipation channel generally implies an exploding dimensionality for the problem at hand. Explicitly considering the motion of substrate atoms was hitherto largely prohibitive within the divide-and-conquer ansatz due to the increasing difficulty that is associated with each added DOF when constructing a continuous PES representation [34]. This situation has just started to change as advanced, high-dimensional interpolation schemes such as the Behler-Parinello atomistic neural network approach [42][43][44] become increasingly applicable in the context of GSD [45,46]. As already conceivable from the simple Baule limit [2,14], however, energy transferred to the lattice in each direct adsorbate-surface collision is still small for light diatomic adsorbates. On this level, one can then tentatively consider the interaction with surface phonons to be more of a small perturbation to the FS-dynamics such that an explicit full-dimensional description of the lattice motion may not be altogether necessary.
This motivates energy sink models that incorporate an effective treatment of surface mobility [47][48][49][50][51][52][53]. Such approaches generally target to reduce the phononic fine structure of the substrate into computationally convenient augmentations of the static-surface interaction in order to include some (approximate) account of energy exchange with the lattice. Motivated by an Einstein-like picture for the phononic system, the surface oscillator (SO) model [47,48] mimics a rigidly moving substrate by a single harmonic oscillator of assigned frequency and mass parameters. As shown in Figure 2, coupling to the adsorbate is then straightforwardly described through a space-rigid shift in the FS expression for Hierarchical approaches to enrich simulations on a FS PES with an account of lattice motion. The SO model first introduces a 3D harmonic oscillator corresponding to a rigid shift of the entire PES in all spatial directions, thus allowing for adsorbate-surface energy transfer following a simple collision model. Building on this, the GLO adds a further ghost oscillator linearly coupled to the SO, which in turn is further coupled to a heat bath within an effective generalized Langevin description. This allows to also include energy dissipation from the SO to the bulk. the PES. Extending on the latter as further depicted in Figure 2, the generalized Langevin oscillator (GLO) additionally incorporates the approximate effect of a bulk thermal bath by coupling the SO reaction zone to an additional so-called 'ghost' oscillator [47,[49][50][51][52][53]. Energy dissipation and thermal fluctuations are consequently accounted for by subjecting the latter to frictional and random forces, respectively, rigorously satisfying the fluctuation-dissipation theorem in the numerical solution of the resulting generalized Langevin equation. While thus incorporating only a very limited number of additional DOFs, these models still allow to lift the FS-approximation by introducing the concept of a surface temperature with concomitant thermal motion [1,48] -even if it is just described in a coarse-grained way.
One of the main practical advantages of the SO/GLO approaches is their convenient application to any pre-existing adsorbate-surface PES that has already been derived on the level of the FS approximation. This also includes, of course, the continuous PES representations obtained from DFT within the aforementioned divide-and-conquer framework. In a first such application to H 2 scattering from Pd(1 1 1) [54], Busnengo and coworkers successfully used the SO model to reveal a channel of dynamic surface trapping that could explain the otherwise puzzling temperature dependence which was observed experimentally in the rotational excitation of scattered molecules [55]. A similar result was reached only shortly afterwards also by Wang and coworkers for D 2 scattering off Cu(1 1 1), but on the basis of quantum-mechanical (QM) wave-packet calculations [56]. Opening also the GLO bulk dissipation channel as an outlet for the accumulating heat could further predict the low-temperature stabilization of a molecularly chemisorbed H 2 species at Pd(1 1 0) [57]. While thus revealing important mechanistic details of the adsorption process, the inclusion of SO/GLO dissipation nevertheless provided only marginal corrections to averaged sticking probabilities for the light H 2 adsorbate [58] (consistent with the previously mentioned good agreement with experiment already at the FS level). A much more dramatic effect was instead demonstrated for the heavier O 2 molecule, during its highly exothermic adsorption at Pd(1 0 0) [59]. Here, even a qualitative change in the shape of the calculated sticking curve was reported in comparison to the corresponding FS results, yielding thus a considerably improved agreement to experimental data [60].
Despite these and many more successes (cf. e.g. Refs. [61][62][63][64][65]), it is nevertheless always prudent to bear in mind the inherent limitations that come with the underlying SO/GLO approximations. Recognizing foremost situations where these are bound to fail are cases where energy uptake cannot be reliably described on the level of only a single phonon excitation. Based mainly on energetic considerations, the frequency of the single SO oscillation is often taken to correspond to the surface (low-frequency, acoustic) Rayleigh modes, as these are generally assumed to be predominantly excited upon first impact with the surface (cf. e.g. Refs. [66][67][68][69] and references therein). The validity extent of this assumption, however, remains to its larger part unknown and could represent a highly dangerous over-simplification [70,71], in particular when going beyond prototypical studies of lightest adsorbates [72][73][74][75]. Even within a dominant one-phonon picture though, the SO model can only go as far as to provide a mechanical coupling of adsorbate-surface momentum exchange, i.e. accounting thus for energy loss to the substrate through lattice recoil and a back transfer to the adsorbate with surface temperature. Thereby considered 'stiff' lattice shifts will thus inherently fail to capture concomitant changes in the FS-PES that are induced by thermal displacements of the metal atoms from their equilibrium positions. This additional effect of vibrational surface motion can play an important role in the adsorption dynamics, as demonstrated prominently for (direct) CH 4 dissociation at several metal surfaces [76]. Here (as also in other reported instances [77,78]), DFT calculations specifically showed significant lowering of the activation barrier as a metal atom is puckered out of the surface at the transition state [79,80]. Effectively accounting for both of the aforementioned effects through two (independent) DFT-derived adsorbatesurface coupling parameters [81], the quantum dynamics simulations of Jackson and coworkers thus arrived at a semi-quantitative agreement to experimental sticking curves (also reproducing the observed mode specificity and bond selectivity) that could not otherwise be achieved on the SO level alone [82]. In conclusion, we note that this so-called barrier height modulation arising from thermal fluctuations of the lattice is expected to depend strongly on both the specific system and surface temperature under investigation, so that its overall relevance remains altogether a priori unclear in the context of other gas-surface reactions, as will be further discussed in the following.

Explicitly resolving surface motion
Direct ab initio molecular dynamics (AIMD) simulations stand fundamentally free of the aforementioned limitations discussed in the context of effective energy sink models. In this method the forces needed to propagate the equations of motion are full-dimensional Hellmann-Feynman forces computed from DFT 'on-the-fly' at each MD time step, so that no fitting or interpolation of the PES is required. Metal atoms are thereby treated explicitly within a fully mobile substrate and no a priori assumptions have to be made with respect to the heat bath Hamiltonian. This further naturally circumvents the extensive PES preevaluation required within the divide-and-conquer ansatz, but in turn comes at a tremendous increase in computational cost for each individual trajectory. With the advent of supercomputers and continuous developments toward improved algorithms for electronic structure calculations, this challenge is being steadily overcome and it is now increasingly possible to compute a meaningful number of AIMD trajectories for many systems that could not even be addressed in static calculations a few years ago [83].
The combined dynamical picture of both adsorbate-surface coupling effects discussed previously (i.e. mechanical coupling of momentum exchange and barrier height modulation due to thermal fluctuations), should in principle be completely accounted for within AIMD simulations. Indeed, for the aforementioned CH 4 dissociation (where both effects are known to be important [76]) Kroes and coworkers predicted satisfying, semi-quantitative agreement to experimental reaction probabilities [84]. In the meantime, similar results have also been reached for H 2 [37,85,86], N 2 [87,88], O 2 [89], and CO 2 [78] adsorbates, overall showing AIMD simulations to provide a reasonable account of surface temperature effects. With reaction probabilities as the target observable, the effective advantage of the detailed AIMD account of adsorbate-surface energy transfer-over corresponding, numerically more attractive SO/GLO modelsultimately depends on the importance of thermal fluctuations (assuming of course the validity of the single-phonon approximation). As such, a simple GLO model has for example been shown to perform on an equal level as AIMD regarding calculated sticking probabilities for the dissociative N 2 adsorption on W(1 1 0) [90], while being entirely inappropriate for other systems such as, e.g. the aforementioned CH 4 on various transition metals [76]. Unfortunately, however, such performance is extremely difficult to predict a priori, especially in view of the fact that both lattice coupling effects will induce qualitatively similar changes to the overall experimental sticking function with surface temperature.
At this point, however, it must be emphasized that the microscopic description of phononic dissipation as provided by present-day AIMD simulations is still far from perfect -a fact which stems (foremost) from practical, rather than conceptual, considerations. As the involved computational cost is largely determined by the substrate size, state-of-the-art AIMD studies are limited in practice to supercell setups involving slab models consisting of typically only a few surface lattice constants and layers [83]. By construction, such models fail to provide the desired (high) resolution when representing the underlying phononic fine structure. For example the description of surface modes -which are assumed to be crucial for the initial energy uptake at the interface [74] -will unavoidably suffer from the limited slab thickness that must be adopted in computationally tractable DFT supercells. Most importantly, the imposed periodic boundaries (which are essential to describing the underlying metallic band structure) will restrict any phonon propagation within the finite extent of the latter, thereby leading to an unphysical confinement of the released chemical energy. Largely exothermic reactions in particular are thus increasingly prone to significant overheating of the metal substrate -an effect which severely compromises the description of the ensuing equilibration process and may even critically modify the actual adsorbate dynamics, as will be discussed in more detail in Section 4.2. While the effect of these shortcomings has not been investigated in detail, they may still be responsible for remaining experiment-theory discrepancies found in literature. A prototypical such example is the dissociation of HCl on Au(1 1 1) for which energy dissipation is believed to play a key role [91]. While (independently performed) AIMD simulations [86,92] have consistently overestimated the experimental probabilities [91], the corresponding reaction dynamics of this system remain an enigma and have become quite a controversy in the field [86].
An explicit treatment of substrate mobility within a (suitably) sizable heat bath thus continues to pose a big challenge to contemporary AIMD simulations and thereby derived accounts of phononic dissipation. This situation is aggravated by the extensive statistical sampling that is hereby required not only for describing the beam-ensemble, but also the thermal distribution of metal atoms at finite surface temperatures. With current computational resources typically providing up to 'only' a few hundred trajectories, AIMD-derived statistics are thus presently still very poor in particular when addressing rare events. In this respect, using existing AIMD configuration sampling as the basis for parameterizing the aforementioned upcoming high-dimensional neural network-based interpolation schemes [42][43][44] represents a promising alternative for the future that will allow for routinely evaluating several orders of magnitude more trajectories as compared to further explicit AIMD simulations [45,46]. In the meantime, while entirely neglected so far, the channel of non-adiabatic energy dissipation adds yet another facet of complexity to GSD and toward which attention is directly turned in the following.

Exciting electron-hole pairs: non-adiabatic effects
Energy dissipation arising from non-adiabatic effects, i.e. the dynamical coupling of the adsorbate nuclear motion to excitations in the metallic continuum of electronic states at the Fermi level, still holds a rather suspect role in gas-surface scattering/adsorption processes [3][4][5][6][7]. One of the most controversial showcases in this respect is the highly exothermic oxygen dissociation at Al(1 1 1) that proceeds through an O 2 triplet-singlet spin transition as the molecule approaches the surface [8,74]. Here, extremely low initial sticking coefficients measured for small O 2 incidence energies are strongly indicative of activated adsorption -a process which goes entirely amiss within the essentially barrierless (adiabatic) PES that is evaluated on the level of the DFT generalized gradient approximation (GGA) [93]. It remains, however, altogether unclear whether this actually indicates a hindered triplet-singlet transition due to spin selection rules that would give rise to corresponding dynamical barriers. Alternatively, this omission could instead be rooted in shortcomings of the employed semilocal exchange-correlation functionals [74,94,95]. Another prominent example vouching for pronounced non-adiabatic effects is the energy loss associated with vibrationally excited NO scattering from Au(1 1 1) [7,[96][97][98]. The non-adiabatic relevance is strongly suggested here by multi-quantum relaxation of the vibrating molecule that is clearly observable on the metal, but completely absent on an insulating LiF surface [97].
Going beyond the Born-Oppenheimer approximation in order to account for such non-adiabatic effects requires in principle the full propagation of a combined, high-dimensional nuclear-electron wave function that would naturally also include the multitude of potential electronic excitations. Such an endeavor, however, is computationally out of reach in the context of GSD for at least decades still to come. Including electronically non-adiabatic effects in computational simulations thus inevitably requires further approximations and simplifications leading to an effective treatment, such as applying a mixed quantum-classical description of the system in terms of a mean field Ehrenfest [99][100][101][102] or surface hopping framework [99,[103][104][105][106]. The latter has, for instance, been successfully invoked to explain the aforementioned NO vibrational de-excitation on Au(1 1 1) through an electron-transfer mechanism [105,107,108]. Further coarse-graining the electron dynamics can be achieved on the basis of time-dependent perturbation theory [109][110][111]. Counting among the disadvantages of all these methods is that they either remain computationally very intense, or require extensive parametrization, or impose specific symmetry constraints on the simulated trajectories such that a general applicability to routinely performed gas-surface calculations is lost.
Presently the only viable solution for large-scale gas-surface simulations thus lies in entirely replacing the electronic degrees of freedom within the concept of electronic friction [112][113][114]. As formally shown by Tully and Head-Gordon [114] -and recently more rigorously also by Dou and Subotnik [115][116][117] -the idea behind this approach is to start with a mixed quantum-classical description of the dynamical (electron-nuclear) system. Following a suitable Miller-Meyer action-angle transformation, the explicit electronic DOFs are subjected to a generalized Langevin replacement that is valid within the weak-coupling approx-imation and results in a generalized Langevin equation for the nuclear DFOs. This, however, is not the form of electronic friction theory that is used in practical MD simulations today. Hence, in a second step the explicit time-dependence of the friction kernel (collectively condensing the dynamics of the electronic DFOs) is removed by invoking a Markov approximation in the electronic subspace implying (infinitely) short electronic coherence times [111]. Thus assuming memoryless eh-pairs that immediately 'forget' about the past yields ultimately a Langevin equation for the nuclear dynamics that combines all non-adiabatic effects in a single electronic friction coefficient. Altogether, the electronic friction approach thus allows to augment the Born-Oppenheimer ground-state dynamics with a dissipative friction force and a temperature-dependent fluctuating force, as manifested through the fluctuation-dissipation theorem. This combination of conceptual simplicity and numerical efficiency, along with the possibility of application to dynamics on a potentially pre-calculated Born-Oppenheimer PES has overall served to the great popularity of this method.
While the idea of including frictional energy losses into nuclear dynamics dates back to more than 30 years ago [112,[118][119][120][121][122] it was not until much more recently that Juaristi and coworkers adopted this concept for routine DFT-based applications within the divide-and-conquer framework [123]. The key ingredient hereupon relied on is an approximate construction of the required molecular friction tensor from independent contributions of individual atoms. As depicted in Figure 3, these are estimated from a simple jellium-based embedding model [118,119,[124][125][126] to finally yield isotropic atomic friction coefficients as a function of an embedding electronic density. Within this so-called local density friction approximation (LDFA) [113,123], the latter is taken as the local electron density of the clean metal surface at the position of the (individual) adsorbate atoms. This simplification consequently allows for a straightforward interpolation of the friction coefficients into a convenient analytical function of the electronic density only. This distinguishes the LDFA from more refined orbital-dependent friction (ODF) methodology [114,122,127,128] whose numerical involvement [128][129][130] (until only very recently [131]) has hindered application to highdimensional models though [132][133][134][135][136][137].
The inherent assumption of independent atoms within the original LDFA approach, however, triggered a controversial discussion [123,138,139], mostly because it ignores the adsorbate's molecular character and thus lacks, for instance, the steep increase of friction coefficients at dissociative transition states predicted from ODF [134,135]. In this regard, it has recently been suggested to introduce molecular information through suitably Hirshfeld-partitioned [140] full system densities to map to the LDFA model system [141] rather than relying on clean metal densities (vide infra). Moreover, one may pragmatically argue that the immanently low velocities in such transition state regions effectively suppress the contribution of the velocity-scaled friction term within the actual Figure 3. Construction of electronic friction coefficients within the LDFA [123]. The interacting molecule-surface system is first approximated through independent atoms being embedded in the electron density of the clean metal surface. The local electronic density at the atomic positions is then used to independently map to an isotropic atomic embedding model system of a spherically symmetric impurity in jellium. Finally, electronic friction coefficients are evaluated from the scattering phase shifts of the Kohn-Sham states at the Fermi-momentum for this model [118,119,[124][125][126]. This ultimately yields, for each element, an electronic friction coefficient as a function of the embedding density and can thus be conveniently evaluated and tabulated prior to dynamical simulations. It has later been suggested to introduce molecular information by constructing the embedding density via a suitable Hirshfeld-partitioning [140] of the full system electronic density [141].
dynamics and thus wash out potential inaccuracies [123] -at least when focusing on reaction probabilities [131].
Retaining the LDFA within the independent atom approximation and its unique numerical efficiency, many divide-and-conquer studies have been revisited over the past years to augment the existing FS-based MD simulations with electronic friction. Alducin and coworkers specifically investigated the dissociative adsorption of several small molecules such as H 2 [123], N 2 [123,142], H 2 O [143] and CH 4 [144] on various metal surfaces. The overall conclusion drawn from these studies is that the non-adiabatic energy dissipation channel produces only small effects on calculated sticking probabilities. Alone a full dimensional account of the adsorbate PES was in fact found to alleviate certain discrepancies with respect to experiment that had previously been erroneously assigned to nonadiabatic effects through a low-dimensional ODF-description. Prominent such examples are N 2 adsorption on Ru(0 0 0 1) [134] and the vibrational de-excitation of H 2 scattering from Cu(1 1 1) [37,38,135,145].
Similar to what has already been discussed in the context of phononic dissipation, however, the overall insensitivity reported for averaged sticking coefficients does not necessarily extend to mechanistic details of the underlying adsorption/scattering process. Conclusions along these lines have, e.g. already been independently reached for the dissociative H 2 adsorption on Ru(0 0 0 1) [146] and Ag(1 1 1) [131]. Here, the authors specifically noted minor nonadiabatic energy losses that may be unimportant to reaction probabilities, yet significantly influence other experimentally accessible observables such as the energy distribution of backscattered molecules. This is of particular interest considering the tensorial character of the electronic friction that is -in contrast to LDFA -directly accessible from ODF [114,134]. In this regard, recent studies by Maurer and coworkers [127,128,131] actually demonstrated that hitherto often neglected contributions are indeed rather pronounced and may, besides the wellknown mode-specific electronic friction [127,128,[147][148][149][150], lead to dynamical steering and vibrational mode coupling effects. As a concluding remark, it must be emphasized though that the popular electronic friction approach may even be altogether inappropriate for generally gauging the relevance of this dissipative channel. This is due to its inability to describe strong non-adiabatic coupling, as has, e.g. been most prominently shown for the multi-quantum vibrational transitions occurring during NO scattering on Au(1 1 1) [151]. Notwithstanding, at present electronic friction is the best shot we have at all to address nonadiabatic energy dissipation in high-dimensional surface simulations.

Combining phonons and electron-hole pair excitations
Combining the two dissipation channels in dynamical simulations promises illuminating insight into the relevant importance of each channel, at least on the specific level of the theory employed. Unfortunately, such applications have hitherto been rather scarce. In the late 1990s Tully and coworkers were the first to combine the electronic friction approach with an explicit description of lattice degrees of freedom for the scattering of CO from Cu(1 0 0) [152]. In this early study, however, the authors relied on a semi-empirical potential for describing the adsorbate-metal interaction and no supporting experimental data were available at the time. A dominating role of the phononic dissipation channel was reported that outperformed the non-adiabatic energy losses by up to a factor of ten (depending on the CO incidence energies). Similar findings were also reported much more recently for the scattering of nitrogen atoms and molecules using a combined electronic friction/GLO model applied to DFTderived continuous PES models [61,62,153,154]. Here, phononic effects were found superior over the non-adiabatic counterpart on essentially every measure that was considered. This said, most qualitative aspects of the dynamics could admittedly already be reached on the level of a purely adiabatic FS model when accounting for all adsorbate DOFs.
Quite in contrast, for somewhat lighter (atomic) adsorbates, an outstanding importance of eh-pair excitations has been reported. Evaluating frictional nonadiabatic energy losses non-self consistently from AIMD trajectories, Kroes and coworkers predicted that the latter exceeds phononic energy losses by factor of 2.5 and 6 for H atom scattering off Cu(1 1 1) and Au(1 1 1), respectively [155]. This picture was later substantiated by Wodtke and coworkers who reported a pronounced disagreement between experimentally measured energy-loss spectra for backscattered H atoms from Au(1 1 1) and adiabatic simulations on a carefully parametrized full-dimensional PES based on effective medium theory including all surface DOFs [156][157][158]. Adding LDFA-electronic friction forces, however, resulted in a spot-on agreement with experimental measurements. Finally, a larger relevance of the non-adiabatic against the phononic dissipation channel was also reported for the recombination of H 2 as compared to N 2 via the Eley-Rideal abstraction mechanism (i.e. involving direct collision and recombination with an adsorbed species from the gas phase) [159,160]. Altogether it nevertheless remains doubtful whether these findings actually reflect an absolutely increased relevance of non-adiabatic energy losses. The increased relative relevance compared to phononic channel could also merely arise from a suppression of efficient phononic dissipation for these lightest adsorbates.

Trapped at the surface: gauging the tools for vibrational damping
Vibrations of molecular adsorbates that are 'trapped' at the metal surface can, in principle, decay via coupling to substrate phonons, and/or through the electronically non-adiabatic excitation of eh-pairs [3]. In addition, vibrational energy may flow to other internal molecular modes, adsorbate-surface vibrations, and even (sufficiently close) neighboring species. In disentangling the contributing role of numerous such mechanisms, it is useful to consider corresponding coupling strengths in terms of smaller or larger mismatch in timescale. High-frequency adsorbate vibrations are, for example, expected to be predominantly relaxed through eh-pair excitations, simply due to their large frequency mismatch with all other vibrations present [150,161]. Advocated by Persson and Persson already more than thirty years ago [122,[162][163][164][165][166][167], this notion was further substantiated by Tully and coworkers by investigating the C-O stretching mode of CO on Cu(1 0 0) [168]. In their study, agreement to the experimentally measured vibrational lifetime of about 2 ps [169] could only be achieved by including dissipative non-adiabatic effects in terms of electronic friction, which lowered the corresponding adiabatic prediction by a whopping six orders of magnitude. More recently, Saalfrank and coworkers reported similar findings for the FT zmode (frustrated translation perpendicular to the surface) of adsorbed hydrogen atoms on Pb using explicit high-level AIMD plus electronic friction (AIMD+EF) simulations [170].
Vibrational damping of high-frequency adsorbate modes thus provides the idealized 'isolated' setting for gauging the accuracy of non-adiabatic theories in GSD. The process is (at least, largely) dominated by a single dissipation channel and atoms are confined to the vicinity of their energy minima, thereby eliminating the need for exploring vast PES regions as required for modeling gas-phase impingement. Most importantly, vibrational lifetimes provide a direct measure for the rate of energy flow and are accurately accessible as benchmark observables from real-time experimental measurements using e.g. pump-probe spectroscopy [150].
Based on this foundation, a plethora of computational studies focused on reproducing experimentally measured vibrational lifetimes for small molecular adsorbates on metal surfaces using different underlying models for the electronic friction coefficients [127,128,141,[147][148][149]171]. In detail, Persson and coworkers [149] as well as Tully and coworkers [127,128,147,148,171] evaluated elements of the ODF friction tensor along certain normal modes for CO [127,128,[147][148][149]171], CN [128,148,149] and NO [128,148] on various coinage and transition metals. Given the numerical challenges involved in adequately sampling the underlying metallic band structure and evaluating the required matrix elements to arrive at ODF friction coefficients [128,129,148], the reported values emerging from various implementations and numerical setups show considerable spread. Notwithstanding, in particular for highest-frequency (and thus anticipated to be 'most non-adiabatic') adsorbate stretch modes, calculated vibrational lifetimes were generally found to agree with experiment to within the same order of magnitude. The computed non-adiabatic decay rates further illustrated a pronounced mode-specificity that, at least so far, eludes any clear, system-transferable trend [127,128,148,171]. Unfortunately, an unambiguous experimental verification of this prediction remains unattainable, in particular for low-frequency modes, as the latter are more likely to couple to surface phonons and thus again open Pandora's box of competing dissipation mechanisms.
Confined therefore to investigating only high-frequency adsorbate stretch modes, a more recent benchmark study [141] focused on assessing the performance of the popular LDFA approach [113,119,123,125]. As mentioned before, electronic structure information here only enters on the level of an atomic embedding density. Given this simplicity, it is rather surprising to note that dynamically evaluated vibrational lifetimes obtained with LDFA-based friction coefficients were found to perform on an equal level as previous ODF-based studies for several experimentally well-characterized systems [141]. A recent study suggested this good performance of the LDFA to arise from finite electronic coherence times, which wash out details of the electronic band structure that are in any case neglected in the LDFA [111]. Further work is, however, required to better assess the full validity and performance of the LDFA. As a pragmatic bottom line, at least for the time being it seems that both the very simple LDFA approach as well as the computationally more demanding ODF model account reasonably well for the non-adiabatic energy losses of high-frequency stretch modes of molecular adsorbates on metal surfaces. Yet, as appealing as it may be to interpret this performance as a justification to rely on these models also for other surface dynamical processes, a respective generalization has to be taken with considerable caution. Benchmarking against vibrational lifetimes avoids delicate situations encountered in simulating, e.g. adsorption processes such as a steep increase of friction coefficients at dissociative transition states [123,131,134,135,138,139] or singularities in the latter at spin transitions [132,133,172]. However, these are situations, that anyhow question the limits and validity of the electronic friction approach in general, not only the LDFA [114][115][116][117]172].

Surface diffusion
Similar to the vibrational motion of adsorbed species, on-surface diffusion is largely controlled by the rate of energy exchange with the underlying substrate [173,174]. However, here no striking argument (such as the frequency mismatch detailed in the previous section) can be invoked in order to safely disregard either one of the available microscopic dissipation channels. Quite on the contrary. Good reasons can be brought forward to argue in favor of both channels: Long adsorbate-surface contact times and the concomitant high embedding densities are suggestive of non-adiabatic coupling within an electronic friction-based description, while the relatively long time scales that are characteristic of diffusive motion should intuitively favor a more efficient phononic dissipation channel.

Thermal diffusion
For a long time the classical notion was that once thermalized with the surface, energy exchange during diffusion of surface species happens predominantly through phononic coupling [173]. This picture was also nurtured from a pioneering study by Tully and coworkers who showed that electronic friction forces do not significantly influence the equilibrium diffusive motion of CO on Cu(1 0 0) [152]. Notwithstanding, even the authors themselves remarked at the time that the absence of a notable influence on the diffusion rate in this study neither necessarily implies a non-existent coupling to eh-pairs nor any generality of their findings. Instead, their description of the phononic bath through harmonic oscillators may have just provided an energy sink effective enough to inhibit correlated multi-jump diffusive motion, such that additional non-adiabatic damping effects could be masked when focusing only on diffusion rates. This is in line with findings by Wahnström. These suggested a dominant role of eh-pair excitations over very inefficient phononic couplings for the diffusion of H on Ni(1 0 0) that is indeed found to proceed via correlated jumps [175].
Similar to the situation for adsorption processes, more sensitive observables than mere diffusion rates thus seem to be necessary to develop further insight into the role of the dissipative channels. Along this line of thinking, Rittmeyer and coworkers recently analyzed Helium-3 spin echo ( 3 He-SE) measurements [176], which provide time-resolved access to the surface (auto-)correlation function.
With the corresponding decay rates very sensitive to the adsorbate-substrate coupling and the underlying diffusion mechanism [177][178][179], the authors evaluated the non-adiabatic contribution to these rates within the LDFA approach. Intriguingly, a high degree of non-adiabaticity suggested a more pronounced role of eh-pair excitations than anticipated by the classical 'textbook notion', at least for this system.

'Hot' diffusion
The prevalent notion of thermal diffusion presupposes the adsorbate's continuous equilibration with the surface -regardless of the detailed origin of the coupling -and correspondingly predicts substrate temperature as the ruling factor. This is contrasted by the concept of hyperthermal transient mobility which can arise from the intrinsic exothermicity of an immediately preceding elementary step like (dissociative) adsorption. Non-instantaneous thermalization of the released chemical energy drives the ensuing 'hot' adsorbate diffusion that is then governed by dissipation to the underlying surface, rather than the substrate's overall temperature. This transient mobility thus intricately couples the elementary reaction steps of adsorption and diffusion -an implication hitherto not considered in the present-day Markovian microkinetic formulations in surface catalysis [9,10]. Concepts embracing such processes are nevertheless becoming increasingly established over the past decades [180]. Ensuing for example the (exothermic) dissociative oxygen adsorption, so-called 'hot' adatom motion has been persistently reported by scanning tunneling microscopy (STM) experiments for several metal surfaces [181][182][183][184][185][186][187]. Working at temperatures that are sufficiently low to suppress thermal diffusion, corresponding STM studies infer hyperthermal diffusion from recorded larger or smaller separation distances between adatom pairs after equilibration -an indirect procedure that has in fact caused quite some controversy [188][189][190].
From a modeling perspective, the accuracy of energy sink models in describing such hypterthermal motion is challenged head on as the entire process hinges on the rate of dissipation to the substrate. Such a sensitive dependence was already demonstrated in the aforementioned seminal work of Tully and coworkers for hyperthermal CO migration on Cu(1 0 0) [152]. Including electronic friction in classical MD simulations was found to significantly quench the molecules' on-surface transient mobility -an effect which did not extend to the case of equilibrium diffusion. A similar conclusion was recently reached for hot diffusion of atomic H on a precovered W(1 1 0) surface [191]. As a result, the Eley-Rideal reaction was predicted as the more likely mechanism for H 2 abstraction over hotatom recombination/desorption. On the other hand, Wahnström and coworkers estimated the eh-pair dissipation mechanism to be relatively unimportant for the hyperthermal adatom motion following the dissociative O 2 adsorption at Al(1 1 1) [175,192]. Even when focusing solely on phononic damping within a Langevin framework, however, corresponding simulations showed no evidence of a high transient mobility that could reconcile the exceptionally large O-O separation distances derived from STM experiments [181]. Lateral displacements were instead found to be primarily limited by the high PES corrugation and concomitant rapid randomization of the adatom motion, while adsorbatephonon coupling was overall reported to be weak. At this point, however, it is worth noting that the semi-empirical potentials used in this work to describe the O-Al interaction could not be extended to modeling the actual (preceding) O 2 dissociation event, thereby prohibiting a realistic sampling of the employed initial conditions. Unless there is a problem in the STM interpretation [188,189], this obvious limitation represents an equally likely cause for the discrepancy to experiment as the approximate treatment of adsorbate-surface energy transfer.
Overcoming this limitation, a pioneering study of Groß employed AIMD to consistently model H 2 dissociation and the subsequent adatom thermalization at Pd(1 0 0) [19]. Based on substrates of about 100-200 metal atoms, this predicted 'hot' H adatoms that transiently diffuse to an average separation distance of three to four surface lattice constants. Adding an LDFA-based electronic friction to this description, Blanco-Rey and coworkers later on arrived at halved H-H separations only [193]. The authors consequently advocated a dominating role for the electronic dissipation channel arguing an efficiency of about five times larger than provided by the phonon bath [193]. Subsequent work in this direction [194,195] specifically highlighted the increasing contribution of ehpair excitations during the later stages of the relaxation process (i.e. for long adsorbate-surface contact times), even if non-adiabatic effects can be neglected during the preceding dissociation event [110,123,134,146,196]. At this point, however, it has to be noted that in present-day combined AIMD+EF simulations, the effective 0 K-Langevin description of the non-adiabatic bath provides an infinite energy sink, thus conflicting with an explicitly described, energy conserving lattice motion.
As already mentioned in Section 2, even self-standing AIMD (and AIMD+EF) simulations remain challenged in providing an altogether satisfactory reference for phononic dissipation. Computationally tractable slab models compromise the description of the phonon band structure and limit their propagation to a finite volume of few substrate atoms in each direction. This gains particular relevance in the present context of hyperthermal diffusion by considering that a (preceding) exothermic surface reaction may easily release several electronvolts of energy. Unphysical phonon reflections at the periodic boundaries of the employed supercell may thereby quickly lead to severe substrate 'overheating', while falsifying the ps-scale equilibration of the actual adsorbate dynamics [197].
In solving this problem, the recent development of the QM/MM embedding approach for metals (originally coined, and hence referred to, as 'QM/Me') [197] represents a big step forward. Here, a DFT-based treatment of the immediate reaction zone (i.e. typically around the adsorbate impingement site) is complemented by an extended substrate that is described on the level of a many-body classical interatomic potential. The latter are generally sufficiently accurate in providing a realistic, material-specific representation of lattice deformations (cf. e.g. Refs. [198,199]), while their numerical efficiency allows for capturing all associated long-range elastic effects. Within a multi-scale modeling philosophy, QM/Me thus augments standard AIMD with a fully quantitative account of phononic dissipation as heat flows from the embedding region and into the macroscopic metal bath, as schematically illustrated in Figure 4.
In a first application to oxygen dissociation over Pd(1 0 0), QM/Me showed phonons dissipating the vast majority of the released chemical energy into the  (100) [190]. A QM description of the immediate reaction zone is based on periodic DFT calculations to yield an accurate description of the short range, adsorbate-induced chemical interactions. This QM-cell is embedded into an Me substrate that is treated at the level of a numerically efficient classical interatomic potential to provide the elastic contribution to the forces due to deformation of the lattice. Released chemical energy is thus dissipated out of the 'hot' reaction zone and into a 'cold' macroscopic heat bath, while atoms in the embedding cell are color-coded according to kinetic energy. bulk, i.e. outside of the QM-cell, already within ca. 1.5 ps after the initial O 2 bond breaking [197]. The observed ps-scale rate of heat transfer to the substrate nevertheless clearly demonstrated that this process is not instantaneous on the timescale of the actual adsorbate dynamics. Furthermore, the predicted equilibrium O-O separation distance of four surface lattice constants could not be reproduced within 'pure' AIMD simulations, thus underscoring the importance of an extended and explicit description of the surface degrees of freedom for microscopic details. Analyzing the role of surface symmetry, a subsequent QM/Me study compared this result for different Pd facets to reveal a striking difference in the transient lifetime of the hyperthermal O state [75,190]. More specifically, a much slower equilibration was found on Pd(1 1 1) that would seemingly contradict the shorter O-O end distances previously reported from STM measurements [185]. Relying on the atomic resolution of the QM/Me heat bath allowed for rationalizing this finding through a mode-specific analysis of corresponding phonon excitations: This identified the dominant dissipation channels as qualitatively different groups of localized surface modes that do not necessarily involve a predominant Rayleigh excitation, as principally assumed for energy sinks in model bath Hamiltonians (cf. e.g. Refs. [68,69,200,201] and references therein). Instead, the complex adsorbate-phonon dynamics give rise to a sensitive dependence on details of the phononic fine structure that may lie on either the high-and/or low-frequency end of the spectrum. This can ultimately lead to intrinsically different rates of dissipation to the bulk that would otherwise go entirely unnoticed from the perspective of experimentally accessible product end distances.

Summary and conclusions
In the present review we set out to provide a comprehensive overview on modeling efforts for describing the conversion of energy at the gas-solid interface. Focusing on technologically relevant metal surfaces, we consistently address energy dissipation through the (competing) excitation of substrate phonons and eh-pairs at the level of individual elementary processes. Principal questions of focus concern the appropriate tools for modeling each of these microscopic dissipation channels, their effective role in driving the dynamics, and how this is ultimately reflected in (experimental) target observables of interest. Throughout the review, we report on numerous studies aimed at answering these questions, but note that we do not even come close to a full reference list for this highly active field. While thus demonstrating the considerable progress that has been made over the past few decades, we specifically highlight the persisting difficulty of gauging the performance of the various theoretical tools that have been proposed or, vice versa, the level of detail that is required to be accounted for. While each of the discussed methods comes with its own limitations and challenges, this problem stems mainly from the largely inconclusive picture drawn on the basis of different experimental target quantities and across different systems. This is, for example, clearly demonstrated for the processes of inelastic scattering/adsorption reviewed in Section 2. Here, the importance and ubiquity of reaction probabilities (i.e. sticking coefficients) as measured from molecular beam experiments have typically provided a key measure to compare against. Attempts to reproduce such sticking data from theory has shown that dissipation from phonons can require in some cases an explicit treatment of substrate motion, while satisfactory agreement is often already achieved on the level of effective energy sink models, or even a complete neglect of energy exchange with the lattice. While some prediction of the phononic relevance can be made based on the adsorbate-substrate mass mismatch, this does not always appear to be accurate [78,86,92]. This uncertainty adds yet another facet of complexity when considering also the contributing role of the non-adiabatic dissipation channel, whose relevance remains to date similarly obscure and inconclusive. It must be noted, however, that our judgment here is likely further clouded by the limited accuracy of approximate exchange-correlation functionals within DFT. With the vast majority of studies in GSD still relying on GGA-type functionals for computational convenience and clearly no universally 'best' functional at sight, these inaccuracies are particularly manifested in the calculation of encountered energy barriers [37,38]. As reaction probabilities very sensitively depend on the latter, it could even be that effective dissipation models are simply 'compensating' for incorrect barrier heights. Experimental sticking coefficients thus overall represent a rather inappropriate measure for gauging the accuracy of dissipation models: They convolute many different effects along the trajectory course of impinging gas-phase species, while their stochastic nature prevents a sensitive response to microscopic details of the dynamics which may even 'wash out' all signatures of phononic and/or non-adiabatic dissipation. Alternative experimental measures such as the, potentially state-resolved, energy distribution of backscattered molecules thus appear as more promising benchmark references [131].
Also in this respect, the vibrational damping studies addressed in Section 3 are of great value. Here, only a limited phase-space region of the PES is encountered during the dynamics, the non-adiabatic dissipation channel provides a clearly dominating contribution and, most importantly, detailed, dissipationspecific experimental support is made available from pump-probe measurements. Thereby derived vibrational lifetimes have thus provided an extremely sensitive measure to non-adiabatic effects and have largely served as benchmark observables in establishing the performance of the numerically efficient electronic friction approach. Notwithstanding, the extent to which this trust can be transferred also to other elementary processes remains mostly unclear, especially when accommodating the adsorbate's description at a (reaction) transition state during the making/breaking of chemical bonds. One crucial requirement in this regard is the validity of the underlying physical picture that implies only weak non-adiabatic coupling. Despite its limitations, however, electronic friction presently remains the only tractable method in terms of computational resources needed for extended systems and times scales also relevant to phononic motion. As such, already this assessment as extracted from studies of vibrational damping, even if not necessarily universal, is of great importance.
Generally establishing such trust is very important also for simulating surface diffusion as discussed in Section 4. Here, long adsorbate-surface contact times call for all the more accuracy in describing energy exchange with an atomresolved substrate. Microscopic details are, in principle, hereby expected to gain increasing importance for phonon effects (showing, e.g. a pronounced dependence on surface symmetry [75,190]) and eh-pair exciatations [195] alike. However, the relevance of the two dissipation channels is again here far from obvious. Unfortunately, only very few first-principles diffusion studies have so far simultaneously included both of these contributions, while existing indications support a more pronounced role for non-adiabatic couplings than hitherto anticipated [173]. In this respect, further comparisons against high-quality, atomically resolved experimental data such as, e.g. provided from STM-derived adsorbate end-distances or 3 He-SE signatures will undoubtedly play a key role in systematic future studies. Largely removing uncertainties with respect to the phononic heat bath, the QM/Me embedding scheme (described in Section 4.2) now represents a promising way forward. Coupling to the efficient electronic friction approach will specifically allow for a statistically meaningful number of trajectories that can better assess the performance of the latter and quantitatively disentangle the contribution of non-adiabatic effects in the presence of a realistic mobile substrate. We further envision corresponding QM/Me+EF studies to establish mechanistic trends for yet more complex systems beyond diatomic adsorbates, while elucidating the potential implications of 'hot' chemistry to many important dynamical processes such as the self-assembly of surface nanostructures, the first steps of oxide nucleation and epitaxial growth, or adsorbate-induced surface reconstructions.