Inclusion of Energy Loss in Models of Laser Irradiated Gold Films via Classical Molecular Dynamics
IInclusion of Energy Loss in Models of Laser Irradiated Gold Films via ClassicalMolecular Dynamics
J. M. Molina ∗ and T. G. White † Physics Department, University of Nevada, Reno, NV, USA (Dated: January 5, 2021)The structural evolution of laser-excited systems of gold has previously been measured throughultrafast MeV electron diffraction. However, there has been a long-standing inability of atomisticsimulations to provide a consistent picture of the melt process, concluding in large discrepanciesbetween the predicted threshold energy density for complete melt, as well as the the transitionbetween heterogeneous and homogeneous melting. We employ two-temperature classical moleculardynamics simulations utilizing a highly optimized embedded-atom-method potential. We are ableto match time-resolved electron diffraction data, and find consistency between atomistic simulationsand experiments, by allowing laser energy to be transported away from the interaction region. Thisadditional energy-loss pathway, which scales strongly with laser fluence, we attribute to fast electronproduction in the target.
I. INTRODUCTION
Ultrafast laser excitation of metals is able to bring ma-terial into a state far from equilibrium [1, 2]. The prefer-ential and rapid heating of one subsystem over the otherleads to a system of highly coupled cold ions immersedin a partially degenerate electron sea [3]. These transientstates occur during the formation of all high energy den-sity plasmas, also known as warm dense matter (WDM),with particular relevance to laser micromachining [4–7]and inertial confinement fusion experiments [8]. In thelaboratory, these transient states serve as a testbed wherequantum mechanical theories of electron-ion interactions,nuclei dynamics, and phase transitions can be validated[9–12].Historically, the response of the electron subsystem hasbeen measured in optical pump-probe experiments [13–17]. However, these model-dependent techniques pro-vided only a surface measurement of the electron proper-ties and limited information on the ionic response. Morerecently, X-ray scattering experiments have measured thebulk electron temperature by observing the electron plas-mon feature [18]. In contrast, the bulk ion temperaturehas only been inferred from the atomic structure, mea-sured through ultrafast electron [10] or X-ray diffraction[19, 20]. When the ionic system is still crystalline, theDebye-Waller factor (i.e., the reduction in the intensityof the Laue diffraction peaks) is a practical, albeit model-dependent, method [21, 22]. However, modeling the de-cay of the Laue diffraction peak, which is ultimately de-pendent on the root-mean-square (RMS) deviation of theatoms from their lattice positions, is dependent on sev-eral physical parameters: • The energy density of the sample - (cid:15) ∗ [email protected] † [email protected] • The electron-ion equilibration rate - g ei ( T e , T i ) • The electronic heat capacity - C e ( T e ) • The ionic heat capacity - C i ( T i ) • The Debye temperature - T D ( T e , T i )where the generally assumed dependence on electrontemperature ( T e ) and ion temperature ( T i ) is given.While the electronic and ionic heat capacities are wellconstrained for gold [23, 24], the remaining three exhibitlarge uncertainties in the literature. In particular, theo-retical predictions of the equilibration rate vary by up toan order of magnitude [24–31].During the last decade, three different experimentshave used ultrafast electron diffraction to measure thedecay of the Laue diffraction peaks in laser-irradiatedgold. They each attempt to elucidate the effects of non-equilibrium species on the interatomic potential - i.e.,measure the non-equilibrium Debye temperature and an-swer the long-standing question surrounding the exis-tence of bond hardening/softening in warm dense gold[10–12, 32]. For comparable fluences, each experimentmeasures similar decays in the Laue diffraction peaks.However, they reach opposing conclusions; this is primar-ily due to different assumptions regarding the behaviorof the electron-ion equilibration rate and the initial en-ergy density deposited into the sample. For example, inthe work of Daraszewicz et al. the assumed initial en-ergy density, E , calculated by taking into considerationpurely the reflected and transmitted light, was reducedby a factor of η [11], (cid:15) = η E , (1)where (cid:15) is the corrected energy density.They utilized a value of η ∼ .
5, measured throughindependent optical absorption methods at lower laserfluences. At higher fluences, they found it necessary totreat η as a free parameter [33]. This extra energy-loss a r X i v : . [ c ond - m a t . m t r l - s c i ] J a n pathway, specific to the target and experimental geom-etry, was attributed to ballistic/fast electron transportaway from the system - this, despite contemporary claimsof the negligibly of η [34].Classical molecular dynamics (MD) is one method ofmodeling these experiments. When coupled to an elec-tronic system through an appropriate thermostat, suchsimulations can - in one dimension - simulate the fulltarget thickness; thus, capturing the necessary stronglycoupled nature of the ions. However, models run at theexperimental conditions, generally assuming η =1, dis-play large discrepancies between experimental results andatomistic simulations [35, 36]. The failure of these sim-ulations to match the experimental results is generallyattributed to inadequacies in the interatomic potential;however, recent efforts using the highly optimized andvalidated potential developed by Sheng et al. [37, 38]display the similar deficiencies [39].In contrast to previous simulation work, we treat boththe electron-ion equilibration rate ( g ei ) and initial energydensity as free parameters in order to form a consistentdescription of melting in warm dense gold. We are ableto obtain excellent agreement between the decay of Lauediffraction peaks obtained experimentally and those cal-culated from synthetic diffraction patterns. We find astrongly energy-dependent electron-ion equilibration ratein addition to considerable differences between the as-sumed energy density, E , and the energy density in thesimulations, (cid:15) . This additional energy loss, which wealso characterize by η , scales strongly with laser inten-sity, consistent with previous work invoking a ballisticelectron description. By introducing the free parameter η , we are able to form a consistent description of the melt-ing process in warm dense gold, resolving long-standingdiscrepancies between MD simulations and experiments.We note that the introduction of η as an additional freeparameter further increases the difficulty of drawing re-liable conclusions regarding the system’s evolution. II. COMPUTATIONAL METHODS
We perform large-scale MD simulations in theLAMMPS software package [40] performed in the mi-crocanonical ensemble. The ions, treated explicitly, arecoupled to an electron sub-system through a Langevinthermostat within the two-temperature model (TTM)[41, 42]. We heat the electrons instantaneously to sim-ulate the laser-matter interaction, and synthetic diffrac-tion patterns are produced in order to compare with ex-perimental electron diffraction data of Mo et al. [10]. Inkeeping with the experimental conditions, we model free-standing 35-nm single-crystalline gold films and compareour Laue peak decay to measured results for the three as-sumed energy densities E =0.18 MJ/kg, E =0.36 MJ/kg,and E =1.17 MJ/kg, referred to here as the low, interme-diate, and high energy density cases.We model the non-equilibrium conditions of the laser- irradiated system via a traditional TTM that describesthe evolution of the electron temperature through, C e ( T e ) ∂T e ∂t = − g ei ( T e − T i ) , (2)where C e ( T e ) is the temperature-dependent electron heatcapacity. All simulations utilize an electronic heat capac-ity derived from ab-initio simulations [23] and validatedin thermal conductivity experiments [17]. Within eachsimulation, the electron-ion equilibration rate is treatedas a constant and we employ a 1 fs timestep throughout.In addition, prior to coupling the two sub-systems, weequilibrate just the ions for 6.5 ps using a velocity-scalingthermostat and Berendsen barostat [43] in order to bringthe system to a temperature of 300 K at atmosphericpressure. This introduces an initial energy density intothe ionic sub-system of ∼ .
05 MJ/kg.The ballistic electrons, accelerated by the optical laser,have a mean free path of ∼
100 nm, approximately threetimes the thickness of the foil; thus, we approximate theheating as isochoric [44, 45]. The high thermal conductiv-ity [46] also allows us to retain a single spatially invariantelectron temperature, and thus avoid complications dueto lattice expansion [39]. This is particularly importantto avoid spurious edge effects which may initiate hetero-geneous melting.Our simulations couple the ions to this single-temperature electron bath via a Langevin thermostatgiven by [47, 48], m ∂ v i ∂t = F i − γ v i + f L ( T t ) , (3)where v i is the velocity of particle i experiencing an in-teraction force of F i , m is the mass of the particles, k b is Boltzmann’s constant, and γ is the friction parame-ter that characterizes the electron-ion equilibration rate.The friction parameter is related to the equilibration ratethrough γ = g ei m/ nk b , where n is the number densityof atoms [41, 42]. The f L ( T t ) term is a stochastic forceterm with a Gaussian distribution and a mean and vari-ance given by, (cid:104) f L ( T t ) (cid:105) = 0 , (4) (cid:104) f L ( T t ) · f L ( T t (cid:48) ) (cid:105) = 2 γk b T t δ ( t − t (cid:48) ) , (5)where T t is the target temperature of the thermostat [41,49].We utilize the highly optimized embedded atommethod (EAM) interatomic potential developed by Sheng et al. [37]. Under equilibrium conditions, the prop-erties of the Sheng potential are closely matched withexperimental values, as shown in Table I. Additionally,as demonstrated in Figure 1, the atomic RMS deviationand corresponding Debye temperature, obtained in X-raydiffraction experiments, are well matched at equilibriumtemperatures up to melt. In non-equilibrium matter, the TABLE I. Experimentally measured properties of Face-Centered Cubic (FCC) gold alongside those reproduced inmolecular dynamics simulations using the interatomic poten-tial developed by Sheng et al.
Property Sheng a Expt.Melt Temperature (K) 1320 1337 b Lattice Constant (˚A) 4.078 4.078 c Liquid Density (g/cm ) 17.1 17.1 d Melting Enthalpy (KJ/mol) 11.1 12.8 ea Reference [37], b Reference [50], c Reference [51], d Reference [52], e Reference [53] applicability of the Sheng potential has been able to suc-cessfully match experimental lattice disassembly timesafter laser irradiation [38].We utilize two simulation geometries, either a paral-lelepiped of (86 a × a × a ), or (86 a × a × a ), where a = 4 .
078 ˚A is the lattice constantof gold predicted by the Sheng potential. The smallervolume contained 3,018,600 atoms, where as the largervolume contained 65,738,400 atoms. In the larger direc-tion, which corresponds to the 35 nm thickness of thegold foils used in the experiment, the simulation box ismuch larger than the size of the target geometry, cre-ating a front and rear surface capable of expansion. Inthe two smaller directions, we utilize periodic boundaryconditions. The smaller geometry is employed when cal-culating the decay of the Laue peaks. The larger volume,which produces the same Laue decay curve as the smallervolume, allows for a higher resolution in reciprocal spaceand is used to create the synthetic spatially-dependentdiffraction patterns.
III. MATCHING EXPERIMENTAL DATA
We use the decay of the (220) and (420) Laue diffrac-tion peaks as a metric to compare the simulation andexperiment. We obtain the intensity of the diffraction
FIG. 1. Comparison of the atomic mean square deviationand Debye temperature, T D , produced by interatomic poten-tial developed by Sheng et al. (Sim.) [37] and the valuesobtained from X-ray diffraction measurements at equilibriumtemperatures up to melt by Syneˇcek et al. (Exp.) [54]. peaks from the static structure factor of the system,which is easily calculated from the square of the FourierTransform of the atomic positions [55, 56]. For 3.2 MeVelectrons, whose wave vector ( ∼ × rad/m) is or-ders of magnitude higher than their scattering vector( ∼ rad/m), the Ewald sphere can be considered flatover a finite region of reciprocal space. As such, we takethe Q x - Q y slice of reciprocal space, where the x and y coordinates run parallel to the surface of the foil. Thescattering intensity is calculated by integrating over aspecific reflection.We considered initial electron temperatures ( T e )which, when converted into energy densities via the elec-tron heat capacity [23], correspond to a range in η of0 ≤ η ≤
1, with a step size of ∆ T e = 100 K. Weconsidered values of g ei in the range of 0 ≤ g ei ≤ × W/m /K with a step size of approximately∆ g ei ∼ . × W/m /K. We quantify the best fit byidentifying which ( T e , g ei ) pair minimizes the RMS dif-ference between the experimental and calculated decayin the intensity of the (220) and (420) Laue diffractionpeaks. For each of the three energy densities investi-gated, this difference is shown in Figs. 2(a-c), with thecorresponding best decay curve and temperature evolu-tion for the sub-systems given in Figs. 2(d-i). The blackdashed line in Figs. 2(a-c) encloses the region in whichthe simulation results lie within the experimental errorbars. In each case, we are able to find an ( T e , g ei ) pairwhich gives excellent agreement with the experimentallymeasured decay curves.For the low energy density case, we find that values of T e = 7800 ±
300 K and g ei = 2 . ± . × W/m /Kproduce the best fit. The obtained value of g ei is ingood agreement with experimental measurements of g ei in room temperature gold [57, 58]. The initial electrontemperature corresponds to a corrected energy densityin the electron subsystem of (cid:15) ∼ . ± .
01 MJ/kg,corresponding to η ∼ .
92. In the intermediate en-ergy density case, we find that T e = 8500 ±
300 K and g ei = 5 . ± . × W/m /K produce the best fit,corresponding to a corrected initial energy density forthe electron subsystem of (cid:15) ∼ . ± .
01 MJ/kg, with η ∼ .
57. Finally, in the high energy density case, we findthat the best fit to experimental data comes from valuesof T e = 9900 ±
300 K, and g ei = 15 . ± . × W/m /K- corresponding to a corrected initial energy density forthe electron subsystem of (cid:15) ∼ . ± .
01 MJ/kg, and thus η ∼ . η and the electron-ion equilibration rateon the assumed energy density, i.e., the energy densitycalculated by taking into consideration purely the re-flected and transmitted light (see Fig. 3). It is difficult toprecisely identify the origin of the η parameter. Previouswork has suggested it may be due to energy dissipationinto the supporting grid by ballistic electrons or electronejection from the target rear [11]. The scaling of the η parameter with assumed energy density, which is pro- FIG. 2. A comparison between experimentally obtained Laue decay curves and those produced through molecular dynamicssimulations. (a-c) The RMS difference between simulated and experimental decay curves for a range of initial electron temper-atures and electron-ion equilibration rates. Here, the RMS difference of the (220) and (420) have been added in quadrature.The area enclosed by the dashed line encompasses the region in which the simulation results lie within the experimental errorbars. (d-f) A comparison, for each energy density, of the experimental and simulated decay of the (220) and (420) Laue peaksfor best fit values of T e and g ei obtained in panels (a-c). (g-i) The corresponding temporal evolution of the electron and iontemperatures. Here, the shaded area in blue represents in the variation in the ion temperature across the sample. portional to laser fluence for constant target thickness,corroborates this conclusion. We note that this interpre-tation makes the η parameter specific to the target andexperimental geometry, and thus complicates the inter-pretation of such experiments. A. Spatially-Resolved Diffraction Patterns
We took the best conditions for each of the three en-ergy density cases and predicted, now using the largersimulation box, fully spatially resolved synthetic diffrac-tion patterns; these were again created by Fourier Trans- forming the atomic positions taken from the moleculardynamics simulations and taking a Q x - Q y slice of recip-rocal space. In each case, we ensured that these largersimulations reproduced the Laue peak decay seen in thesmaller simulations. The diffraction patterns, shown inFigs. 4-6, correspond to the three energy density cases.Also shown, for comparison, are the experimentally ob-tained diffraction patterns.In order for a direct comparison with the measuredspatially-resolved electron scattering data of Mo et al. ,the simulated structure factor was multiplied by the elec-tron scattering form factor, which was calculated usingthe Mott–Bethe formula [59] and utilizing tabulated x- FIG. 3. Temporal evolution of best fit values of η and theelectron-ion equilibration rate, g ei , for the assumed experi-mental energy densities E = 0 .
18 MJ/kg, E = 0 .
36 MJ/kg,and E = 1 .
17 MJ/kg ray elastic scattering cross-sections [60]. We found neg-ligible difference with using either a cold or singularlyionized form factor, which exhibit little difference above ∼ − . In addition, a background term of the form[61, 62],d σ inel dΩ = AQ (cid:34) − B Q ) (cid:35) + mQ + c, (6)was introduced to account for inelastic scattering. Thecoefficients, A = 350 and B = 1 .
3, were kept constantfor all comparisons, whereas small variations were re-quired in the linear component between shots at differ-ent laser intensities. In each case the background term isplotted as a dashed line in the figures and represents asmall component of the overall signal. Finally, the entireimage was convoluted with a pseudo-Voigt profile thatrepresented the spatial profile of the electron beam [63](FWHM
Gauss =17 meV, FWHM
Lorentz =6.5 meV).Initially, we found that the diffraction peaks from thesolid gold were an order of magnitude larger than the liq-uid scattering signal, an effect not seen in the experimen-tal data where the difference is closer to a factor of two.We attribute this discrepancy to misalignment and beamdivergence in the experiment that is not present in thesimulations. Thus, we average over a small ∆ Q z , reduc-ing the intensity of the solid diffraction peaks while leav-ing the liquid diffraction peaks, which are broad in recip-rocal space, unchanged. We choose ∆ Q z , which remainsconstant in our modelling, by matching the spatially-resolved scattering signal at both early (t=0 ps) and late(t=17 ps) times in the high-energy-density case. Withthese changes we are able to match, in each the threeenergy density cases, the angularly-resolved line-outs ob-tained in the experiment. It is important to note thatthe applied modifications to the raw data, used to obtainthe angularly-resolved data, affect only the qualitativecomparisons provided in Figs. 4-6, and not the intensity decay curves used in Fig. 2. B. Heterogenous/Homogenous melting
From the time for the total disappearance of the (200)peak, the three energy density cases were thought torepresent the incomplete, heterogenous, and homogenousmelting processes respectively [10]. Based on this analy-sis, they concluded a threshold for complete melting (as-suming total absorption of the laser energy, η = 1) of ∼ ∼ et al. [35, 64].In our analysis, utilizing the interatomic potential de-veloped by Sheng et al. [37], we are able to match theentire decay of the the experimentally measured Lauepeak. Using the OVITO visualization tool [65] togetherwith polyhedral template matching [66] we can directlyassess the structure and melting process in each of ourcases (c.f. Fig. 7). In contrast to the original work,we find that the low energy case exhibits heterogeneousmelting, with a clear propagation of the melt front (seeFig 7(a-c)) and no nucleation of melt inside the target.The two higher energy densities both exhibit homoge-neous melting, with the intermediate case appearing tobe on the cusp of the two regimes, as determined fromnon-uniform melt regions inside the target. Taking thisinto consideration, and the lower predicted energy den-sity from the η parameter, we find a total energy densitythreshold for complete melting below 0.22 MJ/kg (whichincludes the 0.17 MJ/kg present in the electron subsys-tem and the 0.05 MJ/kg already present in our 300 Kgold sample). This value is in close agreement with theexpected value of 0.22 MJ/kg [67]. For the transition be-tween homogeneous and heterogeneous melting, we finda value that lies between 0.22 MJ/kg and 0.26 MJ/kg.Finally, we note that despite achieving excellent agree-ment with the experimental decay curves, we are not ableto reproduce the low-intensity long-lived diffraction fea-tures present in the experiment (c.f. Fig. 5(b) and (f)).We find that the intermediate energy cases melts within100 ps, whereas the original work concludes that it doesnot melt until 800 ps. To investigate if this was a fault ofthe equilibration rate, we manually adjusted the value of g ei throughout the simulation to best match the experi-mental results; we were unable to match both the earlytime Laue peak decay and late time long-lived diffractionpatterns. In the simulation, this may highlight a defi-ciency in either the TTM model or the the inter-atomicpotential. In the experiment, this could be attributed toinhomogeneous heating. FIG. 4. Comparison of simulated and measured spatially-resolved diffraction pattern data for the low energy density case.(a-c) Spatially-resolved data from the simulation for 100 ps, 1000 ps, and 3000 ps. (e-g) Spatially-resolved data obtained byMo et al. [10]. (d) Angularly-resolved diffraction data obtained from azimuthal integrals of the data shown in the panels (a-c).(h) Angularly-resolved diffraction data obtained from azimuthal integrals of the data shown in the panels (e-g).FIG. 5. Same as Fig. 4, but for the intermediate energy density case and times of 20 ps, 100 ps, and 800 ps)
IV. CONCLUSION
In this work we model laser irritated thin-film gold foilsusing large-scale molecular dynamics making use of thehighly optimized EAM potential of Sheng et al. and atemperature-dependent electronic heat capacity. Usingthe Laue peak decay as a metric, and treating the initialenergy density and the electron-ion equilibration rate asfree parameters, we benchmark the simulations againstexperimental results published by Mo et al. . As expected, we find a strong dependence of the electron-ion couplingrate on the initial energy density, with increased couplingat higher energy densities. However, we also find a strongdependence of the η parameter on initial energy density,with increased energy loss at higher energy densities (orlaser fluence), strongly suggesting an energy loss pathwayassociated with the production of fast electrons.With the introduction of the η parameter, we find thatthe complete melting threshold, and heterogeneous-to-homogeneous melting boundaries, are both lower than FIG. 6. Same as Fig. 4, but for the high energy density case and the times of -2 ps, 7 ps, and 17 ps. previously thought. In addition, the introduction of η asan additional free parameter further increases the diffi-culty in drawing reliable conclusions regarding the evo-lution of the system and the Laue peak decay.Future work should concentrate on two specific inade-quacies of the model. First, an electron and ion tempera-ture dependent equilibration rate, i.e., g ei ( T e , T i ), shouldbe incorporated [25]. Second, although ab-inito simula-tions predict no significant bond hardening at electrontemperatures below 25000 K [11], an electron tempera-ture dependent interatomic potential that incorporatesbond softening can be utilized. Such softening is pre-dicted to occur in gold, at its equilibrium volume, fortemperatures above T e ∼ V. ACKNOWLEDGMENTS
The authors acknowledge M. Z. Mo for providing theexperimental data that was used throughout this work.This material is partially based upon work supported bythe U.S. Department of Energy, Office of Science, Of-fice of Fusion Energy Science under Award Number DE-SC0019268. The authors acknowledge the support of Re-search & Innovation and the Office of Information Tech-nology at the University of Nevada, Reno for comput-ing time on the Pronghorn High Performance Comput-ing Cluster. J. M. Molina acknowledges the funding and support from the Nevada NSF-EPSCoR NEXUS Awardand the McNair Scholars Program at the University ofNevada, Reno. Lo w E ne r g y D en s i t y I n t e r m ed i a t e E ne r g y D en s i t y H i gh E ne r g y D en s i t y (a) (b) (c)(d) (e) (f) (g) (h) (i) FIG. 7. Visualization of gold melting in the low (a-c), in-termediate (d-f), and high (g-i) energy density cases via thetechnique of polyhedral template matching [66] in the OVITOvisualization software [65]. For each energy density case, thetime steps corresponding to the diffraction patterns shownin Figs. 4, 5, 6 are presented. The color green identifies theFCC lattice type, where as the non-green colors identify liquidpresent in the sample.[1] K. H. Bennemann, J. Phys. Condens. Matter , R995(2004). [2] E. G. Gamaly, Phys. Rep. , 91 (2011). [3] S. Ichimaru, Statistical Plasma Physics, Volume II: Con-densed Plasmas (CRC Press, Boca Raton, Florida, 2004).[4] B. Christensen, K. Vestentoft, and P. Balling, App. Surf.Sci. , 6347 (2007).[5] K. Sugioka and Y. Cheng, Light: Sci. & Appl. , e149(2014).[6] X. He, A. Datta, W. Nam, L. M. Traverso, and X. Xu,Sci. Rep. , 35035 (2016).[7] R. R. Gattass and E. Mazur, Nat. Photonics , 219(2008).[8] S. H. Glenzer, B. J. MacGowan, P. Michel, N. B. Meezan,L. J. Suter, S. N. Dixit, J. L. Kline, G. A. Kyrala, D. K.Bradley, D. A. Callahan, et al. , Science , 1228 (2010).[9] M. Mo, S. Murphy, Z. Chen, P. Fossati, R. Li, Y. Wang,X. Wang, and S. Glenzer, Science Advances , eaaw0392(2019).[10] M. Z. Mo, Z. Chen, R. K. Li, M. Dunning, B. B. L.Witte, J. K. Baldwin, L. B. Fletcher, J. B. Kim, A. Ng,R. Redmer et al. , Science , 1451 (2018).[11] S. L. Daraszewicz, Y. Giret, N. Naruse, Y. Murooka, J.Yang, D. M. Duffy, A. L. Shluger, and K. Tanimura,Phys. Rev. B , 184101 (2013).[12] R. Ernstorfer, M. Harb, C. T. Hebeisen, G. Sciaini, T.Dartigalongue, and R. J. D. Miller, Science , 1033(2009).[13] J. Hohlfeld, J. G. Muller, S.-S. Wellershoff, and E.Matthias, Appl. Phys. B: Lasers Opt. , 387 (1997).[14] C. Guo and A. J. Taylor Phys, Rev. B , R11921(R)(200).[15] 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, Sci. Rep. , 1 (2016).[16] N. Jourdain, L. Lecherbourg, V. Recoules, P. Renaudin,and F. Dorchies, Phys. Rev. B , 075148 (2018).[17] Z. Chen, B. Holst, S. E. Kirkwood, V. Sametoglu, M.Reid, Y. Y. Tsui, V. Recoules, and A. Ng, Phys. Rev.Lett. , 135001 (2013).[18] L. B. Fletcher, H. J. Lee, T. D¨oppner, E. Galtier, B.Nagler, P. Heimann, C. Fortmann, S. LePape, T. Ma, M.Millot, A. Pak et al. , Nat. Photonics , 274 (2015).[19] T. G. White, J. Vorberger, C. R. D. Brown, B. J. B.Crowley, P. Davis, S. H. Glenzer, J. W. O. Harris, D.C. Hochhaus, S. Le Pape, T. Ma et al. , Sci. Rep. , 889(2012).[20] T. G. White, N. J. Hartley, B. Borm, B. J. B. Crowley,J. W. O. Harris, D. C. Hochhaus, T. Kaempfer, K. Li,P. Neumayer, L. K. Pattison, and F. Pfeifer, Phys. Rev.Lett. , 145005 (2014).[21] T. G. White, N. J. Hartley, B. Borm, B. J. B. Crowley,J. W. O. Harris, D. C. Hochhaus, T. Kaempfer, K. Li, P.Neumayer, L. K. Pattison, F. Pfeifer et al. , Phys. Rev.Lett. , 145005 (2012).[22] T. G. White, J. Vorberger, C. R. D. Brown, B. J. B.Crowley, P. Davis, S. H. Glenzer, J. W. O. Harris, D.C. Hochhaus, S. Le Pape, T. Ma et al. , Sci. Rep. , 889(2012).[23] B. Holst, V. Recoules, S. Mazevet, M. Torrent, A. Ng, Z.Chen, S. E. Kirkwood, V. Sametoglu, M. Reid, and Y.Y. Tsui, Phys. Rev. B , 035121 (2018).[24] Z. Lin, L. Zhigilei, and V. Celli, Phys. Rev. B , 075133(2008).[25] N. Medvedev and I. Milov, Phys. Rev. B , 064302(2020). [26] A. M. Brown, R. Sundararaman, P. Narang, W. A. God-dard, and H. A. Atwater, Phys. Rev. B , 075120 (2016).[27] N. A. Smirnov, Phys. Rev. B , 094103 (2020).[28] K. P. Migdal, D. K. Il’nitsky, Y. V. Petrov, and N. A.Inogamov, J. Phys.: Conf. Ser. , 012086 (2015).[29] X. Y. Wang, D. M. Riffe, Y.-S. Lee, and M. C. Downer,Phys. Rev. B , 8016 (1994).[30] D. A. Papaconstantopoulos, Handbook of the Band Struc-ture of Elemental Solids (Springer US, Boston, MA,2015).[31] Y. V. Petrov, N. A. Inogamov, and K. P. Migdal, JETPLett. , 20 (2013).[32] V. Recoules, J. Cl´erouin, G. Z´erah, P. M. Anglade, andS. Mazevet, Phys. Rev. Lett. , 055503 (2006).[33] Y. Giret, N. Naruse, S. L. Daraszewicz, Y. Murooka, J.Yang, D. M. Duffy, A. L. Shluger, and K. Tanimura,Appl. Phys. Lett. , 253107 (2013).[34] Z. Chen, V. Sametoglu, Y. Y. Tsui, T. Ao, and A. Ng,Phys. Rev. Lett. , 165001 (2012).[35] Z. Lin and L. V. Zhigilei, Phys. Rev. B , 184113 (2006).[36] S. Mazevet, J. Cl´erouin, V. Recoules, P. M. Anglade, andG. Zerah, Phys. Rev. Lett. , 085002 (2005).[37] H. W. Sheng, M. J. Kramer, A. Cadien, T. Fujita, andM. W. Chen, Phys. Rev. B , 134118 (2011).[38] Z. Chen, M. Mo, L. Soulard, V. Recoules, P. Hering, Y.Y. Tsui, S. H. Glenzer, and A. Ng, Phys. Rev. Lett. ,075002 (2018).[39] Q. Zeng, and J. Dai, Sci. China Phys Mech. & Astr. ,263011 (2020).[40] S. Plimpton, J. Comp. Phys., , 1 (1995).[41] A M Rutherford and D M Duffy, J. Phys.: Condens.Matter , 496201 (2007).[42] D M Duffy and A M Rutherford, J. Phys.: Condens.Matter , 016207 (2006).[43] H. J. C. Berendsen, J. P. M. Postma, W. F. van Gun-steren, A. DiNola, and J. R. Haak, J Chem Phys, ,3684 (1984).[44] J. Hohlfeld, S.-S. Wellershoff, J. G¨udde, U. Conrad, V.J¨ahnke, and E. Matthias, Chem. Phys. , 237 (2000).[45] J. Hohlfeld, J. G. M¨uller, S.-S. Wellershoff, and E.Matthias, Appl. Phys. B , 387 (1997).[46] Y. V. Petrov, N. A. Inogamov, and K. P. Migdal, Two-temperature heat conductivity of gold (PIERS Proceed-ings, sent to Proceedings, 2015).[47] T. Schneider and E. Stoll, Phys. Rev. B , 1302 (1978).[48] P. Mabey, S. Richardson, T. G. White, L. B. Fletcher, S.H. Glenzer, N. J. Hartley, J. Vorberger, D. O. Gericke,and G. Gregori, Nat. Commun. , 14125 (2017).[49] A. Medved, R. Davis, and P. A. Vasquez, Fluids , 40(2020).[50] E. A. Brandes and G. B. Brook, Smithells Metals Ref-erences Book, 7th ed. (Oxford, Butterworth-Heinemann,London, 1998).[51] Y. S. Touloukian, R. K. Kirby, R. E. Taylor, and P. D.Desai,
Thermal Expansion, Metallic Elements and Alloys (Plenum, New York, 1975).[52] F. C. Campbell,
Elements of Metallurgy and EngineeringAlloys (ASM International, 2008).[53] E. A. Brandes,
Smithell’s Metal Reference Book , (But-terworths, London, 1983).[54] V. Syneˇcek, H. Chessin, and M. Simerska, Acta Cryst.
A26 , 108 (1970).[55] B. E. Warren,
X-Ray Diffraction (Dover Publications,Mineola, New York, 1990) [56] S. H. Glenzer, and R. Redmer, Rev. Mod. Phys. , 1625(2009).[57] T. G. White, P. Mabey, D. O. Gericke, N. J. Hartley,H. W. Doyle, D. McGonegle, D. S. Rackstraw, A. Hig-ginbotham, and G. Gregori, Phys. Rev. B , 014305(2014).[58] J. Hohlfeld, J. G. M¨ueller, S.-S. Wellershoff, and E.Matthias, Appl. Phys. B , 387 (1997).[59] M. Vos, R. McEachran, E. Weigold, and R. Bonham,Nucl. Instrum. Methods Phys. Res. B , 62 (2013).[60] P. J. Brown, A. G. Fox, E. N. Maslen, M. A. O’Keefe,and B. T. M. Willis, Int. Tables for Crystallogr. C , 554(2006).[61] Ludwig Reimer and Helmut Kohl, Transmission ElectronMicroscopy - Physics of Image Formation (Springer US,Boston, MA, 2008). [62] R.F. Egerton,
Electron Energy-Loss Spectroscopy in theElectron Microscope (Springer US, Boston, MA, 2014).[63] T. Ida, M. Ando, and H. Toraya, J. Appl. Crystallogr. , 1311 (2000).[64] X. W. Zhou, H. N. G. Wadley, R. A. Johnson, D. J.Larson, N. Tabat, A. Cerezo, A. K. Petford-Long, G. D.W. Smith, P. H. Clifton, R. L. Martens, and T. F. Kelly,Acta Mater. , 4005 (2001).[65] A. Stukowski, Modelling Simul. Mater. Sci. Eng. ,015012 (2010).[66] P. M. Larsen, S. Schmidt, and J. Schiøtz, ModellingSimul. Mater. Sci. Eng. , 055007 (2016).[67] Z. Lin, E. Leveugle, E. M. Bringa, and L. V. Zhigilei, J.Phys. Chem. C114