Generalized collective modes approach : Mode contributions to time correlation functions in liquid lead

Contributions of different collective excitations to the density-density and heat-heat time correlation functions in pure liquids are studied within an approach of generalized collective modes. It is shown, that a kinetic relaxing mode, caused by slow density fluctuations, defines almost completely the shape of density-density time correlation function for wavenumbers close to the main peak position of the static structure factor. Analytical threevariable model is used to explain negative amplitudes of contributions to dynamical structure factor from short-wavelength excitations. Such contributions can appear in the range of wavenumbers close to the main peak position of static structure factor.


Introduction
Density-density time correlation functions F nn (k, t) are the subject of primary interest in studies of collective dynamics in liquids.Being obtained in an analytical theory or in computer experiment the F nn (k, t) can be used to explain the experimental dynamical structure factor S(k, ω), which is, in fact, the spectral function of F nn (k, t) with k and ω being wavenumber and frequency, respectively.It is obvious that any time correlation function, obtained in molecular dynamics (MD) simulations, contains already in its shape all the information about various kinds of collective excitations typical for the liquid considered.The simplest picture of different mode contributions to the time correlation functions can be obtained within the hydrodynamic treatment [1,2] when only the most slow processes, which correspond to fluctuations of conserved dynamical variables, are taken into account.Three main contributions to the density-density time correlation function F nn (k, t) are observed only in a pure liquid for that case [3,4], namely, where c s , Γ, D T , γ, and b are the adiabatic sound velocity, sound attenuation coefficient, thermal diffusivity, ratio of specific heats, and some constant, depending on thermodynamic parameters, respectively.One can easily distinguish the contribution to F nn (k, t) from the purely relaxing thermodiffusive mode d T (k) with when k is small, and two oscillating contributions (symmetric and asymmetric ones) from the propagating sound excitations z ± s (k) with in the hydrodynamic region.The asymmetric term in (1) has the leading order ∼ k and does not contribute to the static value F nn (k, 0) because of an asymmetry with respect to t.Thus, in the hydrodynamic region the collective dynamics of a pure liquid is well described by the set of three hydrodynamic collective excitations, which correspond to the most slow microscopic processes in a liquid.Note also that, according to the Bogolubov's hypothesis about hierarchy of relaxation times, the expression (1) means also, that the other mode contributions to the densitydensity time correlation function F nn (k, t), caused by non-hydrodynamic collective excitations, can be neglected in the hydrodynamic limit (such contributions should be at least higher order of magnitude than k).
Beyond the hydrodynamic region the short-time processes become more important in liquid dynamics.Kinetic collective modes (relaxing and propagating ones) begin to contribute sufficiently to all time correlation functions.Therefore, beyond the hydrodynamic region, the simple expression (1) cannot be applied to the study of dynamical properties of liquids and, in particular, for the estimation of a collective mode spectrum.As examples of kinetic excitations known in the literature, we can mention optic-like excitations in ionic liquids, the so-called 'fast sound' and 'slow sound' modes in binary mixtures, and shear waves in transverse dynamics of liquids.All these kinetic propagating collective modes cannot be described within the standard hydrodynamic treatment.
One of the most efficient methods for the study of collective dynamics of liquids in a wide range of k and ω is the approach of generalized collective modes (GCM) [5,6].This method is based on the concept of generalized collective excitations and allows us to estimate the local coupling effects (at fixed k-value) between different relaxing and propagating collective excitations, either hydrodynamic or kinetic ones.Analytical treatment of simplified models, derived within the GCM approach, allows one to also understand different mechanisms responsible for spectra formation in different k-regions.
Dynamical properties of liquid Pb at temperatures 623 K and 1170 K were studied experimentally using the thermal neutron scattering technique and reported in [7] and [8], respectively.It was shown, that for the low-temperature state the effect of collective excitations was visible on dynamical structure factor S(k, ω) as Brillouin peaks up to k ≈ 1.5 Å−1 , while for high-temperature lead, at 1170 K, no side maxima were observed because of the rather wide central peak.That is why according to [8] "a temperature study at small Q values would give an important information on the possibility for a liquid metal to sustain these modes of motion".Another issue of the role played by the structure factor in dispersion law of collective excitations in liquid Pb was raised in [8].
Recently the experimental and theoretical studies were performed on liquid Ga [9], Li [10] and Al [11] with the purpose to find features of non-hydrodynamic behaviour in collective dynamics of the liquids studied.The non-hydrodynamic effects were explained by an additional non-acoustic branch in the spectrum of collective excitations [9] and by the existence of second relaxation time [10,11].
This study of mode contributions to the collective time correlation functions in liquid lead is the follow-up of our recent study on spectra of collective excitations in Pb [13].The scheme of this study is the following: (i) using the numerical results, obtained for the generalized collective modes in [13], to calculate the corresponding amplitudes of mode contributions to the density-density and heat-heat time correlation functions and (ii) to compare the results of calculations with the data expected from the hydrodynamic treatment.This allows us to conclude about the role of non-hydrodynamic effects when wavenumber k increases.Hence, the goal of this study is to investigate the k-dependence of the mode amplitudes, which describe the separated contributions from the generalized hydrodynamic and kinetic collective excitations to time correlation functions in pure liquids.For this purpose we consider a model of liquid lead at two temperatures.We focus special attention on the role of kinetic collective modes and discuss this problem more in detail.
The paper is organized as follows.In the next section the generalized expressions for time correlation functions, which have a wider range of application in comparison, e.g., with (1), are derived.In section III the amplitudes of mode contributions are studied in detail for the density-density and heat-heat time correlation functions.We also discuss the origin of kinetic purely relaxing mode d 2 (k), found in our recent GCM-studies of pure liquids [12][13][14].The last section contains some concluding remarks.

Time correlation functions within the GCM approach
Within the N v -variable approximation of the GCM approach the solutions for time correlation functions can be written (see, e.g., [15]) in an analytical form via the eigenvalues z α (k) and eigenvectors of the generalized hydrodynamic matrix.In this case any time correlation function, constructed on dynamical variables from the N v -variable basis set, is expressed as a weighted sum of N v mode contributions Each term in the expression ( 4) is associated with the relevant collective excitation.The k-dependent function G α ij (k) is the complex weight coefficient of mode contribution to F ij (k, t), associated with the eigenvalue z α (k).In a three-variable approximation for the hydrodynamic set A H = {n(k, t), J l (k, t), e(k, t)} the expression (1) can easily be derived from (4), assuming k to be small.Hence, one may conclude that in the GCM approach the well-known hydrodynamic expressions for time correlation functions (see, e.g., [2][3][4]) can be properly generalized for the case of an arbitrary number of generalized collective excitations.
Let us rewrite now the expression (4), taking into account some properties of the eigenvalues z α (k) and amplitudes G α ij (k).Assuming that among N v eigenvalues z α (k) there are N p pairs of complex conjugated eigenvalues (propagating modes) and N r purely real ones (relaxing modes) d α (k), α = 1, ..., N r (note that N r = N v − 2N p ), one can rewrite equation (4) in the form The new amplitudes A α ij (k), B α ij (k) and C α ij (k) in ( 5) are simply connected with the functions G α ij (k) and depend only on wavenumber k.In further consideration we will call the terms with the amplitudes B α ij (k) and C α ij (k), caused by the propagating modes z α (k), as symmetric and asymmetric contributions, respectively.Obviously, By taking Fourier transform of equation ( 5) one obtains the expression for a spectral function Fij (k, ω) with the separated mode contributions.It is seen from ( 5) that Fij (k, ω) will contain the contributions from N r central Lorentzians, 2N p non-central Lorentzians (symmetric contributions) at frequencies ±ω α , and 2N p non-Lorentzian corrections (asymmetric contributions), respectively.

Results and discussion
We performed MD simulations in a standard microcanonical ensemble for liquid Pb at two thermodynamic points: a high-temperature state at 1170 K with number density n = 0.0289 Å−3 , and a state above melting temperature at 623 K with number density n = 0.03094 Å−3 .In molecular dynamics we studied a system of 1000 particles interacting through oscillating potential Φ ij (r) at a constant volume V = L 3 .The smallest wavenumbers achieved in MD were k h min = 0.1928 Å−1 and k l min = 0.1973 Å−1 for high-and low-temperature states, respectively.The time evolution of hydrodynamic variables and their time derivatives was observed during the production run over 3•10 5 steps for each temperature.The effective two-body potential was taken in an analytical form from [16].This potential very well reproduced the experimental static structure factor of liquid lead over a wide temperature range [16].

Longitudinal collective excitations in liquid lead
Within the nine-variable approximation [17] of the parameter-free GCM method the basis set of dynamical variables for the case of longitudinal dynamics in pure liquids consists of the following operators: where n(k, t), J l (k, t), e(k, t) are the hydrodynamic variables, namely, the densities of particles' number, longitudinal current and energy, respectively.The dots denote the order of time derivative of a relevant operator.The basis set of dynamical variables is applied to generate the eigenvalue problem from the generalized Langevin equation in Markovian approximation [6,17].
The eigenvalues of the generalized hydrodynamic matrix T(k), generated by the basis set (7), were calculated for liquid Pb at low-and high-temperature states.The results are shown in figure 1a and 1b, respectively.As functions of k, these eigenvalues form the spectra of generalized collective excitations.It is seen in figure 1, that for both temperatures the spectra contain in general four branches of generalized propagating modes: three of them exist over the whole k-range considered, and one branch, denoted as z 1 (k), has the propagating gap in small-k region and appears for k larger than some temperature-dependent value k H (k H ≈ 0.2 Å−1 for T l = 623 K and k H ≈ 0.4 Å−1 for T h = 1170 K, respectively).Imaginary and real parts of the complex-conjugated eigenvalues are shown in figure 1a and 1b by the same symbols connected by spline-interpolated solid lines.Three purely real eigenvalues, marked for convenience as d α (k) with α = 1, 2, 3, are shown in the lower frames of figure 1a and 1b by symbols connected by spline-interpolated dotted lines.One can see, that the purely real eigenvalues d 1 (k) and d 3 (k) exist only for small wavenumbers (inside the propagating gap for z 1 (k)), while at some k-value they merged.At this k-point the two relaxing modes disappear, and the pair of propagating excitations z ± 1 (k) emerges instead.Only one relaxing mode d 2 (k) exists in the whole k-region studied.The physical meaning of this mode will be discussed more in detail below.
From the behaviour of eigenvalues at k → 0 one can establish that the pair of propagating modes z 2 (k) corresponds to the generalized sound excitations with the linear dispersion ω s (k) in hydrodynamic region.The straight dash-dotted lines in the upper frames of figure 1a and 1b allows one to estimate the speed of longitudinal acoustic waves in liquid Pb (see [13]).
The pair of propagating modes z ± 2 (k) together with the thermodiffusive mode d 1 (k) form the set of generalized hydrodynamic collective excitations.All the other eigenvalues correspond to the kinetic modes, the damping coefficients of which in contrast to generalized hydrodynamic ones tend to some finite values when k goes to zero, so that these excitations have the finite time of life and do not contribute into the hydrodynamic long-time behaviour.However, as it is seen in figure 1, the real parts of generalized hydrodynamic and kinetic modes become comparable for intermediate and large values of k.Hence, the role of the kinetic modes increases beyond the hydrodynamic region.
In [13] it was shown, that the propagating modes z ± 1 (k) are caused by heat fluctuations and describe the low-frequency heat waves.The branch z ± 2 (k) is well reproduced within the viscoelastic treatment and corresponds to generalized sound excitations.One can also see that in the region of the main peak location of the static structure factor k ≃ Q p the branch z ± 2 (k) has a minimum, which is much more pronounced at the higher temperature.This implies that for very high temperatures one can expect in liquid Pb a propagation gap for generalized sound excitations with wavenumbers k ≃ Q p .The high-frequency branches z ± 3 (k) and z ± 4 (k) are mainly of a thermal and a viscous origin, but their positions are strongly affected by lower lying modes due to the mode-coupling effects.

Mode contributions to the density-density time correlation function
In figure 2 the leading mode contributions to the density-density time correlation function, calculated for liquid lead at the higher temperature T h = 1170 K, are shown for three values of k.These results are obtained based on our expression (4) for time correlation functions with the separated mode contributions.The MD-derived functions F nn (k, t) are shown by solid lines.Dashed lines, which correspond to the GCM functions (4), almost coincide with MD-functions, which means a very good quality of the nine-variable approximation used in our calculations.Note that the GCM approach does not require any adjustable or fitting parameters.Dotted lines show the total contribution from the generalized sound excitations z 2 (k), i.e., the sum of symmetric and asymmetric terms associated with the sound modes.Shortdash-dotted lines correspond to the heat excitations: inside the propagation gap for the heat waves, i.e., for k < k H , this line represents the contribution, caused by the thermodiffusive relaxing mode d 1 (k), while for k > k H it gives the total contribution from the heat waves z 1 (k) (the sum of symmetric and asymmetric terms).Purely relaxing mode d 2 (k) causes the contributions shown by long-dash-dotted lines.It is well seen in figure 2a, that in complete agreement with predictions of the hydrodynamic theory [1,2], the shape of density-density time correlation function for the smallest k-value, considered in our study, is almost completely determined by the contributions from the hydrodynamic excitations (the pair of sound excitations z 2 (k) and the thermodiffusive mode d 1 (k)).For k-values slightly beyond the propagation gap (k > k H ) the heat waves cause a rather small negative contribution, while the shape of function F nn (k, t) is formed mainly by the sound excitations z 2 (k) and the relaxing kinetic mode d 2 (k) in comparable amounts.When k increases further and gets closer to the main peak of static structure factor k ∼ Q p , the contributions of the low-frequency heat waves and the sound excitations become comparable and much smaller than the term associated with the kinetic relaxing mode d 2 (k), so that the long-time tail of F nn (k, t) in the region k ≃ Q p is completely defined by the kinetic relaxing mode d 2 (k).Similar behaviour was also observed for the lower temperature T l = 623 K.
In order to get a more detailed picture of separated mode contributions to the density-density time correlation function we also calculated the normalized amplitudes, describing the mode contributions to the function F nn (k, t)/F nn (k, 0).The results as functions of wavenumber k are plotted in figure 3 for several low-lying collective modes (see figure 1) at two temperatures considered.Note that the weight amplitudes of the high-frequency kinetic propagating modes z 3 (k) and z 4 (k) are nearly two-or tree orders of magnitude smaller than the ones, caused by the sound excitations and low-frequency heat waves.It is worth to mention that the normalized amplitudes satisfy the following sum rule (compare with ( 6)): It is obvious, that the asymmetric amplitudes do not contribute to static quantities.Note also, that due to the existence of the propagating gap for low-frequency heat waves with a relatively small width k H , one has to distinguish two k-regions with different numbers of relaxing and propagating modes.Three contributions from the relaxing modes (N r = 3) and six contributions (symmetric and asymmetric ones, N p = 3), caused by the propagating excitations, exist for k < k H . Beyond the propagation gap the situation changes and we have in (5) one contribution from the relaxing mode and eight terms associated with the propagating modes (N r = 1 and N p = 4) within the nine-variable approximation (N v = 9).For readers' convenience the normalized amplitudes of relaxing modes are shown in figures 3a and 3b by symbols connected by spline-interpolated dashed lines, while symmetric/asymmetric contributions from the propagating modes are shown by closed/open symbols connected by spline-interpolated solid lines.
In general the following specific features in behaviour of the normalized mode amplitudes can be seen in figure 3: (i) In full agreement with the hydrodynamic theory [1] the normalized symmetric amplitude B2 nn (k), which describes the contribution from generalized sound excitations z 2 (k) (solid line with closed triangles), and the amplitude Ā1 nn (k) of the thermodiffusive mode d 1 (k) (dashed line with open boxes), tend to the values 1/γ and (1 − 1/γ), respectively, when k tends to zero [13].It is especially well seen for higher temperature T h = 1170 K, because for the state with T l = 623 K (just 22 K above the melting temperature) the width of the hydrodynamic k-region is much smaller and the smallest wavenumber k min = 0.1973 Å−1 , reached in MD in this last case, is in fact far beyond the range of hydrodynamic behaviour.The asymmetric contribution C2 nn (k) from the acoustic excitations behaves like a linear function of wavenumber in small-k region for T h = 1170 K, which is also in agreement with the hydrodynamic predictions.We note that the terms associated with the kinetic modes become negligible for small k, so that the theory is in complete consistency with the hydrodynamic picture; (ii) The main contribution to the shape of the density-density time correlation functions in the region of the main peak of the static structure factor k ∼ Q p is caused by the kinetic relaxing mode d 2 (k), which is of a different origin in comparison with the thermodiffusive one d 1 (k).We will discuss this last statement below more in detail; (iii) When k increases from hydrodynamic values up to the edge of propagating gap for the heat waves k H , one can see in figure 3, that the magnitude of mode contributions from the heat relaxing modes d 1 (k) and d 3 (k) increase rapidly with opposite signs and may even diverge when k → k H (in particular, such possibility is discussed in section 3.4).On the other hand, just beyond the propagating gap, the symmetric and asymmetric contributions from the low-frequency heat excitations also take large enough negative and positive values, respectively.This implies, that at the point k = k H , where two relaxing eigenvalues merge and the pair of low-frequency heat waves appears, the corresponding amplitudes of relaxing modes diverge with opposite signs, however, keeping their sum finite due to (8).
(iv) The behaviour of amplitudes for the heat waves in the region of wavenumbers k ∼ 2.7 Å−1 at T h = 1170 K (figure 3b) is in some sense similar to the picture, observed for k ∼ k H .This implies the possibility of the existence of another propagating gap for the heat waves with larger k-values at higher temperatures.This thesis is partially supported by the results for the spectra of collective excitations (see figure 1b), where the branch z 1 (k) displays a rapid decay just in the region of k ∼ 2.7 Å−1 .

Origin of relaxing kinetic mode d 2 (k)
It was shown in previous subsection, that the kinetic relaxing mode d 2 (k) causes the leading contribution to the shape of density-density time correlation function for k ≃ Q p .This relaxing kinetic mode is well reproduced within the viscoelastic subset of dynamical variables Hence, its origin is mainly viscoelastic one.This is in contrast with the hydrodynamic behaviour (1), where the relaxing contribution to density-density time correlation function is only due to thermodiffusive processes.To study the origin of relaxing mode d 2 (k) more in detail let us consider the simplest case of dynamics when only one dynamic variable A (1) = {n(k)}, namely, the particles' density n(k) is taken into account.This immediately gives the solution for the only collective mode d 0 (k), which has a very simple form: The quantity τ nn (k) is the generalized correlation time, associated with the densitydensity time correlation function: In figure 4 we show the results for two eigenvalues d 0 (k) and d 2 (k), obtained within the one-(line with triangles) and nine-variable (closed boxes) approximations of the GCM approach for liquid Pb at T l = 623 K and T h = 1170 K.In the lowtemperature state a quantitative agreement between two eigenvalues for k ∼ Q p is a striking feature.This implies that at low temperatures the slow density fluctuations for k-values being in the region of the main peak's position (that corresponds in fact to the short-range fluctuations with the average interparticle distance a , Q p ≈ 2π/ a ) are well separated in time (see figure 1) from the thermal processes as well as from the fast density fluctuations, which are responsible for sound propagation in this region, and determine almost completely the shape of the density-density time correlation function.For the higher temperature the quantitative agreement between the results, found for d 2 (k) and d 0 (k), is not so perfect due to the stronger coupling with the thermal fluctuations.However, the contribution from the mode d 2 (k) to F nn (k, 0) was still dominant when k is close to Q p (see, e.g., figure 2c for the wavenumber k = 1.6809Å−1 ).In particular, such a behaviour of the kinetic relaxing mode d 2 (k) explains why the mode-coupling theory of freezing was so successful by treating just the density fluctuations nearby the region of the main peak of static structure factor [18,19] and completely ignoring the thermal properties.
Within the simplified one-variable theory one gets the expression for the function F nn (k, t) in the following single-exponential form: Thus, within such a treatment the function τ nn (k) gives in fact the estimation for some specific time of relaxation, which, in particular, for k = Q p has a meaning of the lifetime for a particle in the cage of the nearest neighbors, and, therefore, the value 1/τ nn (Q p ) allows us to estimate the cage diffusion coefficient.
Let us summarize the results obtained for the generalized relaxing kinetic mode d 2 (k): (i) The relaxational behaviour of F nn (k, t) in the region of the main peak of the static structure factor S(k) is mainly determined by a single-mode contribution associated with the relaxing kinetic mode d 2 (k), which gives, in fact, the main mechanism of de Gennes slowing the density fluctuations [20]; (ii) The kinetic relaxing mode d 2 (k) in the region of the main peak of S(k) is defined entirely by the density fluctuations, and the thermal fluctuations have no effects herein.Moreover, the eigenvalue d 2 (k) is well reproduced for k ∼ Q p even within the one-mode approximation; (iii) For wavenumbers k ∼ Q p the kinetic relaxing mode d 2 (k) determines the lifetime of a particle in the cage of the nearest neighbors and is directly connected with the mechanism of cage diffusion [21].Thus, one can call this mode a structural relaxation mode; (iv) In the hydrodynamic range the damping coefficient d 2 (k) tends to nonzero value, and its contribution to the density-density time correlation function F nn (k, t) becomes negligible.However, we point out that the role of relaxing kinetic mode d 2 (k) increases rapidly when k becomes larger.Our recent results, obtained for liquid metallic Cs [14] and a semi-metallic liquid Bi [12], show obviously, that beyond the small-k region the relaxing mode d 2 (k) is the lowest one and this mode makes the leading contribution to the shape of density-density time correlation function for intermediate and large wavenumbers in complete agreement with the results found in this study for liquid Pb.

Analytical treatment of amplitudes: three-variable theory
Our next task is to explain the behaviour of amplitudes describing the contributions from the generalized sound excitations z 2 (k) and relaxing kinetic mode d 2 (k) to the density-density time correlation function.We have derived an analytical expression for density-density time correlation function within GCM approach, based on three dynamical variables A (3) = {n(k, t), J(k, t), J(k, t)}, assuming that this basis set would be more appropriate for the description of leading relaxing and propagating processes in the region k ∼ Q p , than the two-variable 'damped harmonic oscillator' model [22].This means that the leading dynamical processes for k ∼ Q p are assumed to be correctly reflected by one relaxing d(k) ≈ d 2 (k) and two propagating z ± (k) = σ(k) ± iω(k) ≈ z 2 (k) modes.The three-term analytical GCM expression for density-density time correlation function has the following form (see (5)): with the normalized amplitudes Here, ω2 k = k 2 k B T /mS(k) is the second frequency moment of S(k, ω).It is seen that depending on the ratio between d(k) and the damping coefficient σ(k) = Re z(k) one can obtain either positive or negative amplitudes for the symmetric contribution Bnn (k).In particular, this explains why the normalized amplitude Ānn (k) can be larger than unity at k ≈ Q p , which can be seen in figures 3a and 3b.It also follows from the expression for Bnn (k), that for strongly over-damped sound excitations, when σ(k) ≫ d(k), and small values of second frequency moment ω 2 k (heavy atoms and low temperatures), the amplitude Bnn (k) may be negative.Thus, in this case one obtains the negative amplitudes of sound contribution to the density-density time correlation function or dynamic structure factor S(k, ω) within the three-variable model.
In the case of strong damping for the short-wavelength excitations a propagation gap can emerge at k ≈ Q p .Within the propagation gap, the model produces, instead of the pair of propagating modes z ± (k), two purely real eigenvalues d (−) (k) and d (+) (k), so that in this case the density-density time correlation function is represented as a sum of three relaxing contributions with the amplitudes given by the expressions: .
Thus, we can predict that at the inner gap boundary, when d (−) (k) and d (+) (k) tend to the same value σ g , two of the amplitudes diverge and one has However, these two amplitudes are not independent.The zero-th order sum rule requires that the sum of Ā− nn (k) and Ā+ nn (k) should be finite.At the outer gap boundary, when ω(k) → 0 and σ(k) → σ g in ( 12), only the asymmetric amplitude diverges, so that we find These results correctly describe the behaviour of the mode amplitudes observed within a three-variable viscoelastic approximation.In general, they reflect the main properties of amplitudes when two relaxing processes merge into a propagating one.Therefore, one can at least qualitatively understand the divergence of amplitudes from thermal excitations on the edge of propagating gap for low-frequency heat waves shown in figures 3a and 3b.

Mode contributions to the heat-heat time correlation function
Within the three-variable hydrodynamic treatment, the following analytical expression for the heat-heat time correlation function F hh (k, t) can be derived: where a is a thermodynamic parameter.One can see, that the hydrodynamic amplitudes, describing the mode contributions from the thermodiffusive mode d 1 (k) and the propagating sound excitations, are now equal to 1/γ and (1 − 1/γ), respectively.Thus, an opposite picture in the ratio of mode contributions from thermal and acoustic excitations is observed for the function F hh (k, t) in comparison with the density-density time correlation function.
Within the nine-variable approximation of GCM approach, used in our study, one has nine separated contributions to the heat-heat time correlation function, associated with generalized hydrodynamic and kinetic excitations (see (5)).The results, obtained for the leading terms contributing to F hh (k, t) at T h = 1170 K, are presented in figure 5 for three wavenumbers k.Again, as it has been before in the case of the density-density time correlation function, the total contributions from the propagating modes z 1 (k) and z 2 (k) are shown, which are simply the sums of relevant symmetric and asymmetric terms.One can see, that for the smallest wavenumber k = k min , reached in MD, (see figure 5a) the results are in good agreement with the predictions ( 14) of the hydrodynamic theory.In the hydrodynamic region the main contribution to F hh (k, t) is caused by the relaxing thermodiffusive mode d 1 (k), and the relaxing shape is modulated by the contribution from the sound excitations (the relative magnitude of this last contribution is about 20%).Note also, that the value of k = k min is slightly beyond the hydrodynamic k-region.Therefore, one can see in figure 5a a tiny short-time contribution from the kinetic relaxing mode d 3 (k), which is completely formed by the heat fluctuations.The contribution from the kinetic relaxing mode d 2 (k) is extremely small in this region of wavenumbers.
For smaller k we expect that the contribution from the kinetic mode d 3 (k) vanishes, while the amplitudes of mode contributions from the thermodiffusive mode and sound excitations tend to the values 1/γ and (γ − 1)/γ, respectively.When k increases the situation changes, and the low-frequency heat waves begin to play a more significant role.It is seen in figure 5b, that the contributions to F hh (k, t) from the structural relaxing mode d 2 (k) and the generalized sound excitations z 2 (k) are of the order of magnitude 40% and 20%, respectively, and an additional, rather strong, contribution from the low-frequency heat waves z 1 (k) appears.For a larger wavenumber, k = 1.6809Å−1 , (see figure 5c) only two contributions from the propagating modes (the kinetic low-frequency modes z 1 (k) and the generalized sound excitations z 2 (k)) almost completely determine the shape of the heat-heat time correlation function.Note that in this last case the relative magnitude of contribution from the heat waves is about 75%, so that these modes become dominant for large wavenumbers.
To complete the study of mode contributions to the function F hh (k, t), we have performed the calculations of the k-dependent normalized amplitudes, using the scheme, described above.The results are shown in figure 6 and can be summarized as follows: (i) In the hydrodynamic limit, when k → 0, all the contributions from the kinetic modes tend to zero and become negligible, while the amplitudes from the hydrodynamic modes d 1 (k) and z 2 (k) go to the nonzero values 1/γ and (γ − 1)/γ, respectively, predicted by the hydrodynamic theory; (ii) For the time correlation function F hh (k, t) the kinetic relaxing mode d 2 (k) does not play such an important role as it is observed for the density-density time correlation function F nn (k, t).The only region of wavenumbers with a sufficient contribution from d 2 (k) is in the vicinity of the edge of the propagation gap for low-frequency heat waves when k ≃ k H ; (iii) In the case of liquid Pb, the low-frequency heat waves make the leading contribution to the shape of heat-heat time correlation function F hh (k, t) beyond the propagation gap; (iv) The amplitudes of the relaxing modes d 1 (k) and d 3 (k) exhibit divergent-like behaviour (as it is observed in the case of density-density time correlation function) in the vicinity of the edge of propagation gap for the low-frequency heat waves.

Conclusions
The main results of this study are as follows: i) The kinetic relaxing mode d 2 (k) makes the leading contribution to the shape of density-density time correlation function in the region of wavenumbers close to the position k = Q p of the main peak of the static structure factor S(k).This kinetic relaxing mode is caused by slow density fluctuations and determines the lifetime of a particle in the cage of the nearest neighbors, which is connected with the effect of cage diffusion.This lifetime can approximately be estimated from the behavior of the dynamical structure factor S(k, ω) as the ratio S(k = Q p , ω = 0)/S(k = Q p ); ii) The specific feature of liquid Pb is the small value of the propagation gap, describing by k H (k H ≈ 0.2 Å−1 for T l = 623 K and k H ≈ 0.4 Å−1 for T h = 1170 K, respectively).As it has been shown, the asymmetric contribution from the lowfrequency heat waves to the shape of the density-density time correlation function (or dynamical structure factor) increases dramatically when k tends to k H from the region of large wavenumbers.In the vicinity of the edge of propagation gap for lowfrequency heat waves, this non-hydrodynamic effect, we hope, can be the subject of experimental studies; iii) In a wide range of wavenumbers, the low-frequency heat waves do not contribute sufficiently to the density-density time correlation function.However, in the case of the heat-heat time correlation function F hh (k, t) their contribution is dominant for the intermediate and large wavenumbers.This is the main reason why the oscillating behaviour can be observed in time-dependence of F hh (k, t) for larger wavenumbers k comparing with the case of density-density time correlation function F nn (k, t).This result can be especially interesting in connection with new experimental possibilities of impulsive stimulated light scattering scattering technique (see, e.g., [23]), which permits to study the density response to laser-deposited heating of a liquid; iv) The three-variable analytical model is used to explain an emergence of negative symmetric amplitudes for short-wavelength sound excitations with wavenumbers k ≈ Q p .У методі узагальнених колективних мод досліджуються розділені вклади від різного типу колективних збуджень до часових кореляційних функцій "густина-густина" та "теплова густина-теплова густина" у простих рідинах.Показано, що кінетична релаксаційна мода, обумовлена повільними флюктуаціями густини, майже повністю визначає поведінку часової кореляційної функції "густина-густина" для значень хвильового вектора в області головного піку статичного структурного фактора.Аналітична тризмінна модель використана для пояснення негативних амплітуд, що описують вклади від коротко-хвильових збуджень до динамічного структурного фактора при значеннях хвильового вектора в області головного піку статичного структурного фактора.

1 )Figure 1 .
Figure1.Spectra of collective excitations of liquid Pb at 623 K (a) and 1170 K (b), obtained for the nine-variable set A(9) (k, t).Complex and purely real eigenvalues are shown by symbols, connected by spline-interpolated solid and dashed lines, respectively.The same symbols in upper and lower frames correspond to imaginary and real parts of propagating modes, except the case of splitting of complex-conjugated eigenvalue z 1 (k) into two real ones for small k-values.The straight dash-dotted lines in upper frames show the linear hydrodynamic dispersion of sound excitations.

Figure 2 .Figure 3 .
Figure 2. Separated mode contributions to the function F nn (k, t) calculated for three k-values at T h = 1170 K.The MD-derived function and the result of GCM study are shown by solid and dashed lines, respectively.Mode contributions from the low-frequency heat waves z 1 (k), the sound excitations z 2 (k) and the kinetic relaxing mode d 2 (k) are plotted by short dash-dotted, dotted and long dash-dotted lines, respectively.For the smallest k-value the contribution from the thermodiffusive mode d 1 (k) is shown by short dash-dotted line.Time scale is τ h = 2.3935 ps.

Figure 4 .
Figure 4. Behaviour of the relaxing kinetic mode d 2 (k) (symbols) obtained in the nine-variable GCM approximation at the temperatures: (a) T l = 623 K, and (b) T h = 1170 K.The one-variable relaxing mode d 0 (k) is shown by splineinterpolated dashed line with triangles.

Figure 5 .Figure 6 .
Figure 5. Separated mode contributions to the heat-heat time correlation function F hh (k, t) for three k-values at T h = 1170 K.The results of GCM study are shown by dashed line.The mode contributions from the low-frequency heat waves z 1 (k), the sound excitations z 2 (k), and the kinetic relaxing mode d 2 (k) are plotted by short dash-dotted, dotted and long dash-dotted lines, respectively.For the smallest k-value considered, the contribution from the thermodiffusive mode d 1 (k) and the kinetic relaxing mode d 3 (k) are shown by short dash-dotted and solid lines, respectively.Time and energy scales are τ h = 2.3935 ps and ε h = k B T h .