Melting of Titanium by a Shock Wave Generated by an Intense Femtosecond Laser Pulse

Laser shock peening with ultrashort laser pulses has been studied by hydrodynamic and atomistic simulations, as well as experimentally. It has been shown that, in contrast to traditional nanosecond pulses, ultrashort laser pulses allow one to increase the produced pressures by two or three orders of magnitude from 1–10 GPa to 1000 GPa (1 TPa). The physics of phenomena changes fundamentally because shock waves generating pressures exceeding the bulk modulus of a metal melt it. It has been shown for the first time that the shock melting depth at pressures about 1 TPa is an order of magnitude larger than the thickness of the melt layer caused by heat conduction. The appearance, propagation, and damping of a melting shock wave in titanium have been studied. The damping of the shock wave makes it possible to modify the surface layer, where the melting regime changes from a fast one in the shock jump to a slow propagation of the melting front in the unloading tail behind the shock wave. It has been shown experimentally that the ultrafast crystallization of the melt forms a solid layer with a structure strongly different from that before the action. The measured depth of this layer is in good agreement with the calculation.


INTRODUCTION
Laser shock peening (LSP) is an industrial technology for significant surface strengthening of goods. The corresponding orders are carried out by several scientific industrial groups in the People's Republic of China, France, Czech Republic, 1 and Ireland. 2 The complete LSP strengthening of a US military aircraft was reported 3 in January 2022. The prolongation of the exploitation of the aircraft whose cost is $100 million by a factor of 1.5-2 provides a gain of $50-100 million, which is two orders of magnitude larger than the expense of LSP.
There are two variants of LSP differing in the duration of used laser pulses. The first variant is nanosecond LSP [1][2][3][4]. The second, relatively young, variant involves ultrashort (femtosecond and picosecond) laser pulses (fs-ps LSP) [5][6][7][8][9][10]. Although nanosecond LSP is more elaborated and widely used, fs-ps LSP has its advantages. First, the protective film is not required; second, irradiation through water is not necessary; and, third, the amplitude of generated shock waves is several orders of magnitude higher [7,8,10]. Let us discuss these features.
The protective film prevents evaporation and melting. Meanwhile, the melt layer is exceptionally important for the modification of the crystal structure (see [11,Sect. 4.2]). The role of the liquid phase on the surface is discussed in this work. The inertia of water increases the momentum transferred from the nanosecond pulse to the target [1,12,13]. However, the presence of a liquid does not affect the amplitude of the shock wave (SW) in the case of ultrashort laser pulses because the ultrashort triangular compression wave behind the shock jump is a simple Riemann wave. The base of this triangular wave increases as 1 www.beamtech-laser.com, www.alphanov.com, www.hilase.cz. 2 www.lsptechnologies.com. 3 www.naval-technology.com/news/frce-verification-milestonepeened-f-35b. [14], where x SW is the path covered by the SW, and a signal from the liquid cannot reach the shock jump until this jump has the amplitude necessary for laser shock peening. Thus, the amplitude and shape of the triangular wave are the same in the presence and absence of the liquid at the same fluence of ultrashort laser pulses F abs absorbed in the target.

OPTICS AND LASER PHYSICS
The possibility of the sharp increase in the pressure is the most important property of ultrashort laser pulses significantly distinguishing them from nanosecond pulses. The presence of water, the filamentation of the laser beam in water, and the optical breakdown of water strongly limit the intensity of nanosecond pulses incident on the surface of the target. As a result, induced pressures are no more than a few gigapascal. At the same time, a beam of ultrashort laser pulses passing through air or vacuum easily provide pressures of tens of terapascal because of a high optical strength of air and high fluences F abs up to hundreds of joule per centimeter squared [8,10].
The elementary estimate of the pressure is p F abs /d T , where d T is the heating depth. Heating by ultrashort laser pulses is studied in detail in [11]. The heating depth is d T = δ + d EHC , where δ is the thickness of the skin layer, which is 10-20 nm in metals irradiated by optical radiation, and d EHC is the thickness of the layer heated by a supersonic electron thermal wave, which propagates at a supersonic velocity at the two-temperature stage (2T, ) [15]. The thickness d EHC in well-conducting metals (e.g., gold) is a factor of 5-10 larger than the thickness of the skin layer δ. In this work, we study the practically important case of titanium, which, as well as zirconium, is a poorly conducting metal. The electrical and thermal conductivities of titanium at room temperature are 1/19 and 1/15 of the respective values for gold. Our 2T calculations for titanium show that the heating depth is about the thickness of the skin layer δ; the detailed description of the 2T model is beyond the scope of this work. In the 2T hydrodynamic (2T-HD) numerical simulation, the 2T model was used with a wide-range multiphase equation of state of titanium [16,17]. The incident fluence F inc in our experiment was 13.8 J/cm 2 .
We used a terawatt laser based on a chromiumdoped forsterite crystal [18]; the wavelength is nm and the FWHM pulse duration is τ = 110 fs. The absorption coefficient of titanium at this wavelength is = 0.4 (see https://refractiveindex.info/). The absorbed fluence is F abs = AF inc = 5.5 J/cm 2 . The initial pressure estimated by F abs /d T with d T = 30 nm is 2 TPa (20 Mbar). The 2T calculation shown in Fig. 1 gives slightly lower pressures at the 2T stage partially because the electronic Grüneisen parameter is low Furthermore, the initial pressure is approximately halved because of the generation of acoustic waves expanding from the high-pressure layer with the thickness d T with initially immobile matter. The 2T stage ends at 4-5 ps; i.e., the electron pressure becomes low.
Melting plays a decisive role in many laser applications. However, melting caused by the energy absorption and heat transfer through heat conduction is usually analyzed. In this work, we study the main effects caused by melting induced by the dissipation of the kinetic energy in the shock front. Melting SWs were studied in [16,[19][20][21], but strengthening processes were not analyzed. Meanwhile, as shown below, they are important in this problem.
First, melting and subsequent solidification of the surface layer fundamentally changes the crystal structure of this layer. In particular, the ultrafast quenching of the melted surface layer in most of the known cases results in the grinding of the grain structure down to the formation of a nanocrystal state in the case under consideration. As known, this state not only ensures significant strengthening but also allows one to keep the necessary plasticity or even to increase it. This is a fundamentally important achievement for applications because this makes it possible to reach record strengthening in the surface layer of metals.
Second, the thickness of the layer d EHC heated owing to heat conduction in poorly conducting metals (titanium, zirconium) heated by ultrashort laser pulses is small. Therefore, the thickness of the melt layer caused by this effect is also small. Below, we demonstrate that the thickness of the shock-melted layer at high absorbed fluences is an order of magnitude larger. Here, is the ion pressure; see, e.g., [11]. The bulk modulus of titanium is B = 110 GPa. The bulk speed of sound in titanium at room temperature is 4.9 km/s. Third, melting in the SW usually occurs at pressures behind the shock front higher than the shear modulus G, about the bulk modulus B.
Fourth, the amplitude of the SW penetrating deep into the target decreases quite slowly, . Writing , we obtain Thus, the thickness of the melt layer d m can be much larger than the thickness of the layer melted owing to heat conduction (the thickness of this layer does not exceed d T ). The aim of this work is to theoretically calculate the melting depth d m and to compare it to the value measured in our experiment.
We performed two types of calculations with (i) a 2T-HD code and (ii) a hybrid code, where the 2T stage is described by the 2T-HD code and the single-temperature (1T) stage is simulated by the molecular dynamics code. At the 1T stage, the electron and ion temperatures are close to each other, and the electron pressure in the condensed phase at temperatures about the melting temperature T m is low; the reference value for titanium is T m (p = 0) = 1.941 kK.
This work consists of two main parts. In the first part, the results obtained with the hybrid code are described. In the second part, we present our experimental data and compare them to the numerical simulation.

NUMERICAL SIMULATION
Both used hydrodynamic (2T-HD) and molecular dynamics (MD) codes have their advantages and disadvantages. The disadvantages can be reduced by combining these codes in the hybrid 2T-HD/MD code in application to the complex physical problem under consideration. It is necessary to consistently describe the spatiotemporal evolution of coupled processes from the formation of high-temperature states (temperatures about 10 eV) in the heated layer to cold elastoplastic processes deep in the target at depths The error of the MD approximation increases at high temperatures, whereas the plastic 2T-HD code disregards elastic effects in the cold solid phase.
We consider processes in a solid target irradiated by an intense (F abs is two orders of magnitude above the ablation threshold) subpicosecond ps) laser pulse. The absorbed fluences F abs of interest are such that the solid substance in the heated layer with a thickness d T of tens of nanometers is transformed to the 2T state with electron temperatures eV. The hydrodynamics of such 2T flows is not described by the MD code, which involves the interatomic potential. Since this potential is chosen according to data on the solid phase [22], it poorly approximates states where the electron contribution is significant. In our case, these are states with densities from the critical to solid-state value and with temperatures one or two orders of magnitude higher than the critical temperature T cr , i.e., states far above the binodal on the (ρ, T) plane; the binodal is the line of equilibrium coexistence of the gas and condensed phases. At the same time, dense media (crystals, liquids) and gases are described by the potential well up to the temperature T cr [23][24][25].
The hybrid 2T-HD/MD code is based on the spatial separation of the HD and MD regions. Time separation was previously used: the initial 2T stage of the process at t < t eq was simulated in the 2T-HD code. The temperature profile after electron-ion relaxation, i.e., at t > t eq , was transferred to the MD code [26]. However, this approach is applicable at relatively low fluences F abs about the ablation threshold; F abs|abl 0 .1 J/cm 2 for bulk metal targets. Change in the density at such fluences to the time t ~ t eq can be neglected in the first approximation. For this reason, the ion temperature profile T i (x, t ≈ t eq ) is transferred to the MD code and the density of the unperturbed crystal is accepted for the density.
In our case with a very large ratio F abs /F abs|abl , the heated layer d T is thin (titanium with a low thermal conductivity), heating is strong, and the speed of sound in the hot layer d T is increased because of heating. Consequently, the expansion rate is high (see Fig. 2), and the acoustic time scale t s = d T /c s (6 ps for d T = 30 nm and c s = 5 km/s) becomes comparable with the 2T relaxation time in titanium t eq = 3-5 ps. Correspondingly, the approach with the time separation of the codes is no longer applicable.
Indeed, it is difficult to transfer the instantaneous (at the HD → MD switching time ) density profile from the HD to MD code. Figure 2 demonstrates how the density profile differs from an unperturbed constant value of 4.506 g/cm 3 . The MD code starts with initial data in the form of crystal. The lattice parameter is chosen such that the crystal before loading (before coupling to the thermostat) is in the unloaded state; i.e., the stress tensor is zero. At the same time, the temperature profile is transferred easily because the Langevin thermostat approaches the temperature distribution in the MD sample to a given distribution [22,26]. This temperature distribution is established in the thermostat operating time t therm , which is chosen much smaller than t s .
Thus, the approach with time separation is inapplicable at high fluences. For this work, we developed the hybrid 2T-HD/MD code with spatial separation, where the hot left and cold right regions in Figs. 2 and 3 are described with the 2T-HD and MD codes. The main points of this simulation are as follows.
Black circles in Fig. 3 are Lagrangian cells in the 2T-HD calculation. The system of Lagrangian 2T hydrodynamic equations is presented in [11]. The solution of the 2T-HD equations provides the trajectory and velocity dx p (t)/dt (see Fig. 4) of the Lagrangian particle in Fig. 3. This particle is chosen such that it is located certainly to the right of the heated layer . Titanium in this particle (and, thereby, to the right along the x axis) is in the 1T state at a temperature remaining below the critical temperature T cr after the passage of the SW. Therefore, the approximation of the interatomic potential is applicable, and the MD code operates with a high accuracy. The potential of the embedded atom model for titanium was used in our MD calculations [27].
The MD code in the hybrid approach is used in the variable-length segment The thickness of the MD layer was chosen sufficiently large L x = 1000 nm for the unloading wave appearing after the reflection of the SW from the free right boundary to reach the melting front at sufficiently long times. The boundary condition of the equality of the matter velocity in the molecular dynamics calculation to the velocity dx p (t)/dt shown in Fig. 4 was imposed on the left boundary x p (t). Figure 4 demonstrates the passage of the shock wave in the 2T-HD calculation through the Lagrangian cell = 300 nm. The shock wave is formed from a smooth (i.e., without jump) compression wave in Figs. 1 and 2. The compression wave overturns because of a nonlinear effect (focusing of characteristics); i.e., the jump appears (see, e.g., [26]). The remnant of the smooth compression wave travels ahead of the jump. The amplitude of the jump increases rapidly with time because of the absorption of the front and back remnants of the smooth compression wave. The jump moves at a supersonic velocity, and the remnant travels ahead at the speed of sound but has a certain initial margin in the x position of the characteristics. Only a small part of the remnant is ahead of the SW in Fig. 4

because the chosen separation/piston position in Figs. 3 and 4 significantly exceeds the heating depth
A rarefaction wave, where the velocity decreases smoothly, follows the jump (see Fig. 4). Fig. 4. The 2T-HD code describes the entire field in the x coordinate from to ∞, whereas the MD flow is bounded from the left by the piston at = 300 nm. The superscript 0 marks the coordinate of the Lagrangian particle x 0 on the Lagrangian coordinate axis. The Lagrangian mark is independent of the time because the coordinates are integrals of motion. Figure 6 shows the HD and MD stresses in the process of propagation of the laserinduced shock wave in titanium. Titanium in the MD simulation is more rigid; for this reason, the MD stress is higher, and the MD shock front is more right. The damping rate of the SW in the MD simulation is slightly higher than that in the HD calculations-cf. profiles in Fig. 6 at a time of 48 ps.  and are the Euler (mobile) and Lagrange (stationary) coordinates of the piston.

Figures 5 and 6 present the 2T-HD and MD calculations with the piston shown in
The temperatures calculated in the HD and MD codes are presented in Fig. 7, where it is seen that MD temperatures are slightly higher possibly since the heat capacity is slightly underestimated because of the incomplete inclusion of electron degrees of freedom at temperatures of 4-5 kK. The reference melting tem-perature of titanium at low pressures is T m (p = 0) = 1.941 kK. The melting temperature in the MD calculation is related to the interatomic potential and is T m (p = 0) = 1.59 kK. The effect of the boundary condition on the piston is manifested in the region between vertical straight segments 2 and 3 in Fig. 7. = 300 nm of the Lagrangian particle. The distance between vertical straight segments 1 and 2 is the displacement of the piston in a time of 36 ps. As seen, the hydrodynamic and molecular dynamics velocity fields on the piston are in good agreement. The shock front in the molecular dynamics calculation is slightly further because titanium in the molecular dynamics calculation is slightly more rigid. shock wave in titanium. The 2T hydrodynamic code covers the entire flow from the left to the right. The molecular dynamics flow is to the right of the piston, which pushes titanium in the molecular dynamics calculation. The unloaded material is to the left of the coordinate 130 nm, where the pressure gradient is small. Consequently, unloaded Lagrangian particles move by inertia, i.e., with the conservation of the velocity in the Lagrangian particle (see Fig. 5). ≈ x Fig. 7. (Color online) Temperature distribution. It is seen that the region covered by the molecular dynamics calculation is expanded far beyond the hot region with high temperatures. Vertical straight segments 1 and 2 are the same as in Fig. 5 and mark the path of the piston. Vertical straight segments 2-5 in the molecular dynamics calculation are described in the main text.
The piston in the MD simulation is described by a steep repulsive potential. The decrease in the temperature between vertical straight segments 3 and 4 is due to the entropy distribution. Lagrangian particles near vertical straight segment 3 passed through the more intense SW and, therefore, have a higher entropy. This effect exceeds the opposite effect caused by a decrease in the pressure in the rarefication wave in the direction from vertical straight segment 4 to vertical straight segment 3 (see Fig. 6 with the pressure profile). The structure 4-5 in Fig. 7 corresponds to the SW. A nucleus of the elastic precursor is seen.

ON THE CALCULATED THICKNESS OF THE MELT LAYER
The transformation of the laser-induced SW from fast melting in the shock wave to slow melting in the unloading tail behind the SW is illustrated in Fig. 8. The SW melts titanium to about 40 ps, resulting in a sharp transition from the solid to liquid phase immediately behind the shock front. The structure of the shock front in this regime is shown by vertical straight segments 4 and 5 in Fig. 7 (see also the 43.2-ps profile in Fig. 8). Melting and establishment of equilibrium temperature of the melt occur in the region between vertical straight segments 4 and 5 in Fig. 7. Melting is manifested in a decrease in the local atomic order parameter below 0.43. The profile of the parameter has a deep dip in the melt region, as seen in Fig. 8. The definition of and its values for various crystal lattices are given in [28].
The shock and melting fronts are separated after 40 ps. Then, a layer with a width up to 50 nm is formed, where the parameter increases from the value for the liquid phase to the value for the solid phase. The parameter in plastically compressed solid titanium is larger than that in cold (at room temperature) unperturbed titanium ahead of the SW. The reason is that initial uncompressed titanium has a hexagonal close packed lattice corresponding to the minimum parameter among all lattices; hence, averaged values in plastically deformed titanium having numerous defects of the initial lattice are larger.
Thus, the plastically deformed titanium layer after 50 ps shown in Fig. 8 is bounded from two sides by wide transient regions. The right bound is the SW and the left bound is the structure between the liquid and solid phases of titanium. The transient structure at the direct melting stage in the SW is thin (see the 43.2-ps profile in Fig. 8). It is noteworthy that the formation of the metastable supercooled liquid because of mechanical melting in the plastic shock front [20] with subsequent crystallization lasts a certain time, which is manifested in the sharp decrease in the parameter and its recovery immediately behind the plastic shock front (see 52.8-81.6 ps profiles in Fig. 8).
Further, the SW penetrates deep in solid titanium, the melt layer slowly expands, and the dip of the parameter becomes deeper already inside the layer covered by the SW (see the profiles in Fig. 8). The expansion of the melt layer and the deepening of this dip are due primarily to a decrease in the pressure in the rarefaction wave. Correspondingly, the melting temperature T m (p), which monotonically depends on the pressure, decreases. As seen, the melting layer thickness, i.e., the maximum penetration depth of the melting front into titanium, is approximately 550 nm.

EXPERIMENT
The laser radiation parameters (absorbed fluence, wavelength, and pulse duration) are presented in the Introduction. The laser beam was obliquely incident on the surface and the illuminated spot on the target surface had the shape of an ellipse with semiaxes of 35 and 42 μm at the 1/е level, which are much larger than the propagation depth of the SW in the bulk of titanium (see Fig. 8). Consequently, the one-dimensional approximation is applicable in calculations. Figure 9 shows electron microscopy images of the microstructure of the thin foil from the cross cut of the surface layer (lamella) of technically pure ВT1-0 titanium obtained in a Tecnai Osiris transmitting electron microscope at an accelerating voltage of 200 kV. The method of preparation of the lamella is described in [29] and includes the preparation of the cross cut by a focused ion beam in the column of the electron micro- scope and its subsequent thinning to a thin foil (lamella). The protective platinum coating is preliminarily deposited (Platinum coating in Fig. 9).
The dark-field image ( Fig. 9) clearly demonstrates nanocrystal layer B with a thickness of 300-600 nm, which significantly differs in microstructure from the initial recrystallized state marked by letter A. The initial structure consisted of large 35-μm crystallites. The part of the substance shown in Fig. 9 was a single crystal before the laser impact. Layer B has an inhomogeneous plate nanocrystalline grain-subgrain structure with crystallites 5-200 nm wide and plates with lengths from 100 to 600 nm in the direction perpendicular to the sample surface (see thin arrows from circle С in Fig. 9). Crystallites marked by arrows from circle С whose plates in the reflecting position are oriented toward the transmitting electron beam have light contrast; otherwise, dark. Correspondingly, light and dark strips alternate in layer В in Fig. 9.
Crystallites, which are elements of the grain-subgrain structure, in layer B are elongated in the direction perpendicular to the sample surface. Heat dissipation along the temperature gradient is maximal in this direction. The white dashed line marks the boundary between the main material А from the surface layer B with the modified structure and is a large-angle misorientation boundary (grain boundary). Electron microdiffraction data recorded from regions А (crystal structure of the initial material with an average grain size of about 35 μm) and B (with nanocrystalline structure) belong to different crystallographic zones.
The analysis of the electron diffraction patterns of region B shows that mutual misorientations at intergrain boundaries of nanograins can be both small-and large-angle. Consequently, it is appropriate to suggest that this structure is formed after the nanocrystallization of the laser-melted surface layer of titanium with the formation of nanograins elongated along the maximum heat dissipation direction.
The thickness d m of the surface layer melted and then crystallized determined according to the signatures indicated above is in good quantitative agreement with the calculated value.

CONCLUSIONS
Melting and solidification are important because they sharply change the initial polycrystal structure. Melting by a strong SW has been considered likely for the first time in this work in application to the technology of ultrashort laser shock peening. Melting ends when the amplitude of the shock wave decreases to about 1 Mbar.
The melting layer at a high initial pressure is sufficiently extended. Its thickness in the case under consideration is an order of magnitude larger than the thickness of the layer heated by heat conduction.
The thickness of the melting layer numerically calculated is in good agreement with the experimental data obtained in this work.

CONFLICT OF INTEREST
The authors declare that they have no conflicts of interest.