Abstract
High fluence femtosecond laser percussion drilling has a potential to produce micro holes for cooling purpose in turbine blades made of titanium alloys due to its short interaction time with the material. For an in-depth investigation, numerical simulation is carried out using 2D axi-symmetric two-temperature and heat conduction models in tandem, mimicking the laser percussion drilling. Moving mesh approach in finite element based COMSOL Multiphysics software is used to arrive at the crater geometry during the ablation process. The heat conduction model provides the value of surface temperature after each pulse. Taking temperature-dependent thermo-physical and optical properties with intraband absorption effect, the simulations are carried out on a 2 mm thick Ti6Al4V plate for 1–15 pulses with fluence in the range of 0.84 to 8.4 J cm−2. For comparison, the simulation results based on room temperature properties are also included. Validation experiments are carried out and the crater morphology is measured using laser scanning confocal microscope. The overall average absolute deviations of the results for crater diameter and depth are within 20% and 40% respectively. The proposed simulation approach is robust and can be used to investigate multi pulses laser ablation process in any other applications.
Keywords
Introduction
Lasers are widely used in the manufacturing applications such as machining, welding, and additive manufacturing. The selection of laser is playing a key role to obtain the quality and repeatability, as the material properties like absorptivity varies for different laser wavelength from metals to non-metals.1,2 The continuous wave lasers are generally used for welding and additive manufacturing applications.3,4 The pulsed lasers are used for laser cutting and drilling applications. 5 The ultrashort pulse lasers with a pulse duration of less than 10 ps minimizes thermal damages due to heat affected zone and recast layer. 6 Therefore, ablation process with femtosecond laser plays a vital role in producing microfeatures with high accuracy and repeatability. A key advantage of the femtosecond laser is its ability to produce peak intensity in the range of 1014 to 1017 W cm−2 in a short duration, which is less than material thermalization time of 5–10 ps. When the surface temperature reaches the evaporation temperature during the laser-material interaction, the material removal takes place through direct vaporization which is termed as laser ablation. Laser-material interaction has been investigated for gold, silver, aluminum, copper, polycrystalline diamond, molybdenum, etc.7–12 Surface integrity and microstructural changes in metals have also been examined by researchers13–15 A very limited work has been reported on Ti6Al4V which finds a wide range of applications in aerospace, automotive, and medical industries. Micro holes for cooling of titanium alloy turbine blades are laser drilled. Experimental results on laser machining of titanium and its alloys have been reported by a few researchers.16–19
A few attempts have been made to bring out the process mechanism of femtosecond laser-matter interaction on titanium alloy using theoretical models and single pulse ablation experiments. In most of the theoretical investigations, material thickness of generally less than 100 µm is considered. In one research work, a 2 mm thick Ti6Al4V plate has been considered and an axi-symmetric finite element model with constant electron-phonon coupling is used to simulate the single pulse laser ablation. The crater geometry is obtained using element removal method. 20 Measurement of ablated crater geometry also poses challenges, as the dimensions are in nano and micro scale ranges. By using laser scanning confocal microscope, the diameter and depth of the crater could be measured reliably by them.
The thermo-physical properties of the materials change with temperature during laser heating, and hence these effects have to be considered in the simulation. The optical properties of metals under laser heating also undergo change with temperature and these properties can be computed based on Drude’s formula. 21 In multi pulses laser ablation, the absorption of laser energy and the crater geometry vary significantly for subsequent pulses. 22 A few researchers have investigated the incubation effect in multi pulses laser ablation of titanium and its alloys at high fluence regime.23–25 To the best of authors’ knowledge, no attempt has been made to investigate the influence of temperature-dependent thermo-physical and optical properties of Ti6Al4V alloy on single and multi pulses femtosecond laser ablation at high fluence.
In the present work, the mechanism of femtosecond laser ablation in high fluence regime (F > 0.7 J cm−2) has been investigated on the basis of simulation and experimental results. Since all the temperature-dependent thermo-physical properties for Ti6Al4V alloy are not available in the literature, the properties such as electron heat capacity, electron thermal conductivity, electron-electron, and electron-phonon scattering coefficients are derived from the fundamental physics. Also, temperature-dependent optical properties of Ti6Al4V alloy are computed from Drude’s formula to investigate the intraband absorption effects. The 2D axi-symmetric two-temperature model (TTM) is numerically solved by using finite element method in COMSOL Multiphysics software. The ablated crater morphology is predicted based on the moving mesh approach available in the software. The surface temperature after vaporization of material is predicted based on 2D axi-symmetric heat conduction model (HCM) by considering temporal and spatial coordinates. The morphology of laser-ablated crater profiles is evaluated by using laser scanning confocal microscope (LSCM). The destructive technique confirms the suitability of LSCM to measure crater depth. The crater profiles from experimental and simulation results for single and multi pulses ablation are compared and discussed. The details of experimental work, 2D axi-symmetric TTM and HCM, material removal, and simulation results are presented and discussed in this paper.
Laser-material interaction
The laser-matter interaction involves two processes, wherein the photon energy is absorbed by free electrons and heated up first, and then heat energy is transferred to the lattice system through the electron-phonon coupling. 26
The following are the assumptions in the present study.
(i) The laser source is treated as Gaussian laser pulse having full width half maximum intensity.
(ii) The electron and lattice temperature systems are influenced by laser properties and thermo-physical and optical properties of the material.
(iii) The temperature dependent properties play a major role at high laser fluence.
(iv) TTM captures changes in the electron and lattice temperature during the laser-material interaction. While solving TTM for systems of temperature only, the metal removal is not considered.
(v) Though the crater geometry is 3D in nature, 2D axi-symmetric simulation captures the crater geometry in the axial plane without any loss of information.
(vi) Since the crater geometry changes with time, moving mesh method is appropriate to predict the metal removal which is controlled by the evaporation temperature of the metal.
(vii) Once the crater formation is completed, the heat conduction model (HCM) can simulate the lattice cooling with less computational effort.
2D Axi-symmetric TTM
In order to investigate the femtosecond laser ablation process, it is necessary to understand the laser-matter interaction and thermal process leading to metal removal. The TTM is described with two coupled differential equations for solving electron (Te) and lattice (Ti) systems of temperature in space and time. The governing equations of the 2D axi-symmetric TTM are given in generalized form as equations 1(a) and 1(b). First term on RHS of equation 1(a) represents heat flux in electron system, while second term signifies electron-phonon coupling. The third is laser source term. Equation 1(b) contains heat flux in the lattice system as the first term.
where Ce – electron heat capacity, ke – electron thermal conductivity, G – electron-phonon coupling factor, Ci – heat capacity of lattice and ki – thermal conductivity of lattice represent constant temperature (300 K) properties. When temperature-dependent properties are considered, the relevant terms become Ce(Te), ke(Te, Ti), G(Te, Ti) and so on, where Te is electron temperature and Ti is lattice temperature.
Considering Gaussian laser pulse having full width half maximum intensity, pulse duration of tp (s) and laser fluence F (J cm−2), the laser source Q is given by equation 1(c). The Gaussian distribution of the laser pulse is described with the analytical function of r and t in spatial and temporal coordinates respectively.
where z – axial coordinate, R – reflectivity, r(z) – laser spot radius as a function of axial coordinate, and zeff – effective penetration depth as the sum of thermal (zt) and optical (zo) penetration depths. The Beer-Lambert law is used to describe the absorption of laser energy into the metal surface as a function of z for multi-pulses laser ablation, as local target surface changes after each pulse. 27 In percussion drilling, the laser beam propagates with spot radius as a function of axial coordinate z depending on beam quality factor M2 as given by equation 1(d).28,29
where ro – laser spot radius at focus point (µm), ZR – Rayleigh range (µm) and M2 – beam quality factor (1.45).
Thermo-physical and optical properties
The thermo-physical and optical properties of Ti6Al4V alloy at constant temperature of 300 K are given in Table 1. Whenever certain values are not available for Ti6Al4V alloy, the relevant properties are taken for pure titanium.
Properties of Ti6Al4V.
Pure Titanium.
Thermo-physical properties
where Be– electron heat capacity coefficient (314 J m−3 K−2). 41 The room temperature value Ce is 9.42 (104) J m−3 K−1 as given in Table 1(a).
where νF– Fermi velocity (1.79 (106) m s−1). 42 Temperature dependent electron relaxation time τe(Te, Ti) is given by equation (4). 31
where Ae – electron-electron scattering coefficient in s−1 K−2 and Bi – electron-phonon scattering coefficient in s−1 K−1 are given by equations (4b) and (4c). 32
where EF – Fermi energy (8.85 eV), 41 Be – electron heat capacity (314 J m3 K−2) 41 and We-e – thermal resistivity (5 (10−4) m K W−1) of titanium. 40 m* is effective electron mass.
where ρer – electrical resistivity (1.78 (10−6) Ω-m), 37 q – electron charge (1.60 (10−19) C) 38 and ne – electron number density (5.57 (1028) m−3). 31 Since the effective electron mass m* for use in equations (4b) and (4c) is not available for titanium, it is computed as 1.14 (10−30) kg from equation (4d). 43
where ωp – plasma frequency (4 (1015) Hz), 45 nc – conduction electron number density (0.63 (1028) m−3) 33 and ε0 – permittivity in free space (8.85 (10−12) C 2 N−1m−2). 38
The electron-electron scattering coefficient Ae and electron-phonon scattering coefficient Bi are 4.42 (108) s−1 K−2 and 7.57 (1012) s−1 K−1 respectively (Table 1(a)). The room temperature electron relaxation time τe of 0.4 (10−15) s is computed using equation (4a) with Te = Ti = 300 K, and this value is used in equation (3) to compute electron thermal conductivity ke as 43.5 W m−1 K−1.
where G- electron-phonon coupling at room temperature is 40 W m−3 K−1. 34
where τp-p – phonon-phonon scattering time (s) and vi – phonon velocity (ms−1) are computed from equations (7a) and (7b). 36
where ar – atomic radius (0.147 (10−9) m),
38
ni – ion number density (m−3) is calculated from relation ni = Zne taking averaged charge state Z as 2.
39
The phonon velocity vi(
where kB – Boltzmann’s constant (1.38 (10−23) J K−1), 38 vs – speed of sound (4.987 (103) m s−1) 37 and mi – ion mass (1.672 (10−27) kg). 38
where τep(Te, Ti ) – electron-phonon coupling time (s) is given by equation (9). 34
The thermal penetration depth zt at room temperature is 73.3 nm.
Optical properties
The energy absorption in the laser-matter interaction depends on the optical properties of material such as reflectivity and absorption coefficient. In present work, the optical properties are computed from Drude’s formula.
The temperature- dependent refractive index and extinction coefficients are calculated from equations (11a) and (11b). 21
The optical components such as real and imaginary are calculated form equation (12) based on Drude’s theory. 43
where ωp – plasma frequency (4(1015) Hz) 41 and ω – laser frequency (Hz) from equation (13a).
where c – speed of light (3(108) m s−1) 38 and λ – wave length (800 (10−9) m).
The electron collision frequency is expressed as inverse of election relaxation time as given in equation (13b). 43
Taking the refractive index n and the extinction coefficient κ as 1.37 and 2.38 respectively, the room temperature R is computed as 0.52 (Table 1(b)).
where α (Te, Ti) – absorption coefficient (m−1), κ(Te, Ti) – temperature dependent extinction coefficient and λ – wave length (800 (10−9) m). The room temperature value of α is 37.38 (106) m−1 as given in Table 1(b).
The room temperature value zo is 26.74 (10−9) m.
Numerical simulations
In the present work, the 2D axi-symmetric TTM governing equations are numerically solved in space and time domains using the finite element module in COMSOL Multiphysics software. The axi-symmetric two-dimensional geometry with 1000 µm in the radial direction and 2000 µm height in the axial direction is considered in rz-plane, as shown in Figure 1. In order to have proper comparison with the experimental work, the present simulation is carried out with a laser spot radius at focus point ro (µm) as estimated from equation (16).20,34
where f – focal length and do – initial diameter of laser beam. In the present work, ro is calculated as 6.15 (10−6) mm by considering do = 12 mm, M 2 = 1.45, f = 100 mm and λ = 800 (10−6) mm.

2D axi-symmetric geometry in rz-coordinates.
In the present work, a 2D quadratic element order is selected with Lagrangian shape function type. The mesh setting sequence type is selected as user-controlled mesh option. The 2D axi-symmetric geometry is discretized with a fine mapped meshing. In femtosecond laser heating process, the penetration depth is in the range of few nanometers hence a finer mesh is needed at the top surface of geometry. Therefore, in both r and z directions, a mesh refinement is applied by using predefined distribution mesh. The number of elements in radial direction is defined as 200 with element ratio of 4000, while it is 400 with element ratio of 4000 in axial direction. The element growth rate is defined by selecting the geometrical sequence option for exponential element distribution. The predefined element size is calibrated with plasma module available in COMSOL. Convergence tests carried out also confirm the selection of element sizes, resulting in a total number of 80,000 quad elements.
In COMSOL Multiphysics software, appropriate discretization scheme is applied to the time-dependent partial differential equations (PDEs) to obtain system of ordinary differential equation (ODEs) which are then solved by suitable iterative method. 46 In view of proprietary nature of the software, it is difficult to know the complete details of the algorithms used. However, the user is provided with certain options. One such option is to select a value for relative tolerance, whose default value is 0.01. With a view to capture the pulse on time step of 100 fs at each pulse, a relative tolerance value of 10−4 is selected in present simulation work.
The coupled PDE governing equations are solved using the COMSOL. After checking the preliminary results for 10 fs time-step, an initial time-step of 100 fs followed by the adaptive time steps up to 0.1 ms is chosen. In the present simulation work, the temperature-dependent material thermo-physical and optical properties based on Drude’s formula are combined into a single system of dependent variables and are made available during each step during the simulation. For comparison, the results are obtained with constant temperature properties also.
Single pulse laser ablation
The numerical simulations are carried out with laser fluence in the range of 0.84 to 8.4 J cm−2. The temporal distributions of electron and lattice temperature for single-pulse laser ablation are predicted. When the lattice temperature reaches the evaporation temperature, the material directly evaporates forming a gaseous phase. As a result of the material removal from the surface, there is a change of geometry. The instantaneous geometry after each time step of FE simulation must be captured in terms of the surface node location change and hence a moving mesh approach is implemented by using the Deformed Geometry module which is available in COMSOL Multiphysics software and taking temperature-dependent heat transfer coefficient, density of material, and heat of evaporation. 47 In order to handle changing geometry during the simulation process, the hyper-elastic mesh smoothing technique is used. While a default value of 0.5 is provided for the smoothing tuning parameter, a value of 0.1 is selected in the present simulation work.
The crater profiles are predicted for a single pulse with laser fluence ranging from 0.84 to 8.4 J cm−2. The material displacement for each iteration is traced by using domain point probe at center location in z coordinate. The crater formation time (tc) is obtained with the help of a stop condition by using a logical function on the material displacement. When the crater depth reaches the maximum value, the stopping condition is fulfilled. After complete formation of crater, the material cools down from the evaporation temperature.
Since the cooling can be considered as a heat conduction problem, an attempt is made to use 2D axi-symmetric heat conduction (one-step) model as in equation (17). 48
where ρ – density of material, Cp – specific heat capacity and k – thermal conductivity. The heat source term Q (r, z, t) considers the time period from tc to 0.1 ms. The thermal properties of Ti6Al4V alloy for a wide range of temperature are obtained from the literature and are given in Table 2. The 2D axi-symmetric heat conduction model is numerically solved using COMSOL Multi-physics software. In the present simulation, two sperate time dependent solvers are used for TTM and HCM respectively.
Temperature-dependent thermo-physical properties of Ti6Al4V alloy for heat conduction model.
Multi pulses laser ablation
The hybrid TTM-HCM is further extended to simulate multi pulses laser ablation for 5, 10, and 15 pulses with a separation time of 0.1 ms between successive pulses (repetition rate of 10 kHz). The surface temperature at the end of each pulse is taken as the initial condition for the next pulse mimicking the laser percussion drilling. In order to capture more accurate and efficient pulse train, the periodic load (on and off) condition is applied through an explicit event interface module, which is available in COMSOL Multiphysics software. The time-dependent solvers are used taking respective solver setting and parameters based on an automation program in COMSOL software. For 15 pulses, 30 time-dependent solvers are needed to complete the simulation run. The results for each pulse can be extracted automatically in a post-processing program.
Programs were run on Lenovo, Thinkstation P310 signature edition, with 64-bit OS, x64 based processor having Intel(R) Core (TM) i7-6700 CPU 3.40 GHz 3.41 GHz and 32 GB RAM and it was connected to main server at the institute computer center having COMSOL Multiphysics (Version 5.4) Academic licensed version.
Experimental work
Single and multi-pulses laser ablation experiments have been carried out on a 2 mm thick Ti6Al4V plate using the femtosecond laser micromachining system. The laser source (Make: Spectra physics, USA: Model: Spitfire Ace) consists of solid-state Ti: Sapphire laser source with regenerative amplifier having a central wavelength of 800 nm, pulse duration <120 fs and pulse energy of 1.2 mJ at a high repetition rate of 10 kHz. The laser unit is integrated with a computer-controlled workstation consisting of a galvoscanner and motion controller (Make: Aerotech, USA). The workstation is synchronized with a high-resolution in-line CCD camera to focus the work surface. The motion controller consists of X-, Y- linear stages and Z-stage to control the laser beam focal point through the 50x objective lens with a focal length of 100 mm. A fast-mechanical shutter is used to control the number of pulses for multi pulses ablation process.
A series of experiments were carried out at laser fluence in the range of 0.84–8.4 J cm−2 for 1, 5, 10, and 15 pulses with laser pulse duration of 100 fs. All the ablated surfaces were cleaned by placing the sample in an ultrasonic bath with a 0.5 molar sodium hydroxide (NaOH) solution (a very weak etchant) for about 20 min.
Measurement of laser-ablated crater morphology
Characterization of femtosecond laser ablated crater profile is a very challenging task, as the dimensions are in the sub-micron range.18,23 In the present work, laser ablated crater profile are measured using a laser scanning confocal microscope (LSCM). LSCM (Make: Olympus, Japan; Model: LEXT 4000) has a resolution of 120 nm in the lateral direction and 10 nm in the vertical direction. The Gwyddion 2.48 (open source) analysis software is used to analyze confocal scanned data. 50 The crater profile is extracted from scanned data by moving the cursor along with the deepest point.
The surface morphology of the ablated crater was also studied using scanning electron microscope (FEI, Quanta 400). The laser-ablated specimen was then sectioned, and hand polished to reveal the cross-section of the crater. The estimated values from destructive technique and LSCM measurement techniques were found to be in good agreement. Therefore, the crater profiles for all experimental conditions were measured by LSCM.
Results and discussion
The simulation results are presented and discussed first, followed by the comparison with the experimental values. The simulation is quite challenging as time scales involved are fs (pulse duration), ps (crater formation), ns (pulse interval) and ms (in case of multi-pulses).
Numerical simulations
The convergence tests were run for the selection of element size, taking more rapidly varying electron and lattice temperatures as responses. Only typical results are presented here for single pulse laser ablation at F: 8.4 J cm−2 with temperature dependent material properties and without material removal. For predefined element numbers of 150 (radial) and 300 (axial) using calibration with plasma module, the maximum electron and lattice temperatures reached were 29,433 K and 23,880 K with changes of 0.36% and 0.24% respectively with reference to previously selected element numbers of 100 × 200. Taking higher element numbers of 200 × 400, minimum and maximum element sizes were obtained as 0.01 and 40 µm in both directions resulting in a high-density mesh at the irradiation zone. Subsequently, all the simulations were run with 200 × 400 elements resulting in a total number of 80,000 elements, in order to capture the crater geometry properly. The simulation results are dealt in three separate sub-headings, namely laser-material interaction, single pulse laser ablation and multi pulses laser ablation, to provide insights into the percussion drilling process.
Laser-material interaction
During the laser-material interaction, high-intensity pulse is absorbed by free electrons on the top layer. The electrons are excited to a highly non-eq27 state, and the electron temperature rises rapidly as in Figure 2(a) showing the simulation results with temperature dependent properties. The heat energy is transferred to the lattice system through electron-phonon coupling. Even though the pulse duration is 100 fs, the electron temperature rises to a peak value of 15,590 K at 0.3 ps. The lattice temperature also rises to 10,154 K and both temperature systems reach an equilibrium value at teq of 10.5 ps. When the laser fluence is increased to 8.4 J cm−2, a higher peak electron temperature of 29,183 K is reached (Figure 2(b)), whereas lattice temperature reaches an equilibrium value of 23,244 K with a slight delay (teq = 11.48 ps). This shows that the laser-material interaction is more intense as the laser fluence is increased. For comparison, the temperature distributions for simulation with constant temperature (300 K) properties are superimposed in Figure 2. In this case, higher values of peak electron and equilibrium temperatures are reached in shorter times.

Predicted temperature distributions with temperature dependent and constant temperature material properties for single pulse laser ablation (Te – Electron temperature, Ti – Lattice temperature) (a) F: 0.84 J cm−2 and (b) F: 8.4 J cm−2.
Single pulse laser ablation
In simulation of laser ablation process, the material removal must be considered along with laser-material interaction. The temperature distributions obtained in simulation with moving mesh method are shown in Figure 3 for F: 8.4 J cm−2. Insets in the figures show that the temperatures obtained with 10 fs time-step within the pulse duration have the same trend as the results obtained with an initial time-step of 100 fs. Therefore, temperature distributions obtained with the time-step of 100 fs followed by adaptive time-steps up to 0.1 ms are considered for discussion. Considering temperature dependent properties, it is seen that the electron temperature reaches a peak value of 29,969 K which is closer to 29,183 K in Figure 2(b). The lattice temperature reaches the evaporation temperature (3560 K) well ahead of the equilibrium point (teq = 7.68 ps) as shown in Figure 3(b). The result in Figure 3(a) for simulation with constant temperature properties shows a peak value of 738,560 K using moving mesh method in comparison with 2853,400 K in Figure 2(b) without material removal. This is due to the high energy associated with electron system, as less energy is transferred to the lattice system.

Predicted temperature distributions with moving mesh method for single laser pulse ablation at F: 8.4 J cm−2. Inset showing results for 10 fs time-step. (Te – Electron temperature, Ti – Lattice temperature): (a) constant temperature properties and (b) temperature dependent properties.
The element removal method adapted by Kiran et al. 20 can be used to predict the crater geometry in single-pulse ablation, but cannot be extended to percussion drilling with multi pulses. Therefore, the moving mesh approach is used in the present work. When the material reaches the evaporation temperature, it directly evaporates forming gaseous phase. The dimensions of fully formed crater and the time taken (tc) are given in Table 3 for simulation with constant temperature and temperature dependent properties using moving mesh method. The crater dimensions are slightly higher with temperature dependent properties and time taken are also longer than those with constant temperature properties. Figure 4 shows simulated crater geometries after they are fully formed along with the respective times taken at laser fluence of 8.4 J cm−2. With temperature dependent properties, it takes 8.72 ps (after reaching equilibrium at teq = 7.68 ps) for formation of crater with diameter of 18.20 µm and depth of 0.407 µm (Figure 4(b). While considering constant temperature properties (Figure 4(a)), a crater of 16.85 µm diameter and 0.393 µm depth is fully formed at a shorter duration of 1.95 ps (after reaching teq at 1.86 ps). The formation of smaller crater confirms the fact that energy transferred to the lattice system is lesser in the latter case.
Simulation results for single pulse ablation based on moving mesh approach.

Contour images in moving mesh approach for single pulse laser ablation at F: 8.4 J cm−2 (a) constant temperature properties, time (tc): 1.95 ps and (b) temperature dependent properties, time (tc): 8.72 ps.
The simulation of single laser ablation using TTM was continued even beyond the complete formation of crater till 0.1 ms (pulse interval) and the CPU times were recorded for all cases and included in Table 3. The CPU times increased from nearly 1.5 to 2.5 h in case of simulations with constant temperature properties when laser fluence was varied from 0.84 to 8.4 J cm−2. With temperature dependent properties, the CPU times increased from nearly 2.5 to 5.5 h. The worst-case estimate was that the simulation of laser ablation with 15 pulses at a fluence of 8.4 J cm−2 would require nearly 82.5 h. Therefore, it was decided to use HCM after complete formation of crater (time tc). The CPU time could be reduced from 5.5 to 2.25 h for a single pulse ablation at F = 8.4 J cm−2. The temperatures at the end of each pulse from TTM and hybrid TTM-HCM did not show much difference.
Multi pulses laser ablation
Figure 5 shows surface temperatures for laser ablation with 5 pulses at F: 8.4 J cm−2. The final temperature at the end of one pulse becomes initial temperature for the next pulse. Also, the crater surface created in the previous pulse is the initial surface for the subsequent pulse. The present work therefore mimics the practical condition in the percussion drilling. In this simulation, the program was run using hybrid TTM-HCM for 15 pulses seamlessly taking a specific value of fluence. For F: 8.4 J cm−2, the surface temperature rises to 305.4 K at the end of 5 pulses and after 15 pulses it reaches 308.0 K. The relevant results were extracted from the final output at the end of 1st, 5th, 10th, and 15th pulses. The crater dimensions obtained for different number of pulses at different fluence are listed in Table 4.

Predicted surface temperature for multi-pulses laser ablation for 5 pulses with temperature-dependent material properties based on heat conduction model, F: 8.4 J cm−2; nP: 5 pulses.
Comparison of experimental values with simulation results based on hybrid TTM-HCM using temperature dependent properties and moving mesh method.
Deviation (%) = (dsim−dexp)/dexp or (hsim−hexp)/hexp.
Comparison with experimental results
The experimentally measured crater topography is compared with simulation results for single pulse and multi pulses laser ablation process.
Crater diameter and depth
SEM image of laser-ablated crater surface at laser fluence of 8.4 J cm−2 for 15 pulses is shown in Figure 6(a) along with its cross-sectional image in Figure 6(b). The crater depth is estimated as 2.025 µm at the deepest point. The crater profile was extracted from the scanned data of LSCM by moving the cursor along the deepest point. The deepest crater profile dimensions were estimated to be 16 µm diameter and 2.26 µm depth at a laser fluence of 8.4 J cm−2 after 15 pulses. As LSCM was able to capture the deepest point of the crater profile, the cross-sectional profiles of the craters obtained for 1, 5, 10, and 15 pulses at laser fluence of 8.40 J cm−2 were measured using LSCM as shown in Figure 7. A second-order polynomial curve was fitted over the crater profiles and the crater dimensions were estimated. Average of five samples at each laser fluence ranging from 0.84 to 8.40 J cm−2 for 1, 5, 10, and 15 pulses are included in Table 4.

SEM image of laser ablated surface, F: 8.4 J cm−2; np:15 pulses: (a) surface morphology of crater and (b) cross-sectional view.

Crater-sectional profiles of femtosecond laser ablated craters measured using laser scanning confocal microscope (LSCM), F: 8.4 J cm−2; np: 1, 5, 10, and 15 pulses.
Single pulse laser ablation
The experimental (solid black line) values and simulation results based on the temperature-dependent properties (blue dashed line) for single pulse laser ablation are shown in Figure 8(a) as laser fluence varies from 0.84 to 8.4 J cm−2. The simulation with temperature dependent properties over-predicts crater diameter, while crater depths are under-predicted in comparison to the experimental results. The same trend is observed for values predicted with constant temperature properties (red center line). However, the predicted crater diameters are closer to the experimental values to some extent than the values predicted using temperature dependent properties (blue dashed line). It is observed that average absolute deviation between the results predicted with temperature-dependent properties and measured crater diameters is 33.91 %, while it is 28.50% for crater depths (Table 4).

Crater diameter and depth in laser ablation for different number of pulses: (a) number of pulses, nP: 1, (b) number of pulses, nP: 5, (c) number of pulses, nP: 10 and (d) number of pulses, nP: 15.
Exact comparison with the results reported in the literature is not straightforward as there are always some differences. One research work closer to the present investigation is on the ablation of Ti6Al4V using 10 ps laser. For a quick comparison, the values extracted from the plots included in a published paper 43 are cited here. The experimental values of crater depth were in the range of 26–314 nm and the values predicted by the critical point phase separation model lie between 20 and 205 nm for laser fluence in the range of 2–10 J cm−2 exhibiting deviation ranging from 13% to 35%. One must note that the laser used was different, and the simulation was done using a one-dimensional model. There is another work in which the results of single pulse laser ablation are reported for Ti6Al4V using femtosecond laser with fluence varying from 0.84 to 8.4 J cm−2. 20 In this simulation based on 2D axi-symmetric model, the electron heat capacity had a linear dependence on electron temperature and the electron-phonon coupling at room temperature was taken. The metal removal was achieved by removing elements which were above the evaporation temperature for Ti6Al4V (3560 K). The maximum measured crater depth of 1004 nm was observed and the deviation from the simulated value was 61% at a laser fluence of 8.4 J cm−2.
Multi pulses laser ablation
The experimentally measured values of crater diameter and depth from the cross-sectional profiles traced using LSCM for different number of pulses (5, 10, and 15) are also included in Table 4. The predicted crater dimensions (dashed blue line) and experiments (solid black line) are plotted for 5, 10, and 15 pulses for different laser fluences in Figure 8(b) to (d) respectively. The predicted dimensions with constant temperature properties show mixed trend for different number of pulses and, in general, the deviations from the experimental values are higher at higher fluence. Average absolute deviation of predicted values from the measured values increases as the number of pulses increases. The deviation reaches a maximum value of 76% for depth of crater formed at 15 pulses and fluence of 8.4 J cm−2. The overall average deviation values considering np = 1, 5, 10, and 15 and fluence varying from 0.84 to 8.4 J cm−2, are 19.70% and 39.57% respectively for the crater diameter and crater depth.
The deviation of simulation results from the experimental data is low for crater diameter in comparison with that of crater depth in multi pulses ablation. This is because experimentally obtained craters exhibit irregular edges with bumps at the edges, making the assessment of crater depth somewhat difficult. Such larger deviations between simulation and experimental results are already reported in laser ablation. 14 Even though the changes in the optical properties with temperature have been incorporated in the present work, the dynamic changes in the optical properties under multi-pulses laser ablation still play a significant role in the prediction of crater depth. 8
The CPU time taken was maximum for the 15 pulses simulation at F: 8.4 for, with a total time (h:min:s) of 21:10:12 (TTM: 20:20:42 plus HCM: 0:49:30). It is clear that the program manager in COMSOL selects appropriate system parameters to minimize the CPU time.
Comparison of crater volumes
The measured crater diameter deviates from the predicted values to a larger extent in comparison with the crater depth in laser ablation process. Since laser machining comes under toolless machining category, it is prudent to compare the volume of metal removed. Also, a single parameter allows effective comparison of machining performance and hence the crater volume may be a preferred parameter. It can also serve as a bench mark for validating different simulation approaches and the assumptions made in their development. Therefore, an attempt is made in this work to estimate the crater volume taking the crater geometry to be a paraboloid using equation (18).
where Vc – crater volume, d – crater diameter and h – crater depth.
The comparison of experimental and simulation results of ablated crater volume for single and multi-pulses at laser fluence of 0.84 and 8.4 J cm−2 is shown in Figure 9. It is observed that the results for the single pulse match closely and, as expected, the volume of metal removed by ablation is higher with high fluence laser. With increase in number of pulses, the volume of metal removed shows increasing trend and the increase is steeper with high fluence than low fluence. Like measured crater depth, measured volume also shows wider range of scatter. Interestingly, the predicted values are higher than the experimental values and the gap between them widens with increasing number of pulses. The lower experimental values are attributed to the incubation effect which leaves behind a layer of irradiated material with changes in the properties to certain extent. This incubation effect seems to increase with the number of pulses and the laser fluence. The researchers have attempted to account for the incubation effect using empirical relations for crater diameter and depth separately. 23 The incubation effect based on crater volume may be a unique one and further investigations can be carried out to arrive at the value from the basic principles.

Comparison of crater volume for multi-pulses laser ablation, F: 0.84 J cm−2and F: 8.4 J cm−2.
Conclusions
In the present work, numerical simulations are carried out with temperature-dependent thermo-physical and optical properties based on Drude’s model to investigate the single and multi-pulses femtosecond laser ablation at high laser fluence. It is found that the electron and lattice temperatures reach equilibrium temperature after certain period and equilibrium time is longer with increased laser fluence, due to strong energy coupling between the electron and the lattice system.
The laser ablated crater morphology is predicted based on the moving mesh approach in COMSOL Multiphysics software for laser fluence ranging from 0.84 to 8.4 J cm−2. It is shown that the material removal starts as soon as the surface temperature reaches around 3560 K and continues even beyond the equilibrium point till complete formation of crater.
The laser ablated crater morphology is evaluated experimentally and compared with laser ablated crater morphology predicted in the simulation. It is observed that the average absolute deviations between experimental and simulation results with temperature-dependent thermo-physical and optical properties are slightly higher for crater diameter than that for crater depth in case of single pulse laser ablation.
The surface temperature rise due to previous pulse in multi-pulses laser ablation process is obtained from 2D axi-symmetric heat conduction model. Even though the surface temperature rise predicted at laser fluence of 8.4 J cm−2 is only 8 K for 15 pulses, it can become significant when the number of pulses is in the order of hundreds.
From the results for multi-pulses ablation, it is seen that the average absolute deviation between simulation results and experimentally measured crater depth is much higher than that for the crater diameter by nearly two-fold for pulses in the range of 5–15. The crater volume proposed in the present work shows promising results.
Since the proposed simulation approach has been validated by carrying out the experiments systematically, it is useful for the selection of process parameters such as repetition rate and laser energy in percussion drilling of Ti6Al4V alloy.
The present work considers only evaporation temperature as the criterion for determining the material removal. Future extension of the work is possible by incorporating phase explosion and other mechanisms to arrive at the material removal.
Footnotes
Appendix
Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
Funding
The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The authors would like to thank Dr. Ravi Bathe for his help in conducting experiments at ARCI, Hyderabad. The authors are also thankful to the sponsoring agencies DST (Grant No. SR/S3/MERC-68/2004 dated 08-06-2007) and Aeronautical Research and Development Board (ARDB), Government of India (Project number: ARDB/01/2031768/M/I dated 10 Aug 2015) for providing necessary facilities at the Manufacturing Engineering Section, Department of Mechanical Engineering, IIT Madras.
