Ultrafast dynamics of optically-induced heat gratings in metals -- more complicated than expected

Diffusion of heat in metals is a fundamental process which is crucial for a variety of applications of metal nanostructures. Surprisingly, however, {\em ultrafast} heat diffusion received only limited attention so far. Here, we show that heat diffusion can be made faster than $e-ph$ energy transfer rate, in which case, it dominates the spatio-temporal dynamics of the temperature. This enables the metals to overcome the conventional limitations of the nonlinear optical response of materials - it can be simultaneously fast and strong. As a specific example, we identify the underlying (femtosecond and few picosecond) time scales responsible for the generation and erasure of optically-induced transient Bragg gratings in thin metal films. Further, we show that heat diffusion gives rise to a significant nonlocal thermo-optic nonlinearity - it affects also the nonlinear optical response such that the overall change of the permittivity (hence, reflectivity of the transient grating) has a significant dependence also on the illumination period rather than only on the illumination intensity.


Keywords
Ultrafast heat dynamics, Plasmonics, Transient gratings Intro Heat generation and dynamics is central to our understanding light-matter interactions in metals, 1-3 specifically for the purpose of separating thermal effects from (non-thermal) electronic effects. 4,5 It is also critical for many applications 6,7 and specifically, for the correct interpretation of the role of heat in nanophotonic applications such as surface femtochemistry 8,9 and plasmon-assisted photocatalysis. [10][11][12][13][14][15][16] The standard model of heat generation and dynamics in metals is the well-known Two Temperature Model (TTM) 17 in which a distinction is made between the electron subsystem and the phonon subsystem (each assigned by a temperature, T e and T ph , respectively); Further, the electron system is assumed to become thermalized instantaneously, hence, before significant energy transfer between the electrons and the phonons occurs. Under these conditions, the dynamic equations for the electron and phonon temperatures are coupled via a simple linear term, namely, ∂T e (t) ∂t = −G e−ph (T e − T ph ) + p abs ( r, t), Here, C e and C ph represent the heat capacity of the electrons and phonons, respectively; G e−ph is the electron-phonon coupling factor, representing the rate of energy exchange between the electrons and the lattice; p abs , is the density of absorbed photon power. These equations describe the decay of the electron temperature following an initial heating stage (which sometimes for convenience is even skipped altogether) due to energy transfer to the phonons; in conjunction, it describes the phonon heating. The TTM was originally in-troduced in Ref. [17] and derived in detail from the semi-quantum mechanical Boltzmann equation for the few picosecond regime in Ref. [18]. Unfortunately, as well known, the assumption underlying the TTM is never strictly valid -the thermalization of the electron subsystem extends roughly over the first picosecond after the excitation pulse exited the metal. [18][19][20][21] As a remedy, using a unique model that ensures energy conservation, we derived in Ref. [ 4,5] an extended version of the TTM (referred to below as the eTTM) from the semi-quantum mechanical Boltzmann Equation.
In this model, the early stages of the thermalization of the electron subsystem are accounted for via the exact electron distribution and the total energy of the non-thermal electrons; the latter is characterized by a fast (pulse-duration limited) rise time and slow decay time, corresponding to the thermalization of the electron subsystem.
This derivation confirmed the phenomenological models presented much earlier in Ref. [ 19,22] as well as the classical derivations in Ref. [23,24]. In contrast to the TTM, the eTTM allows the electron subsystem to be non-thermal, and only assumes (e.g., by adopting the relaxation time approximation, or adopting a rigorous many-body formulation 18,25,26 ) that the electron subsystem can be described by some temperature -the one to which the electron system would have relaxed if it was isolated from the photons and phonons, see discussion in Ref. [4,5]. In that sense, within the validity conditions of the Boltzmann description, the eTTM is an exact coarse-grained description of the energy dynamics; it is valid at all times (in particular, also before all the non-thermal energy is depleted) and the temperature varies instantaneously upon absorption of photons as it reflects the total energy of the electron system.
Most research on the heat dynamics in metals to date focussed on the derivation of the thermal properties, the details of the e − ph coupling and the dynamics of the temperatures and electric fields. However, quite peculiarly, early studies regularly did not account for heat diffusion. This approach suited most of the configurations studied experimentally. 20,27 Indeed, on one hand, nanometric metal particles or thin films are characterized by a uniform electric field and hence, uniform temperatures. In fact, the temperature is uniform even for metal nanostructures extending to several tens or even hundreds of nanometers, due to the strong heat diffusion in metals. 3,28 Furthermore, whenever a large beam was used for illumination, diffusion was a slow process (see Eq. (2) below) occurring only on the edges of the illumination spot, away from where the processes of interest occurred. [18][19][20] However, on the other hand, it is easy to appreciate that there are several scenarios in which heat diffusion cannot be neglected. These include, for example, thick metal layers into which light penetration is far from being complete, 1,27,29,30 metal-dielectric composites, 31,32 diffusion of a localized heat spot in a thin film (studied in Ref. [33]), highly non-uniform illumination of micron-scale metal objects etc.. In such cases, the heat diffusion is expected to cause the strongly illuminated regions to reach lower maximal temperatures (compared to the diffusion-free case) and the weakly-illuminated regimes to initially get hotter (before cooling down with the rest of the system due to heat transfer to the environment), see e.g. Ref. [ 32]; heat diffusion may also affect the overall time scales for the dynamics in a non-trivial way (e.g., to increase the electron temperature rise time, 34 to tailor the slower e − ph energy transfer rate, 35 or more generally, to affect the overall dynamics due to the relevant modal response 32 ).
In this manuscript, we focus on the specific case of the diffusion of a periodic heat pattern (aka transient Bragg gratings, TBGs) in a thin metal film. In this scenario (depicted schematically in Fig. 1), two pump pulses are interfered onto a metal film to create a periodic pattern of absorption; in turn, the generated heat and the temperature dependence of the metal permittivity gives rise to a periodic permittivity modulation that outlives the pulse, but self-erases due to heat diffusion. Such a pattern can enable transient reflectivity of guided modes in the metal film that can be turned on and off on a sub-picosecond time scale, as for free-carrier gratings. 36 This property is appealing for ultrafast switching applications, as it may enable switching speeds significantly faster than what is commonly achieved with free-carrier generation in semiconductors. 37 To the best of our knowledge, TBGs in metal films were studied before only in the context of the slow (nanosecond and microsecond) phonon dynamics (see, e.g., Ref. [38][39][40][41]), where the electron thermalization and sub-picosecond temperature dynamics can be conveniently ignored. Studies of faster dynamics were performed in Ref. [42] for time scales of several hundreds of picoseconds and in Ref. [30] on a scale of few tens of picoseconds. The neglect of the subpicosecond dynamics might have also originated from the un-availability of a proper theoretical tool valid in that temporal range. An exception is the study of Ivanov et al. 43 which was, however, focussed on employing the periodic illumination pattern for patterning the surface of a metal film.
Here, D e and K e are the diffusion coefficient and thermal conductivity of the metal electrons; d is a characteristic length scale associated with the initial temperature non-uniformity; here, it is simply the illumination period. The magnitude of the diffusion constant D e of metal electrons is approximately D e ≈ 100 − 300 cm 2 s −1 (as e.g., demonstrated experimentally 33 ).
With sufficiently strong inhomogeneities of the illumination (specifically, small periods, in the current context), the heat diffusion can become as fast as a few hundreds of femtoseconds.
We show that due to the non-trivial interplay between the temporally-and spatiallynonlocal nature of the thermo-optic response of metals (namely, thermalization and e − ph energy transfer vs. heat diffusion), heat diffusion dominates the dynamics when it is faster than e−ph heat transfer (∼ 1.25ps for Au). This physical scenario explains why the standard trade-off of the nonlinear optical response of materials between speed and strength 44 can be surpassed by metals (see discussion in Ref. [ 45]), enabling the thermo-optic nonlinearity in metals to be both fast and strong. When the diffusion time is slower, the dynamics is dominated by the spatially-local but temporally-nonlocal e − ph energy transfer. In this regime, the gratings do not get totally erased, but rather persist for many picoseconds. This means that, contrary to our initial expectations, TBGs in metals are not ideal as self-erasing TBGs in ultrafast switching applications, except maybe in the extreme spatially-nonlocal limit.
Finally, we also show that the spatial nonlocal nature of the nonlinear thermo-optic response of the bulk metal (not to be confused with the structural nonlocality in metaldielectric composites, see Ref. [46]) affects the overall magnitude of the nonlinear response.
This newly observed effect might have been playing a role in the multiplicity of reported values of the intensity-dependent optical nonlinearity of metals and their composites (see e.g., Ref. [ 47,48]) and provides yet another motivation to avoid describing the nonlinear optical response by a nonlinear susceptibility, a concept which was derived in the context of spatially-and temporally-local (i.e., Kerr) electronic nonlinearities.
While the results of the current manuscript rely on a rather generic "minimal" model applied to metals only, in its last part we describe several aspects of the model that can Here, Γ e represents the rate of exchange of energy between the non-thermal electrons and the thermal ones, or in other words, it represents the thermalization rate of the electron distribution as a whole. 21 It was computed and measured to be on the order of several hundreds of femtoseconds (see e.g., Ref. [20,22,24,49] to name just a few of the studies of the thermalization). In Ref. [19] it was found experimentally that Γ e is independent of the laser fluence in the range of 2.5 − 200 µJ cm 2 , corresponding to an electron temperature rise of up to about 200K. Similarly, Γ ph represents the rate of energy transfer between the non-thermal electrons and the phonons. Due the similarity in the physical origin, it is set to the same value chosen for the energy transfer rate between the thermal electrons and the phonons (Γ T ph ), see below. Together they form the total decay rate of the NT energy (also referred to as the thermalization rate) denoted by The second term on the RHS, p abs , describes the time-averaged power density of absorbed photons. Using the Poynting theorem for dispersive media and assuming that all the absorbed photon energy is converted into heat, together with the assumption of a slowly varying envelope approximation in time for the electric field, it can be written as where ǫ m is the metal permittivity (and ǫ ′′ m is its imaginary part), ω pump is the pump frequency, E is the local electric field and stands for time-averaging over the period; note that p abs ( r, t) maintains the time dependence of its envelope.
The equations for the temperatures in the original TTM are In these coupled heat equations, U N T e serves as the heat source, and C e and C ph are the heat capacities of the electrons and the lattice; G e−ph (mentioned above) represents the rate of energy exchange between the electrons and the lattice which occurs on time scales defined as 20 Γ Te ≡ G e−ph /C e and Γ T ph ≡ G e−ph /C ph , respectively; Since Γ Te ≫ Γ T ph , the latter is essentially negligible. These equations neglect heat coupling to the environment and consequent cooling, as it occurs typically on time scales much longer than considered in the current work.
In order to account for heat diffusion, we add to the standard equations above a standard heat diffusion term. We emphasize that as there is currently no self-consistent derivation of Here, K e and K ph are the thermal conductivities of the electrons and the lattice. Note, however, that due to the relative smallness of the phonon diffusion K ph , it is henceforth neglected.
The model described above is used below to characterize the ultrafast heat diffusion in the periodically-illuminated thin metal film. We emphasize that in order to maintain simplicity and without loss of generality, inclusion of more advanced modelling (e.g., accounting more accurately for thermalization, 21,50 for the temperature dependence of the optical and thermal parameters 51,52 (see also discussion below), for potential anisotropy of the parameters, for differences between different metals, for quantum mechanical effects etc.) is intentionally avoided in order to keep the generality and simplicity of the discussion. Those effects are discussed briefly at the end of manuscript and would be included in future studies.

Numerical results
We now solve Eqs. (7) for the case in which the metal film is optically thin (compared to the optical penetration (skin) depth into the metal) and the illumination is periodic.
Accordingly, we assume that all quantities depend only on x and t and look at the dynamics within a single period d. In addition, we set the pump pulse to be the shortest time-scale in the system (i.e., τ pump ≪ τ dif f , 1/Γ N T , 1/Γ Te , ...); the latter choice is motivated a-posteriori by the weak dependence on the pump duration observed in further simulations (not shown).
The parameters used in the solution of the eTTM (7) are given in Table 1, suitable for Au; the parameters for Ag and Al are quite similar. We assume that the pump illuminates the sample near its plasmon resonance with a local field of ∼ 30MV /m. 2 Importantly, the comparisons below are performed for the same local field (hence, absorbed power density) within a single unit-cell 3 .
Solving the complete eTTM requires knowledge of the temperature dependence of all parameters. In the current manuscript, we neglect the temperature-dependence of all parameters (most notable of which is of the heat capacity), as appropriate for low intensities whereby the temperature rise is modest. 52 Due to the brief duration during which the system temperatures are more than several hundreds of degrees above room temperature, any changes associated with such temperature dependence make only modest quantitative changes to the results shown below.
Under these conditions, the solution for the non-thermal energy dynamics was shown in Ref. [53] to be where (following Ref. [36]) This shows that the grating contrast decays exponentially at the diffusion time scale without modifying the average spatial absorption profile.
For a short pump pulse, U N T e ∼ δ(t)e −Γ NT t , such that the pump pulse and NT energy peak together. For a longer pump pulse, the NT energy becomes maximal later than the peak of the pump pulse. Furthermore, the non-instantaneous thermalization gives rise to a temporally smeared source for the electron temperature compared to the source appearing in the TTM 17 (which is just p abs (5)).  (no illumination). In all cases, T e (x max ) initially builds up due to thermalization (i.e., due to energy transfer from the NT electrons to thermal ones); This occurs on a time scale Γ −1 e which is longer than the pump duration, hence, this is referred to as a temporally-nonlocal effect. 56 The maximal temperature rise with respect to room temperature is ∼ 85%, ∼ 60%, ∼ 45% for the three cases we simulated. The temperature at x max then gradually decreases due to both e − ph energy transfer and electron heat diffusion. The latter effect causes T e (x min ) to build up as well, despite not being initially illuminated (hence, this is a spatiallynonlocal effect). Figs. 2(g)-(i) show the electron temperature grating contrast, defined as Further simulations show that as the pump duration is increased (at a fixed pulse energy), the temporal maximum of T e (x max ) drops only slightly and the temperature at the minimum remains roughly the same, so that the electron temperature grating contrast decreases only slightly as well. Accordingly, we do not pursue the study of longer pulses. In order to quantify the decay dynamics of the grating contrast, we plot in Fig. 3(a) the temporal maximum of T e (x max ) and T e (x min ) for different values of τ dif f . One can see that the dependence of the maximal values is strongest for short diffusion where T e (x max ) drops and T e (x min ) increases. Accordingly, the grating contrast δT g e decreases for shorter diffusion times (Fig. 3(b)). It however does not completely disappear, because of the delayed heating of the unilluminated regions -a transient temperature contrast builds up before the temperature at the non-illuminated regions.
We now define τ e g as the characteristic time scale in which the electron temperature contrast is erased. To determine τ e g , we fit the decay stage of the grating contrast to an exponential function (see Fig. 2(g)-(i)). Potentially unexpectedly, one can see in Fig. 3(c) that when the diffusion time is long, the grating is erased much faster than the diffusion time (2); this happens due to the spatially-local but temporally-nonlocal effect of e − ph energy transfer which occurs roughly on a 1/Γ Te time scale. When the diffusion is faster, then, the grating erasure time approaches τ dif f , as one might expect 4 . We now look again at the later stages of the dynamics shown in Fig. 2(g)-(i). Peculiarly, one can see that the rapid grating contrast decay slows down significantly. Additional numerical simulations of the complete eTTM (not shown) show that the diffusion proceeds at the phonon diffusion rate (as already shown in Ref. [33,42]); clearly, this occurs due to the coupling of the electrons to the phonons (the grating erasure indeed proceeds much more 4 Mathematically, it is indeed easy to see that in the extreme spatial nonlocal limit (d → 0, hence, τ dif f → 0), the e − ph coupling term is negligible with respect to the diffusion term, such that the dynamics of the electron temperature is dictated by the convolution between the temporally smeared source (8) and the impulse response of the system (namely, the Green's function of the single temperature heat equation).
rapidly if the electrons are artificially decoupled from the phonons). Within the few picosecond range, the residual contrast is 1%, specifically, it is a few degrees for the local case but, as one might expect, it is much lower, a fraction of a degree, in the nonlocal case; naturally, the (residual) contrast decreases for a stronger electron thermal conductivity.
In order to complete the description of the heat dynamics, we also show the spatiotemporal dynamics of the phonon temperature. One can observe a behaviour reminiscent of that of the electron temperature in the local and nonlocal limits. Importantly, the actual phonon diffusion is faster than its intrinsic value because of the fast electron heat diffusion and the spatially-local e − ph energy transfer. We also observe the non-zero residual phonon temperature grating contrast in the few picosecond regime; this residual phonon grating contrast is, however, much stronger than for the electron temperature -it is 85%, 45%, 20% of the maximal phonon temperature rise, respectively, for the simulations in Fig. 2, specifically, at a few degree level.

Permittivity dynamics
We use the temperature-dependent permittivity model described in Ref. [33] to study the spatio-temporal dynamics of the metal permittivity ǫ m (T e , T ph ) induced by the change in the temperatures. This dynamics thus determines the thermo-optic nonlinear response of the metal to the incoming light. Specifically, the probe pulse is centered at 850nm for which the contribution of interband transitions to the permittivity can be lumped into a single (temperature-independent) constant ǫ ∞ . Accordingly, we refer below only to ǫ intra m (T e , T ph ). Fig. 5 shows an initial rise of the permittivity due to the T e rise and a corresponding relative permittivity change below 1% for the cases plotted in Fig. 2 (compare to a 90% relative change in the electron temperature). At later stages, the sensitivity to the phonon temperature becomes dominant such that for the 1% change of T ph , a ∼ 0.1% change of the permittivity is obtained. This implies on much higher sensitivity of the permittivity to the phonon temperature, and a permittivity thermoderivative of ∼ 0.1/K; such values are comparable to those found in the ellipsometry measurements described in Ref. [57].
It is interesting to note that the Au permittivity change occurs mostly due to the change in the imaginary part of the permittivity. This is similar to the theoretical and experimental observations in the CW case reported in Ref. [28,52,58]. For Ag, the changes of the real part of the permittivity may play a non-negligible role. 52,59 These well-understood aspects of the problem have important implications to the spatiotemporal dynamics of the metal permittivity. First, one can see that heat diffusion reduces effectively the maximal change of the permittivity -this is a result of the non-locality mediated heating at the illuminated regions. In contrast, the permittivity change at the few picosecond time scale is nearly unaffected by diffusion -this is because the permittivity change is dominated by the phonon temperature (see Fig. 4(d)-(f)). Second, while it seems from Fig. 5(c) and (f) that the permittivity grating contrast (defined as δ∆ǫ m = ǫ m (t, x max ) − ǫ m (t, x min )) vanishes in the spatial non-local limit, a more careful look (Fig. 5(i)) shows that the permittivity grating contrast is reduced, but still significant (a few tens of percent) -this is due to the finite time required for heat to diffuse to the unilluminated regions. In that sense, the rise and decay time scales of the grating contrast are faster than those of the maximal permittivity change itself, a fact that could be exploited for ultrafast switching applications 5 . Nevertheless, the maximal permittivity contrast is ∼ 35% of the maximal permittivity change in the non-local case (Fig. 5(i)), whereas it is only ∼ 6% at the few picosecond regime. This level is far higher than the relative residual electron temperature grating (∼ 0.1%, compare to Fig. 2(f) and (i)). Furthermore, although the residual permittivity contrast is small, it persists for much longer than the short transient peak. As a result, the reflectivity from the permittivity grating will persist beyond the initial picosecond, thus, putting the applicability of TBGs in metal films for ultrafast switching configurations based on self-erasure (see e.g., in Ref. [36]) in question.
Overall, the above results show that the permittivity change depends on the interference patten (i.e., on the period) even though the local field hence absorbed power density are the same for all cases. More generally, the nonlinearity at short times (and to a lesser extent, also in the picosecond time scale) changes with the spatial frequency contents of the illumination (by more than a factor 2). This implies that in addition to the temporal nonlocal nature of the metal response (see Ref. [56]), the nonlinear thermo-optic response of metals has also a spatial nonlocal aspect, an aspect which may contribute to the multiplicity of values assigned for the nonlinearity across the literature (see, e.g., Ref. [47,48]). (T e , T ph ); permittivity range plotted is 0 to 6 · 10 −3 in all plots.

Discussion
Our simulations show that the heat diffusion in metals determines the grating erasure time for τ dif f ≪ 1/Γ Te (∼ 1.25ps in Au). In particular, diffusion determines directly the build up time of the heat in the minimum illumination points and mediates the heating level of the strongly illuminated regimes. This enables metal to overcome conventional limitations on the nonlinear response 45 -the response is both fast and strong. When the diffusion time is longer, the grating erasure time is dominated by spatially-local but temporally-nonlocal effects (Γ Te ).
Our analysis, and especially the identification of the nonlocal nature of the permittivity change, is also essential for understanding the nonlinear thermo-optic response of metals. In particular, it shows that the nonlinear thermo-optic response depends not only on the pump and probe wavelengths, the illumination intensity and duration (see e.g. Ref. [47]), but also on the spatial distribution of the pump, namely, its spatial (in)homogeneity.
From the applied physics perspective, our results show that most of the grating contrast disappears on a sub-picosecond time scale. However, even then, there is a residual grating contrast that persists for many tens of picosecond even for very strong heat diffusion (nonlocal case), i.e., for very short-scale thermal inhomogeneities (the induced grating period in the current case). Therefore, TBG applications based on metals might be more constrained compared to transient dynamics of free carriers in semiconductors. 36 The results shown in this manuscript are robust to the various approximations made in our model. Specifically, the phonon diffusion does not change the dynamics described above (indeed, it is a much weaker effect compared with the electron-induced phonon diffusion); there is also only weak sensitivity to the pump duration and to the metal layer thickness.
Indeed, previous studies (e.g., Ref. [33,34,42] showed that accounting for the finite thickness gives rise only to modest quantitative changes in the temperature dynamics. Yet, our result might be sensitive to anisotropy in the diffusion and to structuring. 34,46 Finally, as the inclusion of heat diffusion in our model was done phenomenologically, in a manner compatible with a classical description of the underling physics, it is clear that a self-consistent model has to be derived from a more fundamental point of view, even going beyond the Boltzmann equation to a Density-Matrix formulation. Such a model might reveal some subtleties that the current model did not capture, and may even provide first theoretical predictions for ballistic (i.e., faster-than-diffusion) transport. These effects were envisioned theoretically, but so far studied with phenomenological models that do not account for the non-thermal energy dynamics (Eq. (3); see e.g., Ref. [1,11,43]). Experimentally, there have been some studies of vertical heat transport based on a configuration in which the pump and probe are incident on opposing interfaces of a metal film (see Ref. [1,27,30]), however, for lateral heat transport, there has been no evidence for ballistic diffusion even when using femtosecond temporal resolution, see Ref. [33].
However, thinner structures, combined with advanced ultrafast techniques applied to plasmonic systems (see, e.g., ), might allow one to observe ballistic transport.
Convincing models and measurements of heat transport in this regime will also be useful for resolving arguments regarding charge transfer on the femtosecond scale which are relevant for non-thermal (the so-called "hot") electron studies, and their application for photocatalysis and photodetection 4,5,60 and should eventually be combined with models of additional transport effects associated with interband transitions (see e.g., Ref. [61,62]).
Finally, our findings motivate the study of heat and charge diffusion in other systems, e.g., for materials with nearly-instantaneous thermalization such as graphite or graphene 63,64 and other 2D materials 65,66 on one hand, as well as for semiconductors such as GaAs for which the thermalization is slower. [67][68][69]