Dynamics of ruthenium mirror under action of soft x-ray ultrashort laser pulse

We investigate behavior of the ruthenium targets under x-ray. In this paper, the two-temperature equation of state for Ru is developed. Electronic spectrum of ruthenium is calculated using density functional theory. We have defined an electron–ion coupling parameter of ruthenium and its two-temperature thermal conductivity. With this input we run our two-temperature hydrodynamic code.


Introduction
In the present paper, the oblique illumination of ruthenium (Ru) film by soft x-rays is considered. Ruthenium is interesting because it is widely used in x-ray optics as a focusing mirror or a protective pellicle of the laminates [1]. Protection helps to prevent contamination of laminates. Ruthenium is chemically very stable metal, e.g. it is resistive against oxidization. Ruthenium has high reflectivity for x-rays at grazing angles and is a promising material for free electron laser beamlines.
Properties of the refractory metals like W, Mo, Ru are purely known in a warm dense matter regime. In this regime the densities are comparable to the solid state densities while temperatures are ∼ 1 eV. For considered here cases of ultrashort durations of the laser pulse τ L ∼ 0.1 ps the two-temperature (2T) effects are significant [2][3][4][5][6]. Electron temperature T e is higher than ion temperature T i in the 2T states created by ultrashort illumination.
The straight-through description from 2T to 1T (one-temperature) warm dense matter states is necessary for quantitative calculations because pulse duration is shorter than duration t eq of the 2T to 1T transit while the electron-ion temperature equilibration time t eq is shorter than acoustic time scale t s = d T /c s defining the rate of hydrodynamics processes [2,[4][5][6]; here d T is thickness of a heat affected zone (HAZ), c s is speed of sound; 1T relates to the onetemperature situation when T e ≈ T i . The HAZ is formed during the 2T stage. This straightthrough description of thermodynamics and kinetics of 2T and 1T states is presented in the paper below. Let us say briefly about heating and coupling between thermal and mechanical effects. The heat affected zone d T is created at the 2T stage lasting during the time interval t eq because at this interval the absorbed laser energy expands fast from an attenuation depth d att with enhanced (relative to the 1T value χ 1T ) electron thermal diffusion coefficient χ 2T , see [2,[7][8][9] and [10].
Later in time the 1T stage begins, the coefficient χ significantly decreases, and the heated zone expands more slowly, thus strongly delaying in its expansion relative to the interplay of acoustic waves; it should be mentioned that at t > t s the heat expansion weakly depends on the current acoustic situation; but in situations above an ablation threshold the redistribution of thermal energy depends on nucleation and evolution of phase composition because heat conductance of the layer with two-phase mixture is low, see papers describing thermal properties of foam [11][12][13][14]. Expansion of heat causes decrease of temperature. Temperature of the heated zone gradually decreases as c(T − T 0 ) ∼ F abs / √ πχ 1T t until thickness of the heated zone (which includes HAZ formed at the 2T stage and addition formed during the 1T stage) equal to ∼ d T + 2 √ χ 1T t is less than thickness d f of a film; in the estimate T − T 0 ∼ 1/ √ χ 1T t we neglect appearance of a foamed layer; here T 0 is initial temperature of a target, F abs is absorbed energy, t is time elapsed after an ultrashort pulse. Decrease of temperature is slowed down when the heated zone achieves the contact with dielectric or semiconductor substrate if substrate thermal conductivity is low; we illuminate a film from its vacuum side; then the contact is the rear-side of a film relative to the illuminated boundary. In the paper below we consider the case of a thick film. This means that thickness of a film d f is significantly larger than thickness of the heat affected zone d T ; in the considered case d f is 50 nm, while d T is ∼ 15-30 nm in the case of ruthenium. Energetics of laser heating analyzed below is rather high, thus absorbed energy is around an ablation threshold or higher. In all metallic targets considered previously [15,16] and in the case of a ruthenium film the ablation threshold is few times higher than the melting threshold. Then thickness of a molten layer d m is of the order of thickness of the heated zone for such energetics. Increase and decrease of thickness of a molten layer d m (t) in time follows thermal history of a target. First, it increases fast (supersonically) during 2T stage and is d m (t ∼ t eq ) ∼ d T at the end of the 2T stage. Second, the molten zone d m (t) expands slowly (subsonically) at the 1T stage. In the cases considered the energetics of heating is not enough to melt the ruthenium film through, along its entire thickness d f = 50 nm. Third, the molten layer begins to shrink thanks to cooling and finally solidifies.
The freezing history and structure formation are one of final results of the paper. But to come to description of them we have to begin and pass a long way with presenting physical model of thermodynamic states and kinetic properties of hot ruthenium and with presenting of a simulation approach.
Ultrashort (τ L ∼ 0.1 ps), soft x-ray (hν = 92 eV) irradiation begins with absorption of x-ray photons mainly by 4s 2 and 4p 6 shells of ruthenium 4d 7 5s 1 shown in figure 1. These shells are approximately 43 eV (4p) and 73 eV (4s) below a Fermi level of cold Ru. Electrons ionized from these shells have energies of the order of few tens of eV above the Fermi level of cold ruthenium.
We use method of projected augmented waves (PAW) and the exchange-correlation functional in the Perdew-Burke-Ernzerhof (PBE) form [17] to describe electron behavior. In calculations of thermodynamic properties and ionic configurations we use density functional theory (DFT) code VASP [18,19]. We use the cutoff energy 400 eV, 32 empty levels per atom, and the magnitude of the error at convergence of electronic iterations is no more than a 10 −5 eV for the total energy of the whole system. We use 528 k-points located in the first Brillouin zone of the hcp primitive cell with two atoms.
Providing full potential-linear augmented plane waves (FP-LAPW) calculations with Elk package [20] we employ the following conditions for convergence: we use 21×21×21 Monkhorst-  Pack grid of wave vectors and 20 empty electron levels per atom. The multiplication of muffin-tin radius on maximum electron wavenumber is equal to 8.5.
We can neglect energy loses connected with photo ionized electrons escaping from ruthenium in vacuum. Indeed, to stop and return back a photo ionized electron at a spatial scale d esc of its flight in vacuum we need a positive surface charge σ e and corresponding number of positive holes at a surface of ruthenium σ e /e equal to σ e /e = ε esc /(e 2 d esc ), where e is charge of electron, ε esc ∼ 30 eV is initial kinetic energy of an escaping primary electron when it intersects surface of metal. Velocity of the 30 eV electrons is 3000 km/s. During the pulse τ L ∼ 0.1 ps these electrons will fly out from a surface at the distance 0.3 µm. It is plausible to take d esc ∼ 1 µm. Then, first, electron flight time out and back to surface is ∼ 1 ps and, second, this flight distance is less than an illuminated radius R L ∼ 10 µm.
In this case a negative charge of electrons in vacuum holds equal charge of the positive holes inside an irradiated spot-the double layer preventing spreading of a positive charge out from the area of the irradiated spot exists. Configuration of electric field similar to that in capacitor is formed. The capacitor plates are the electrons in vacuum and holes at the surface of a spot. Taking given above values into account we obtain number of holes σ e /e ∼ 10 12 cm −2 . While if we put energy 100 mJ/cm 2 into holes we obtain σ e /e ∼ 10 16 cm −2 . Thus very small part of this amount will prevent escaping of energy in the form of photo ionized electrons escaping into vacuum.
ln experiments [1] considered here, ruthenium (Ru) is irradiated by soft x-rays (with 92 eV photon energy, 13.5 nm wavelength) under glancing angle of incident plane wave equal to 20 • (angle between the wavevector direction and target surface). In this conditions attenuation depth is extremely small (figure 2); data for this plot are taken from the Henke's tables [21].  The electron heat conductivity of Ru is rather small while electron-ion coupling parameter is large. Therefore the heat affected zone in Ru is thin. In this paper we will focus our attention on the thermodynamic and kinetic coefficients of ruthenium. It is especially important that they are of interest to us in a two-temperature state, the appearance of which is characteristic for the action of femtosecond laser pulses onto metals.

Thermodynamic functions of ruthenium
At the normal conditions ruthenium is a transition metal with hexagonal close-packed (hcp) structure with an edge of the base and height of hexagonal cell equal to a = 2.706Å and c = 4.282Å respectively. We have calculated the electronic spectra of solid ruthenium and its internal energy and pressure in dependence on the density and electron temperature T e at zero ion temperature T i . Program VASP based on the density functional theory was used for calculations. When calculated the electronic properties of hcp ruthenium, the variation of hexagonal ratio c/a has been made (figure 3).
In figure 4 the electron density of states is shown at low electron temperature T e = 0.01 eV and high temperature T e = 6 eV. In this density of states two domains are clearly visible, corresponding to s-band (more longer) and d-band electrons. In our subsequent analysis of the kinetic characteristics of ruthenium, we will use parabolic approximations of the electronic spectra in these bands, so that the energy of electron as a function of its wave vector k in s-band is  and in d-band Here ε s0 is a bottom of s-band, ε d1 is a bottom of d-band (unlike the top of d-band ε d2 ) and m s , m d are effective masses of s-and d-electrons. When measuring from the Fermi level at T e = 0, we obtained ε s0 = −8.0 eV, ε d2 = 1.0 eV. Pressure at zero electron temperature (cold pressure) calculated using the VASP program as a function of volume is shown in figure 5 by circles. These data correspond to zero pressure equal at the density 12.26 g/cm 3 . Calculated pressure is in a good agreement with that one calculated for the compressed ruthenium in [22]. We extended our calculations to significant expansion. To satisfy the experimentally more accurate value 12.47 g/cm 3 of density corresponding to zero cold pressure and maintain the minimum pressure value, cold pressure curve can be approximated by the expression with the volume per atom v 0 and parameters A = 3.81 a.u., a = 1.5886, b = 1.3333.
Here we have introduced the reduced density x = ρ/ρ 0 with the equilibrium density at zero temperature ρ 0 = 12.47 g/cm 3 .
Respectively cold internal energy per atom is presented as Hartree units (a.u.) for A were used.  In the framework of Mie-Grüneisen approach we have for the ion internal energy per atom at the ion temperature T i Here, G(x) = d ln θ/d ln x is the Grüneisen parameter with θ(x) being the Debye temperature. Debye temperature can be presented as with s(x) being the speed of sound. Further and from (2) θ( p 0 _(GPa) Figure 5. Cold pressure as a function of density. Solid circles are results of VASP-calculations. Solid curve is a pressure according to the formula (1). Also the results of calculations [22] are shown (empty circles).
Here we designate To avoid the discontinuity of function y(x) at decreasing values of x we change it onto the close function Functions y(x) and y 0 (x) and their derivatives coincide at x = 1 and the asymptotics of two functions coincide for large x. Then Grüneisen parameter as a function of density is shown in figure 6.

Electron thermal conductivity
Because of small mobility of d-electrons we consider a contribution of only s-electrons into the electron thermal conductivity. Relaxation of electrons at thermal conductivity process is due to the electron-ion and electron-electron scattering. Then calculating partial thermal conductivity κ ei because of electron-ion relaxation and thermal conductivity κ ee due to electron-electron relaxation we can calculate resulting thermal conductivity as [23,24] First we consider ei scattering. Thermal conductivity due to it can be presented as  Figure 6. Grüneisen parameter as a function of density.
with the electron heat capacity per unit volume C e (T e , x), mean electron velocity v(T e , x) and electron mean free path λ ei (T i , x). Electron heat capacity of s-electrons per unit volume as a function of T e and x can be represented by an expression of the type C e = n e k B f 1 (k B T e /ε F (x)), where n e is the electron concentration, ε F (x) is the Fermi energy, k B is the Boltzmann constant and f 1 is a functiom of k B T e /ε F (x) only. Fermi energy is a difference between the chemical potential at T e = 0 and the bottom of s-band. Mean electron velocity v(T e , x) can be written where v F is a Fermi velocity. Calculated dependence of the Fermi energy on the compression is consistent with that of ε F ∝ x 2/3 , m s ≃ const, v F ∝ x 1/3 . So the product C e v can be reproduced as an expression Here we introduced dimensionless parameter τ = 6k B T e /(ε F0 x 2/3 ). Electron mean free path λ ei can be written as where n is the ion concentration, σ is the electron-ion scattering cross section. Wherein with u 2 0 , u 2 t being the mean squared amplitude of zero-point and thermal lattice vibrations, θ is the Debye frequency. Here we make no distinction between acoustical and optical vibrational modes at ion temperatures T i under consideration, exceeding the Debye temperature θ (hcp lattice of ruthenium has three acoustical and three optical vibrational modes). Taking into account that Substituting dependence of C e v and λ on T e , T i , x into the electron thermal conductivity due to the electron-ion scattering (4), we obtain the dependence of κ ei on these three parameters: Then introducing the electron thermal conductivity due to the electron-ion scattering κ r at normal room conditions T i = T r = 298 K, x = x r = 12.41/12.47 (12.41 g/cc and 12.47 g/cc being densities at normal conditions and respectively at zero temperature and pressure), we can write κ ei (T e , T i , x) as Now consider κ ee therm in (3). At not too high electron temperatures it can be written in the form Taking into account that in this case ν ee ∝ (ε F / )(k B T e /ε F ) 2 , we have We suppose that at arbitrary temperatures κ ee depends on k B T e /ε F in a more general form: Function Q(τ ) was calculated at x = 1 at various electron temperatures T e with the use of electron-electron collision frequencies calculated as in [23]. Results can be well described by expression where γ 1 = 2.8544, γ 2 = 2.9883, γ 3 = 2.5245, γ 4 = 31.5261, γ 5 = 8.9447. Then using (5) and (6) we obtain κ(T e , T i , x) from (3). In order to reproduce the value of thermal conductivity under normal conditions 117 W/(m K), the value of k r was taken equal to 125.5 W/(m K). Electron thermal conductivity in two-temperature situation in dependence on the electron temperature for different ion temperatures is shown in figure 7.
When calculating the equilibrium one-temperature heat conductivity at ρ = 12.41 g/cm 3 , we obtain results shown in figure 8 together with the experimental points (corresponding to the constant pressure conditions). Difference between two curves can be partially explained by the difference between the isochoric and isobaric conditions.

Electron-ion energy exchange in ruthenium
Laser radiation in metals is first absorbed by electrons. Energy, transferred from electrons to ions per unit time and unit volume, can be written aṡ HereṄ q is the rate of change of concentration of longitudinal phonons N q with momentum q and energy ω q (in this paper we consider interaction of electrons with the acoustical phonons). The rate of change of phonon concentrationṄ q can be represented in terms of spontaneous and induced phonon emission and absorption of phonons: Here V is the volume of a solid, W q δ(ε p−q + ω q − ε p ) is the probability of electron transition from momentum p to momentum p − q. Considering electron interaction with longitudinal acoustic phonons, we can write W q as with U (q) = 4πz i e 2 q 2 ξ(q) n. Here z i to be the ion charge, n to be the concentration of atoms. Matrix element of screened Coulomb electron-ion interaction U(q) depends upon the electron dielectric permittivity ξ(q). Supposing two-temperature situation, we have separately equilibrium electron (10) and phonon distribution functions with electron and phonon temperatures T e and T i . Statistical factor Here ω = ω q , ε = ε p . Then integration over p in (8) is reduced to the evaluation of integral  12 To evauate (13) we introduce polar coordinates with z axis directed along the vector q and polar angle θ. Denoting t = − cos θ and integrating over azimutal angle from 0 to 2π, we obtain dp (2π ) 3 = 2πp 2 dp dt (2π ) 3 .
The argument of δ-function in (13) can be written as when introducing the effective mass of electron m. It gives dt = m pq dβ and integral I reduces to Taking into account that within the effective mass approach pdp = mdε we obtain Integration over ε now gives We introduce longitudinal phonon dispersion law to be simply ω(q) = sq with a sound velocity s. It corresponds to Debye approach. Now we consider energy transfer from d-electrons to ions. d-band will be approximated by the parabolic dispersion law ε(p) = ε 1 + p 2 /(2m d ) with the upper boundary ε 2 . Here m d is the effective mass of d-electrons. We introduce momentum p d = 2m d (ε 2 − ε 1 ). Now we must take into account that for d-electrons p p d and also (p − q) 2 p 2 d . This inequality can be written as Depending on the point t 0 position with respect to the interval −1 t 1 two cases arise.
In the case I t 0 1. It gives In the case I because of point t 0 lies out of −1 t 1 interval, integration over t includes the hole interval −1 t 1. Then limits of integration over β are We have β 2 > 0 for all values of p. To obtain different from zero result of integration over β, we must have β 1 < 0.This inequality takes the form It gives for p: in addition to earlier obtained (18). These inequalties define two-dimensional region of integration in p-q-plane.
In the case II t 0 point lies within −1 t 1 interval: Now integration over t is restricted by the interval −1 t t 0 . It gives limits for the integration over β: Integration of δ-function leads to the condition  (20) and (21) together with the restriction p p d define the region of integration in p-q-plane. Within the Debye approach phonon momentum module is restricted by the Debye momentum q D = (6π 2 n) 1/3 . Now we can write expression forṄ q : When considering ion temperatures exceeding the Debye temperature, we can expand the terms exp( ω/[k B T i ]) and exp( ω/[k B T e ]) into series and obtain more simple expressioṅ Substituting this expression to (7), we obtaiṅ It can be written aṡ  15 where the electron-ion coupling parameter g(T e ) is Now consider expression (9) for the matrix element of electron-ion interaction. It contains the dielectric permittivity function ε(q). We calculate it within the Lindhard approach Here first term in brakets results from screening by d-electrons and the second term is due to the screening by s-electrons. First term exceeds the second one. Its value multiplied by 4π is close to the squared screening constant, obtained within the Thomas-Fermi approach, if when calculating the screening effect we consider d-electrons as independent from s-electrons and having own chemical potential. In this case for low temperatures for example and this value does not change significantly with electron temperature increase. For q inside the Debye sphere Lindhard and Thomas-Fermi dielectric permittivities are close to each other if d-electrons are screened by d-electrons. When considering the energy transfer from s-electrons to ions, we obtain the formula, analogous to (22): It can be reduced tȯ  The electron-ion coupling parameter in this case is The electron-ion coupling parameter g as a function of the electron temperature T e is shown in figure 9.
The dielectric function in matrix element of screened Coulomb electron-ion interaction was taken in Lindhard approach [25,26].

Two-temperature hydrodynamics
We use equation of state and kinetic parameters presented above to run our 2T-HD code. At the early stage (up to picosecond) thermal processes dominate. Profiles of electron temperature T e (x, t = 0) are shown in figure 10. Time is elapsed from maximum of laser pulse. Durations τ L of x-ray pulse is equal to 100 fs; this is the e-folding time: ∝ exp(−t 2 /τ 2 L ) Electron heat conduction sharply enhanced in the 2T states (see figure 7) is important from very early stage. Figure 11 presents the electron temperature profiles together with an exponential law corresponding to photon absorption. In simulations absorption depths is equal to an attenuation depth 3.5 nm for 92 eV soft x-rays (see figure 2). Electron heat conduction gradually levels shallow and deep penetrations of electromagnetic energy. Late stage evolution of a temperature field is shown in figure 12. Internal energy fluxes from absorption depth and fills all thickness of a film. In this set of simulations we neglect heat conduction in support. x (nm) (K) Figure 10. Distribution of electron temperature for x-ray pulse at the maximum of the heating pulse (half of total absorbed energy is delivered to the ruthenium film). x (nm) Figure 11. Expansion of the energy absorbed in electrons from the heat absorption depth (solid line; distribution of absorption powers is shown by dashed line) thanks to strong electron heat conduction.
Dynamic oscillations of a film significantly influence average temperature, compare the instants t = 10 and 50 ps in figure 12. x (nm)  Figure 13. Temperature relaxation after x-ray excitation.
The late (during the first acoustic oscillation of a film) stage for x-ray laser is given in figure  12. Spreading of absorbed heat along thickness of a film smeared out the initial which was defined by absorption of laser energy (see figure 11). The electron-ion coupling coefficient is rather high in Ru. Therefore electron-ion relaxation proceeds rather fast. Equilibration of electron and ion temperatures is shown in figure 13. Equilibration takes approximately one picosecond. Temperature relaxation proceeds in narrow layer in the case of oblique illumination by a x-ray laser because the attenuation depth d att is extremely thin (3.5 nm).
Absorption of energy 100 mJ/cm 2 heats a 50 nm ruthenium film above 3 kK, see figure 12. We do not include melting in our approach developed above. Dividing a latent heat of melting for ruthenium 38.59 kJ/mol to its heat capacity 24.06 J/(mol K) we obtain approximate value of the ruthenium heat of fusion in degrees 1.6 kK. Melting temperature of ruthenium is 2.607 kK. Thus ruthenium should be molten in our calculations if it is heated above approximately 4.2 kK. Threshold for total melting of a film (if it has spatially homogeneous temperature) is approximately 60 mJ/cm 2 .

Conclusion
We have calculated thermodynamic functions and kinetic coefficients of the ruthenium, required for use in the hydrodynamic code, describing the dynamics of transformation of ruthenium sample under the action of femtosecond laser pulse. Of these thermodynamic functions are pressure and internal energy and kinetic coefficients under consideration are electron thermal conductivity and electron-ion coupling parameter responsible for the electron-ion energy exchange. All characteristics are calculated in two-temperature situation when the electron temperature exceeds the ion temperature. Using obtained information about equation of state and kinetic properties we have studied thermal and hydrodynamics processes triggered by a soft x-ray pulse in a ruthenium film.