Nature of Non-Adiabatic Electron-Ion Forces in Liquid Metals
SSupplemental Material for “Nature of Non-Adiabatic Electron-Ion Forces in LiquidMetals”
Jacopo Simoni ∗ and J´erome Daligault Theoretical Division, Los Alamos National Laboratory, Los Alamos, NM 87545
PACS numbers: 75.75.+a, 73.63.Rt, 75.60.Jk, 72.70.+m
I. DETAILS ON FIG. 1 OF THE MAIN PAPER
The proton is displaced along the red and green line segments shown in the cartoon and defined as follows:panel (a): FCC structure with proton position along the Miller direction d = [
12 12
0] starting at R o = (0 , a c / , d = [100], R o = (0 , a c / , d = [100], R o = (0 , a c / , a c /
2) (dashed lines) and R o = (0 , a c / ,
0) (full lines).panel (d): BCC structure, d = [100] with R o = (0 , a c / , xx = yy and xz = yz ; in (b), xz = yz ; in(c), yy = zz when R o = (0 , a c / , a c / x in units of the cell size a c = 7 .
65 (FCC), 4 .
82 (cubic) and 6 . II. THE FRICTION TENSOR CALCULATION
In order to compute the friction tensor γ ax , by for the systems described in the main paper we combine a finitetemperature Density Functional Theory (FT-DFT) approach that allows us to estimate the ground state of theelectronic system at any given atomic configuration together with classical Molecular Dynamics (MD) simulations totemporally evolve the ionic positions. The ground state of the electronic system is obtained by solving the followingset of Kohn-Sham (KS) equations. ˆ H KS | Ψ n, k i = (cid:15) n ( k ) | Ψ n, k i , (1)where ˆ H KS is the KS Hamiltonian, | Ψ n, k i is the single particle Kohn-Sham state with energy eigenvalue (cid:15) n ( k ). TheKS Hamiltonian written on a spatial grid acquires the following form H KS ( r ) = − ~ ∇ r m e + v H ( r ) + v xc ( r ) + v ext ( r ) , (2)where the effective KS potential of the system is given by the sum of the Hartree component, v H ( r ), the exchange-correlation component, v xc ( r ), and the external ionic potential v ext ( r ).In order to obtain the friction tensor we numerically evaluate the following quantity γ ax , by ( t ) = − M Re (cid:28) X n = m X k ∈ IBZ W k p n ( k ) − p m ( k ) (cid:15) n ( k ) − (cid:15) m ( k ) f ax nm ( k ) f by mn ( k ) cos (cid:18) (cid:15) n ( k ) − (cid:15) m ( k ) ~ t (cid:19)(cid:29) , (3) γ ax , by = Z ∞ dt γ ax , by ( t ) , (4)that is formally equivalent to Eq. (2) of the main paper, with the advantage to be easier to implement and computenumerically. In Eq. (3) < . . . > is a thermal average over the ionic degrees of freedom, while the first summationis computed over all the possible transitions between the KS bands n and m with n = m and the second over allthe k -points belonging to the Irreducible Brillouin Zone (IBZ). p n ( k ) = 2(1 + e − β e ( µ − (cid:15) n ( k )) ) − is the Fermi-Diracoccupation for the spin unpolarized KS state | Ψ n, k i and W k defines the k -point integration weights.The force matrix elements f ax nm ( k ) associated to atom a are obtained from the following integral in real space f ax nm ( k ) = ˆ e x · Z Ω d r Ψ n k ( r ) ∗ f a ( r )Ψ m k ( r ) , (5)where Ω is the system’s volume and f a ( r ) is the effective electron-ion forces resulting from the electronic shielding ofthe bare electron-ion force centered on atom a . In all the calculations for both aluminum and hydrogen plasmas weemploy local pseudo potentials. a r X i v : . [ c ond - m a t . d i s - nn ] J u l t (fs) g (t)( f s − ) × g(t) = X a , b X x , y γ ax , by (t)g(t) = X a X x γ ax , ax (t) g(t) = X a X x , y x γ ax , ay (t)g(t) = X a , b a X x , y γ ax , by (t) t (fs) Z t d s g ( s )( f s − ) × FIG. 1: (Color online) Analysis of different components of the tensor’s time correlation function (3), upper panel, and of theircumulative sum (4), lower panel. The system considered is aluminum at ρ = 2 . g/cm and T i = T e = 0 . eV .TABLE I: Typical parameters used during the self-consistent KS DFT calculations. The number of bands ( N b ), the numberof k -points ( N k ), the cut-off energy ( E cut ), the number of atoms ( N i ) in the simulation box and the number of selectedconfigurations ( N c ). T e ( eV ) ρ ( g/cm ) N b N k E cut ( Ry ) N i N c H 0.1 1.0 100 64 150.0 128 15H 1.0 5.0 400 64 150.0 128 5H 10.0 1.0 1300 8 150.0 128 10Al 0.1 2.35 250 8 150.0 64 18Al 0.5 2.35 350 8 150.0 64 18Al 1.0 1.0 500 8 150.0 64 7
III. DETAILS OF THE CALCULATIONS
All the QMD calculations presented in the main paper were performed by using the QUANTUM ESPRESSO 5 . , a typical calculation always consists of two main parts. The first part is a standard QMD simulationwhere the atoms evolve according to the Born-Oppenheimer dynamics. An Andersen thermostat is employed to ensurethat the ionic temperature does not change during the temporal evolution, the system is first equilibrated and thenlet evolved for a sufficient amount of time (few picoseconds) allowing the accumulation of a number of well separatedatomic configurations. The set of KS equations (1) are solved until convergence in the ground state electron densityis reached. Then we can compute the forces acting on each ion and update the atomic positions at the successive MDstep.In the second part of the calculation we compute Eq. (4) for the friction tensor by averaging over several atomicconfigurations collected during the QMD run. For a given selected configuration a refined electronic structure calcu-lation is performed where the number of bands, the energy cut-off E cut of the plane wave expansion and the numberof k -points are increased in order to achieve full convergence. Table (I) shows the set of chosen parameters in thecase of different systems analyzed in the main paper (see Fig. (2)). The number of bands N b used depends stronglyon the electronic temperature of the system, the higher the electronic temperature is, the higher the number of bandsrequired in the calculation is in order to converge. The number of atoms used in the simulation box is also an impor-tant parameter, in the case of aluminum plasmas we generally use 64 atoms in the periodic box, while for hydrogenplasmas more atoms are usually necessary at a given temperature in order to generate a richer manifold of states andachieve a better convergence of the friction coefficients. For the same reason at low temperatures, in particular in thehydrogen case, a higher number of k -points is also required.We used in all the cases the Perdew-Zunger Local Density Approximation (LDA) to compute the exchange correla-tion potential, v xc ( r ). However, the friction coefficients have a very weak dependence on the choice of the exchange-correlation functional. IV. THE EXACT SUM RULE
From Fig. (1), lower panel, it is easy to observe that the exact sum rule derived in Ref. 2lim t →∞ G ( t ) = X ax , by γ ax , by = 0 , (6)is satisfied by our tensor to a high degree of precision (see the black solid line in the figure). Here we show onlythe case of aluminum at liquid density and T i = T e = 0 . eV , but the validity of Eq. (6) was verified in all thecalculations presented in the main paper. The quantity P a P x,y = x γ ax,ay is very low in magnitude compared to theother contributions to the tensor, this means that the validity of the sum rule is due to a perfect cancellation ofthe diagonal components of the tensor, P a P x γ ax,ax (red dashed line), and of the remaining out of diagonal terms P a , b =a P x,y γ ax,by (green dotted curve). This behaviour is common to all the other aluminum and hydrogen plasmacases that has been considered. ∗ Contact email address: [email protected] P. Giannozzi, O. Andreussi, T. Brumme, O. Bunau, M. Buongiorno Nardelli, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli,M. Cococcioni, N. Colonna, I. Carnimeo, A. Dal Corso, S. de Gironcoli, P. Delugas, R.A. Di Stasio Jr., A. Ferretti, A. Floris,G. Fratesi, G. Fugallo, R. Gebauer, U. Gerstmann, F. Giustino, T. Gorni, J. Jia, M. Kawamura, H.-Y. Ko, A. Kokalj,E. Kkbenli, M. Lazzeri, M. Marsili, N. Marzari, F. Mauri, N.L. Nguyen, H.-V. Nguyen, A. Otero-de-la-Roza, L. Paulatto,S. Ponc, D. Rocca, R. Sabatini, B. Santra, M. Schlipf, A.P. Seitsonen, A. Smogunov, I. Timrov, T. Thonhauser, P. Umari,N. Vast, X. Wu, S. Baroni,
J. Phys.: Condens. Matter , 465901 (2017). J. Daligault and J. Simoni,
Phys. Rev. E , 043201 (2019). ature of Non-Adiabatic Electron-Ion Forces in Liquid Metals
Jacopo Simoni ∗ and J´erˆome Daligault † Los Alamos National Laboratory, Los Alamos, NM 87545, USA (Dated: July 15, 2020)An accurate description of electron-ion interactions in materials is crucial for our understanding oftheir equilibrium and non-equilibrium properties. Here, we assess the properties of frictional forcesexperienced by ions in non-crystalline metallic systems, including liquid metals and warm denseplasmas, that arise from electronic excitations driven by the nuclear motion due to the presence ofa continuum of low-lying electronic states. To this end, we perform detailed ab-initio calculationsof the full friction tensor that characterizes the set of friction forces. The nonadiabatic electron-ion interactions introduce hydrodynamic couplings between the ionic degrees of freedom, which aresizeable between nearest neigbors. The friction tensor is generally inhomogeneous, anisotropic andnon-diagonal, especially at lower densities.
The study of interactions between electrons and ionsin metals is a cornerstone of condensed matter physics.Despite the small electron to ion mass ratio, the inter-actions are never strictly adiabatic as a continous spec-trum of electronic excitations of arbitrarily low energyis available at the Fermi level to couple with the nuclearmotions [1–3]. Similar couplings influence a host of phys-ical and chemical processes at metal surfaces [1, 4, 5],which has generated a great deal of experimental andtheoretical interest [6–12]. The non-adiabatic transitionsresult in exchanges of small amounts of energy that main-tain thermal equilibrium between electrons and ions, anddrive the irreversible evolution towards thermal equilib-rium from a non-equilibrium state [13–18]. In solid met-als, the non-adiabatic interactions are well understoodin terms of electron-phonon interactions [19]. In metal-lic systems where the electron-phonon picture no longerholds because ions have the ability to travel throughoutthe system, such as in liquid metals and in warm denseplasmas created in various matter under extreme condi-tions experiments, the basic properties of these nonadia-batic electron-ion interactions remain largely unexplored.In general, a detailed description of nonadiabaticcouplings with first principles simulations remains aformidable challenge. Fortunately, a simplified coarse-grained description that avoids the explicit propagationof the electron dynamics is possible [22, 23]. Indeed, 1)the small electron to ion mass ratio – or, more accurately,the large electron to ion velocity ratio – and 2) the pres-ence of a continuum of electronic states imply the exis-tence of two distinct time scales t γ and τ e , respectively re-lated to the slow relaxation of ionic momenta induced byelectronic frictional forces and the fast electronic densityfluctuations, such that the variations of the ion velocitiesover an interval of time δt with τ e (cid:28) δt (cid:28) t γ satisfy theLangevin equations [22] M ¨ R = −∇ R V ii − X n p n ∂E n ∂R − M γ ↔ ( R ) · ˙ R + ξ ( R ) , (1)while the electron dynamics is described by a masterequation ˙ p n = P m [ W nm p m − W mn p n ] for the popula- tions p n of the adiabatic electronic states E n ( R ) [22].By way of illustration, for aluminum at 2 .
35 g . cm − and0 . τ e < t γ ’ R = { R a } = { R ax } denotes the set of Cartesian positionsof all the ions, V ii ( R ) is the interaction energy betweenions, N and M are the total number of ions and theion mass, and γ ↔ = { γ ax,by } is the 3 N × N -dimensionalelectronic friction tensor. The second term in the rhs ofEq.(1) is the adiabatic force; in thermal equilibrium, itreduces to the usual Born-Oppenheimer (BO) force de-scribing the interaction between the ions and the densitythe electrons would have if they were in thermal equilib-rium with the instantaneous ionic configuration R . Theother force terms account for the fact that electrons donot adjust instantaneously to ionic motions. The fric-tional forces − M γ ax,by ( R ) ˙ R by arise from the electronicexcitations induced by ion a ’s own motion ( b = a ), or bythe motion of all other ions ( b = a ) and mediated to a bythe conduction electrons. Finally, ξ is a white-noise ran-dom force caused by the rapidly varying electronic den-sity fluctuations. In thermal equilibrium, p n is given bythe normalized Boltzmann factor e − E n /k B T / Z , and thefrictional and random forces completely determine eachother via a fluctuation-dissipation relation. Each γ ax,by is related to the correlation function of the fluctuatingelectron-ion forces on a and b such as [22, 24] γ ax,by ( R ) = 12 M k B T × (2)Re Z ∞ dt Z Z d r d r ∇ ax V ei ( r ) (cid:10) δ ˆ n e ( r , t ) δ ˆ n e ( r , (cid:11) R ∇ by V ei ( r )where V ei ( r ) = P a v ei ( r − R a ) is the total electron-ion potential energy and the correlation function (cid:10) δ ˆ n e ( r , t ) δ ˆ n e ( r , (cid:11) R describes the dynamics of the elec-tron density fluctuations in the ionic configuration R .In this work, we use first-principles simulations to mea-sure the strength of electronic friction and assess the im-portance of its tensorial properties and of its dependenceon the instantaneous ionic positions R . Following themethod discussed in Ref.[15], the friction coeffcients (2)are calculated from the electronic and ionic structures a r X i v : . [ c ond - m a t . d i s - nn ] J u l obtained with Density Functional Theory (DFT) basedquantum molecular dynamics (QMD) simulations. Fulldetails on the simulations are given in the SupplementalMaterial (SM) [26]. To help comprehend the data, wecompare them with two limiting models. Because liq-uid metals and plasmas are isotropic and homogeneousat large scales, the tensor γ ↔ ( R ) greatly simplifies whenaveraged over a thermal ensemble. Thus, the canoni-cal average of the “self” terms ( a = b ) over all ionicconfigurations with a fixed satisify h γ ax,ay ( R ) i = γ d δ x,y where γ d is independent of R a . Similarly, the canoni-cal average of the “cross” terms ( a = b ) over all config-urations with a and b fixed are diagonal in the coordi-nate system where the x axis is directed along the in-terparticle direction R ab = R a − R b with h γ ax,by ( R ) i = γ k ( R ab ) δ x,x + γ ⊥ ( R ab ) [ δ y,y + δ z,z ]. In addition, we com-pare our results with the model e γ ax,by = − M Z d k (2 π ) (cid:12)(cid:12)(cid:12)(cid:12) v ei ( k ) (cid:15) ( k, (cid:12)(cid:12)(cid:12)(cid:12) ∂ Im χ ( k, ∂ω k x k y e i k · R ab , (3)obtained by approximating Eq.(2) to second order in theelectron-ion potential or, equivalently, by substituting inEq.(2) the density correlation function of the homoge-neous electron gas (jellium) model. Here, (cid:15) ( k, ω ) is thejellium dielectric function and χ ( k, ω ) the ideal gas re-sponse function [24]. At this order, γ ↔ has the same sym-metry properties as h γ ↔ i , and the diagonal e γ d reduces toa celebrated model for the energy loss by slow ions in ajellium [28].We begin with simple illustrative calculations to fa-miliarize oneself with the friction tensor. Figure 1 showsthe self components γ px,py ( γ xy for short) felt by a pro-ton p immersed in perfect crystalline structures of Al atnormal density 2 . . cm − . Three stuctures are consid-ered, namely FCC (panels a and b), simple cubic (c) andBCC (d). The proton position R p is varied in a unit cellalong the rectilinear segments illustrated in the cartoon.At each R p , the thermal electronic structure is calcu-lated assuming an electronic temperature of 0 . γ xy ( R p ). The results high-light important properties of the instantaneous frictiontensor. Firstly, each γ xy generally depends on the protonposition in relation to the crystal ions as it feels a differ-ent electronic environment at different positions and itsinteraction with electrons changes. γ xy ( R p ) grows or de-creases between extrema, whose locations correlate withthe high-symmetry sites. The amplitude of variationsdepend on the spatial directions and can be significant:e.g., in (a) and (b), γ zz varies by 40 %, while γ xx isnearly insensitive to the position; in (c), γ xx changes bya factor 3 . { γ xy } is generallyanisotropic. Differences in the diagonal elements alongdifferent crystallographic directions can be significant: in(d), γ yy and γ zz get up to 1 . γ xx ; in(c), γ xx is up to 3 . γ zz . Thirdly, the xyz cubic bcc fcc F r i c t i o n i n p s (a) fcc red trajectory xxyyzz xyxzyz (b) fcc green trajectory x/a c F r i c t i o n i n p s (c) cubic dashes:green trajectory full: red trajectory x/a c (d) bcc FIG. 1: (Color online) Diagonal and off-diagonal friction co-efficients felt by a proton immersed in perfect crystals of alu-minum at 2 . . cm − and 0 . x in units of the cell size a c friction tensor γ xy is generally non-diagonal: in (a) and(b), γ xy develops sizeable nonzero off-diagonal values atpositions between the high-symmetry sites (about 30% of γ xx at x/a c = 0 . x/a c = 0 .
5) or in (a), when it isequidistant from the six nearest neighors ( x/a c = 0 and1) .We now analyze the electronic friction tensor in warmdense H and Al plasmas under the conditions listed in ta-ble I. In each case, we have calculated all the elements ofthe tensor for N c = 10 ionic configurations equally spacedin time along a 5 ps-long QMD simulation with N = 64or 128 [26]. Figures 2a and 2b show the histograms of ρ T Θ Γ sc κ sc e γ d γ d σ d γ od σ od ρ in g . cm − and temperatures T ineV considered in this work. Θ = E F /k B T with E F theFermi energy. κ − sc is the finite-T Thomas-Fermi screen-ing length in units of the average inter-ionic distance r WS ,Γ s = ( Ze ) e − κ sc /r WS the effective coupling parameter of ions.The five rightmost columns are in ps − . the self diagonal elements ( γ ax,ax with a = 1 to N and x = 1 , ,
3, open symbols) and of the self off-diagonal el-ements ( γ ax,ay with a = 1 to N and x, y = 1 , , x = y ,solid symbols), including all the N c configurations. Asexpected from isotropy and homogeneity, the distribu-tions are independent of the Cartesian coordinate systemand their mean values represent estimates of the canon-ical average discussed above. In all cases, each quantityshows a single-peak distribution with mean and standarddeviation given in table I (four rightmost columns). Forcomparison, the full lines show the normalized Gaussiandistributions obtained with these mean values. For H,the distributions are nearly Gaussian at 5 g . cm − butdepart from a Gaussian law at 1 g . cm − (e.g., the diag-onal components have a right-skewed distribution). ForAl, they are all very nearly Gaussian.In all cases, we find that, unlike h γ ↔ i , the self part of γ ↔ is both anistropic and non-diagonal. The diagonalcomponents spread over a sizeable fraction of the meanvalue. For instance, in Al, the full width at half maximum( ’ . σ d ) is 46% of the mean value at 2 .
35 g . cm − and0 . . . cm − and 1 eV. The increased dispersion at lower den-sity is also seen in H, where, in addition, the tails ofthe skewed distributions at 1 g . cm − extend to over twotimes the mean value. The off-diagonal elements are typ-ically much smaller than the diagonal components. Yet,at lower density, they reach values comparable to the di-agonal elements; e.g., in H at 1g . cm − and 0 . γ d increases by 18% between 0 . . . . cm − , and by only 5% between 0 . . cm − . These variations generally strongly depend onthe details of the electronic DOS, and different behaviorscan be expected in other metals [15].To further understand the distributions, Fig. 3 showsthe variations of γ ax,ay ( R ( t )) for a randomly chosen Alion a as it travels through the Al plasma. The coeffi-cients are plotted versus the distance traveled from an ini-tial position and measured at equally spaced time steps.We find that, over the course of a quite short trajectory(the maximum distance traveled is 1 . r WS ), the dispersion of the frictioncoefficients in Fig. 3 is similar to the dispersion of thedistributions in Fig. 2b (blue symbols). Like in Fig. 1,the variations correlate with the spatial variations of theelectronic fluid along the ion trajectory, which, in theliquid-like state under consideration, consists in a suc-cession of localized oscillations in the transient potentialenergy cages formed by neighbors followed by the passageinto another cage [29].We now make simple comparisons to assess thestrength of frictional forces. Firstly, it is interesting tocompare the diffusion time scale t D = r WS / D , where D is the self-diffusion coefficient, with the typical ve- locity relaxation time scale t γ = 1 /γ d . For liquid den-sity Al at 0 . t D = 0 .
46 ps and t γ = 5 . t γ = ∆ R / D with ∆ R = 3 . r WS [30]. As ex-pected, electronic friciton is weak, yet finite, correspond-ing to a theoretically difficult regime beyond the limitsof standard methods (e.g., Smoluchowski equation [33]).Secondly, the inset in Fig. 2b shows the histogram ofBO forces || F BO a || /M v th measured for the same set ofionic configurations as in the main frames, in units ofthe mean frictional force on an ion with thermal veloc-ity v th = p k B T /M . In H, the distributions are singlepeaked around 50 at 1 g . cm − , 10 eV and around 35 at5 g . cm − , 1 eV. In Al, the distribution is peaked around275 at 2 .
35 g . cm − , 0 . . cm − , 1eV. Thus, albeit small, the frictional forces are not neg-ligibly small compared to the BO forces. As discussed in[15, 24], the small nonadiabatic couplings are responsiblefor the irreversible evolution toward equilibrium of thenon-equilibrium states typically created in experiments[20, 21], in the limit of thermal equilibrium they act as athermostat to keep electrons and ions at the same tem-perature. Although further work would be needed to as-sess their impact on the material properties, one cannotrule out at this stage that they are not strong enough,e.g., to affect the fluctuations that allow the potentialbarrier crossing events underlying particle diffusion andnucleation.We now move on to discussing the “cross” terms γ ax,by that couple the motions of two distinct ions a and b . Thedetailed statistical analysis is challenging because eachpair ( a, b ) has a different orientation in the Cartesiancoordinate system of the simulations and the distinctionbetween diagonal (e.g., xx ) and off-diagonal terms ( xy ) ismeaningless. We thus limit ourselves to an analysis of theensemble averaged data. For each pair ( a, b ), we considertheir coupling in a coordinate system where the x axisis directed along R ab = R a − R b , and denote by γ k a,b , γ ⊥ a,b and γ ⊥ a,b the diagonal elements in this coordinatesystem. As expected by isotropy, upon averaging over allpairs of ions and over several ionic configurations, γ k a,b and ( γ ⊥ a,b + γ ⊥ a,b ) / γ k and of γ ⊥ . We verified that upon averaging the latter dependsonly on the separation distance R ab , that γ ⊥ a,b and γ ⊥ a,b become equal, and that the off-diagonal elements ( k⊥ ,etc.) vanish. Figure 4 shows γ k (solid symbols) and γ ⊥ (open symbols) as a function of R ab /r WS for H and Alsystems; the value at R ab = 0 is set to the mean selfdiagonal friction γ d of Fig. 2. The data are compared tothe predictions (full lines) of the model (3), which yields e γ k ( ⊥ ) a,b = − M Z ∞ dk π k (cid:12)(cid:12)(cid:12)(cid:12) v ei ( k ) (cid:15) ( k ) (cid:12)(cid:12)(cid:12)(cid:12) ∂ Im χ ( k, ∂ω f k ( ⊥ ) ( kR ab ) , with f k ( r ) = sin( r ) r − f ⊥ ( r ) and f ⊥ ( r ) = sin( r ) r − cos( r ) r . As R ab increases, γ k first reaches negative values at distancescorresponding to the first layer of neighbors (the firstpeak of the pair-distribution function g ( r ) (not shown)is at r = 1 . r WS ). The absolute magnitude of the firstminimum is a significant fraction of γ d (65% for H, 45%for Al). The negative values mean that the total fric-tional force on a is reduced (increased) when a and b move in the same (opposite) direction along R ab . Be-yond the first layer, γ k slowly decays with R ab in anoscillatory manner around zero. As for γ ⊥ , it rapidlydecays to values significantly smaller than γ d . Thus,the “hydrodynamic” couplings between ions mediated byelectrons are mainly directed along the direction of sep-aration, is sizeable between closest neighbors, and neg-ligibly small with all the other ions. Regarding the lat-ter, one should nevertheless keep in mind the exact sumrule P a P x,y γ ax,ay = − P a = b P x,y γ ax,by [24], whichcouples all coefficients and results from momentum con-servation. The model remarkably reproduces these fea-tures even beyond the closest neighbors (see inset), whichshows that the strength of the coupling by electronic fric-tion is first a property of the electron gas that mediatesit.In summary, we have presented first-principle calcu-lations of the electronic friction tensor γ ↔ ( R ) in warmdense H and Al to characterize the frictional and randomforces that affect the dynamics of ions in non-crystallinemetallic systems due to their non-adiabatic interactionswith electrons. We have shown that, unlike the ther-mally averaged tensor and independently of the frameof reference, the instantaneous tensor is generally inho-mogeneous, anisotropic and non-diagonal, and that theseeffects are stronger at lower density when electronic den-sity variations are larger. We have found that the nona-diabatic interactions introduce “hydrodynamic” couplingeffects between the different ionic degrees of freedom,which is particularly sizable between nearest neighbors.The model (3) gives a satisfactory description of the ther-mally averaged frictions and could be incorporated intoclassical molecular dynamics simulations.This work was performed under the auspices ofthe U.S. Department of Energy under Contract No.89233218CNA000001 and was supported in part bythe U.S. Department of Energy LDRD program atLos Alamos National Laboratory through the grantNo.20200074ER. ∗ Electronic address: [email protected] † Electronic address: [email protected][1] W. Dou, and J.E. Subotnik,
J.. Chem. Phys. , 230901(2018).[2] Y. Li and G. Wahnstr¨om,
Phys. Rev. Lett. , 3444(1992);[3] Y. Li and G. Wahnstr¨om, Phys. Rev. B , 12233 (1995).[4] A.M. Wodtke, J.C. Tully and D.J. Auerbach, Int. Rev. Phys. Chem. , 513 (2004).[5] J.C. Tully, J. Chem. Phys. , 22A301 (2012).[6] M. Askerka, R.J. Maurer, V.S. Batista, and J.C. Tully,
Phys. Rev. Lett. , 217601 (2016).[7] R.J. Maurer, M. Askerka, V.S. Batista, and J.C. Tully,
Phys. Rev. B , 115432 (2016).[8] S.P. Rittmeyer, D.J. Ward, P. G¨utlein, J. Ellis, W. Alli-son and K. Reuter, Phys. Rev. Lett. , 196001 (2016).[9] D. Novko, M. Blanco-Rey, J.I. Juaristi and M. Alducin,
Phys. Rev. B , 201411(R) (2015).[10] J.D. White, J. Chen, D. Matsiev, D.J. Auerbach andA.M. Wodtke, Nature , 503 (2005).[11] N.H. Nahler, J.D. White, J. LaRue, D.J. Auerbach andA.M. Wodtke
Science , 1191 (2008).[12] N. Shenvi, S. Roy and J.C. Tully,
Science , 829(2009).[13] C.P. Race, D.R. Mason, M.W. Finnis, W.M.C. Foulkes,A.P. Horsfield, and A.P. Sutton,
Rep. Prog. Phys. ,116501 (2010).[14] A. Tamm, M. Caro, A. Caro, G. Samolyuk, M. Klin-tenberg, and A.A. Correa, Phys. Rev. Lett. , 185501(2018).[15] J. Simoni and J. Daligault,
Phys. Rev. Lett. , 205001(2019).[16] R.J. Magyar, L. Shulenburger, A.D. Baczewski.
Contrib.Plasma Phys. , 459 (2016).[17] V. U. Nazarov, J. M. Pitarke, Y. Takada, G. Vignale,and Y.-C. Chang, Phys. Rev. B , 205103 (2007).[18] V. U. Nazarov, J. M. Pitarke, Y. Takada, G. Vignale, andY.-C. Chang, Int. J. of Mod. Phys. B , 3813 (2008).[19] F. Giustino, Rev. Mod. Phys. , 015003 (2017).[20] P. M. Leguay, A. L´evy, B. Chimier, F. Deneuville, D.Descamps, C. Fourment, C. Goyon, S. Hulin, S. Petit, O.Peyrusse, J. J. Santos, P. Combis, B. Holst, V. Recoules,P. Renaudin, L. Videau, and F. Dorchies, Phys. Rev. Lett. , 245004 (2013).[21] B. I. Cho, T. Ogitsu, K. Engelhorn, A. A. Correa, Y.Ping, J. W. Lee, L. J. Bae, D. Prendergast, R. W. Fal-cone, and P. A. Heimann,
Scientific Reports , 18843(2016).[22] J. Daligault and D. Mozyrsky, Phys. Rev. B , 205120(2018).[23] J. Daligault and D. Mozyrsky, Phys. Rev. E , 026402(2007).[24] J. Daligault and J. Simoni, Phys. Rev. E , 043201(2019).[25] J. Simoni and J. Daligault,
Phys. Rev. E , 013205(2020).[26] See Supplemental Material at xxxx for the details of thequantum molecular dynamics simulations, which includesRefs. [24, 27].[27] P. Giannozzi, O. Andreussi, T. Brumme, O. Bunau,M. Buongiorno Nardelli, M. Calandra, R. Car, C. Cavaz-zoni, D. Ceresoli, M. Cococcioni, N. Colonna, I. Carn-imeo, A. Dal Corso, S. de Gironcoli, P. Delugas,R.A. Di Stasio Jr., A. Ferretti, A. Floris, G. Fratesi,G. Fugallo, R. Gebauer, U. Gerstmann, F. Giustino,T. Gorni, J. Jia, M. Kawamura, H.-Y. Ko, A. Kokalj,E. Kkbenli, M. Lazzeri, M. Marsili, N. Marzari,F. Mauri, N.L. Nguyen, H.-V. Nguyen, A. Otero-de-la-Roza, L. Paulatto, S. Ponc, D. Rocca, R. Sabatini,B. Santra, M. Schlipf, A.P. Seitsonen, A. Smogunov,I. Timrov, T. Thonhauser, P. Umari, N. Vast, X. Wu,S. Baroni,
J. Phys.: Condens. Matter , 465901 (2017). [28] T.L. Ferrell and R.H. Ritchie, Phys. Rev. B , 115(1977).[29] J. Daligault, Phys. Rev. Lett. , 065003 (2006).[30] We used the experimental value D ’ .
25 ˚A / ps [31],LDA and GGA QMD calculations give 1 . and .
85 ˚A / ps[32].[31] F. Demmel, D. Szubrin, W.-C. Pilgrim, and C. Morkel, Phys. Rev. B , 014307 (2011).[32] T. Sjostrom and J. Daligault, Phys. Rev. E , 063304(2015).[33] D.L. Ermack and J.A. McCammon, J. Chem. Phys. ,1352 (1978). ax , ay in ps d a t a f r e q u e n c y (a)Hydrogen1 g/cc , 0.1 eV1 g/cc , 10 eV5 g/cc , 1 eV0.10 0.050.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 ax , ay in ps d a t a f r e q u e n c y (b) Aluminum F B / Mv th D a t a f r e q u e n c y H, 5g/cc, 1eVH, 1g/cc, 10eVAl, 2.35g/cc, 0.1eVAl, 1g/cc, 1eV
FIG. 2: (Color online) (a) Distributions of self diagonal (opensymbols) and self off-diagonal (solid symbols) friction coef-ficients in warm dense H. The dashed lines show Gaussianfunctions calculated using the means and standard deviationsof simulation data. (b) Same as in (a) for Al. The inset showsthe histogram of BO forces in units of the mean friction force.
Distance traveled R(t)/r WS F r i c t i o n i n p s (c) xxyyzz xyxzyz FIG. 3: (Color online) Self diagonal (open symbols) and selfoff-diagonal (solid) friction coefficients felt by an ion in Al at2 .
35 g . cm − and 0 . p s Hydrogen, 1 g/cc, 0.1 eV R ab / r WS p s Aluminum, 2.35 g/cc, 0.1 eV