Research and Optimization of the Influence of Process Parameters on Ti Alloys Surface Roughness Using Femtosecond Laser Texturing Technology

Abstract

One of the major disadvantages of Ti alloys is their poor wear resistance. To increase their wear resistance, before applying a wear-resistant layer, the surface of the substrate should be carefully prepared to ensure the required coating adhesion. Femtosecond laser (fs) texturing is a technology that can be used for surface texturing of Ti alloys because it enables a controlled heat input on a small surface area. The process of laser texturing is very sensitive to the choice of input parameters, such as the number of passes (P) and laser power (W), the choice of which may significantly influence the ultimate surface roughness values (R_a). It cannot be expected that by using the fs process a given default R_a value will be achieved, but it is assumed that the obtained roughness values will be within the given interval. As a result of this research with a significance level of 95% using a design of experiments (DOE) and Monte Carlo simulations, a general linear model of R_a = f (P, W) and optimal input parameter intervals (P and W) of laser texturing were obtained both for the given interval as well as for the default surface roughness value (R_a). Considering that an industrial process is involved here, a process performance capability index (C_pk) has been also defined, which shows that optimal process parameter intervals give roughness values for the given interval or given default roughness value.

1. Introduction

The service life of all moving parts used in transportation, production, and energy generation highly depends on their surface quality. Therefore, to extend their durability, it is important to research how to reduce friction and wear. Friction and wear depend on the choice of tribopairs and contact surface roughness. The study of surface and surface preparation is the most active field in tribology. In the field of surface preparation, various techniques of surface texturing have been developed, mainly in conventional procedures [1,2,3]. Application of advanced technologies to improve the wear resistance of metallic materials has become very interesting in recent years [4,5,6,7,8,9]. Titanium and its alloys are used in aeronautical, marine, and medical applications because of good corrosion resistance, high strength, non-magnetic properties, and biocompatibility. One of the main problems of Ti alloys is their poor wear resistance and this requires careful preparation of the surface, that is, the surface of the tribopair should be as less-rough as possible [10]. Advanced surface modification methods such as friction stir processing (FSP), laser surface modification, and micro-arc oxidation (MAO) offer better surface properties of Ti and its alloys [11]. Laser processing techniques are used to improve surface topography [12,13] because they enable a controlled heat input on a small surface area [14]. Research has proved that this technology can modify the micro- and nano-structure of the material surface [15,16], increase nano hardness [17,18], improve wear resistance [19,20], and modify surface wettability [21,22]. In [23], the author showed that laser microscale textures of sliding surfaces of the tribopair reduce friction and wear. Laser texturing of tribosurfaces influences the behaviour of the surface wetting angle and the authors showed that such surfaces can be used as an alternative to the use of high-performance lubricants with complex formulations. Maalouf et al. [24] investigated that small-scale texturing of titanium surfaces makes it possible to couple-enhanced osseointegration with potential antibacterial properties that are more sensitive at the nanoscale than the microscale. In Ref. [25], the authors examined the influence of laser peening with varying fluence on the fretting behavior of Ti-6Al-4V samples. The results showed a decrease in the wear of samples. Despite better wear results, the roughness of laser-textured surfaces produced higher wear rates than needed. Since wear resistance is a weakness of Ti alloys, they need to be coated. TiN, TiCN, and diamond-like carbon (DLC) coatings are coatings with excellent tribological properties and are often used to coat titanium alloys [26,27]. A satisfactory degree of surface roughness is important for good adhesion of the coating and for coating thickness. Femtosecond (fs) laser texturing of surfaces is used to examine the influence of process parameters on the roughness of metallic materials. Femtosecond lasers have a higher degree of precision in surface texturing which means that the heat affected zone is much smaller than those formed with nanosecond and microsecond laser pulses for the same laser spot size [28,29]. They enable the formation of different shapes on metal surfaces such as grooves, dimples, and waves [30,31,32,33]. Bonse et al. in [34] presented a review of the tribological properties of steel and titanium alloy surface morphology (ripples, grooves, and spikes). By changing the femtosecond laser parameters such as the speed, power, and size of the laser beam, it is possible to produce surface roughness that can be hydrophobic or hydrophilic [35,36].

Femtosecond lasers are increasingly used in the surface texturing of metal alloys to increase the wear resistance of coating substrate systems. However, the correlation between the obtained surface roughness of the alloy substrate (R_a) depending on the change in fs laser parameters, primarily the number of passes (P) and power (W) values, has not been sufficiently investigated. The texturing process is sensitive to the selected process parameter values (P and W), which to a greater or lesser extent influence the changes in the surface roughness values (R_a). It cannot be expected that by using the fs process a given default R_a value will be achieved, but it is assumed that the obtained roughness values will be within a given roughness interval. The same is true when choosing the process parameter values, which will also have interval values. The objective of this research is to define with a significance level of 95% an interval with acceptable roughness values for the given process parameters interval values. An overview of research steps and methods is shown in Figure 1.

According to the above and Figure 1, the main benefits of this research are as follows:

By applying the Monte Carlo method and the performance capability index (C_pk), an optimized robust regression model was obtained, which at the 95% significance level gives the following answers important for laser setup texturing process:

• For a given interval of roughness R_a, the intervals of process parameters P and W are obtained.
• For given intervals of the process parameters P and W, the roughness interval R_a is obtained.

2. Materials and Methods

2.1. Preliminary Research

Wrought Ti-6Al-4V alloy samples were prepared in dimensions 17 mm × 17 mm × 7 mm. Before proceeding with laser treatment, all samples were wet and fine abraded with silicon carbide sandpapers (gradations P320, P1000, P2000, and P4000). Then, all samples were mechanically polished on the Struers DAP-V device in two phases, first, with a finer diamond paste of 9 μm and after that, with a liquid silica particle granulation of 0.03 μm.

The surface roughness measurements were carried out according to norm ISO 4287:2009 on pre-treated samples and the laser-treated samples in four repetitions using a stylus device, SurfTest SJ-500 (Mitutoyo, Kanagawa, Japan). The Gauss filter was used for data filtering, a cut-off value, (λ_c = 0.25 mm), evaluation length l_n = 1.25 mm, traverse speed 0.1 mm/s, and a tip of radius r = 2 µm (60°). After grinding and polishing samples, the average maximum height of samples was R_a = 0.038 ± 0.0025 μm, the average maximum height of the profile was R_z = 0.277 ± 0.0111 μm and the maximum roughness depth was R_max = 0.315 ± 0.026 μm (Figure 2).

A Ti:Sapphire laser system (SpectraPhysics Spitfire seeded by Tsunami, M 3960, Spectra-Physics, Inc., Mountain View, CA, USA) operating at a 1 kHz pulse repetition rate was used in the experiments. The laser system delivered 120 fs pulses at 800 nm with a maximum laser energy per pulse of 0.8 mJ. The energy stability was within 2%. In typical measurements, less energy was needed, and the original beam was attenuated by a set of absorptive neutral density filters. The laser pulses were focused on the sample surface to a 170 μm radius spot. The size of the spot was measured by strongly attenuating the beam and then sending the focused beam to a CCD camera.

A sample was fixed to an XYZ computer-controlled translational stage, with the Z-axis being the laser beam direction and the surface of the sample being in the XY plane. The texturing was performed by scanning the sample in the XY plane, line by line, at a constant line velocity of v_x = 10 mm/s and with an inter-line (y) spacing of h = 250 μm.

For preliminary and DOE experiments, laser fluences of 0.057, 0.139, 0.205, and 0.297 J/cm² were used, which correspond to average laser powers of 52, 126, 186, and 270 mW, respectively. We made 3, 5, or 7 consecutive layers for each of the aforementioned fluences. In the rest of the paper, we refer to average laser powers for easier calculations.

The distance between two consecutive spot centers (d) is expressed as d = v_x/f, where v_x is the scan velocity, and f is the repetition rate of the laser. The effective number of pulses per unit area, N_eff, can be calculated as the ratio of the single spot area and the area of the single cell (d × h) [37]. The single cell in a 2D scanning is the area defined by two consecutive pulses (separated by d) in the row and two corresponding pulses in the next row, which is separated by h,

$$N_{\mathrm{eff}}=\frac{{\pi \omega }_0^2}{d h}$$

(1)

where ω₀ is the radius of the focal spot. In our measurements, N_eff = 36.3 for each layer.

Scanning electron microscopy analyses were performed using a TESCAN VEGA 5136 (TESCAN Brno, Brno, Czech Republic) mm system to analyse the laser-treated surfaces. The roughness profile and the SEM image of the surface of the reference sample-RS, the mechanically pre-treated sample, are shown in Figure 3.

To initially and roughly determine the influence of process parameters (number of passes (P) and power (W)) on the surface roughness (R_a) of Ti alloys, a preliminary experiment with 27 samples was performed. The process parameter values were intuitively determined: the number of passes P_o = [3,5,7], and power (mW) W_o = [80,100,128]. After the preliminary experiment, a 2-factor interaction DOE experiment was performed, where the process parameter values were no longer intuitively chosen, but rather based on randomly chosen values from a predefined distribution for each process parameter. Considering that the research also foresees the determination of the general functionality of R_a = f (P, W), a general roughness model was determined using Monte Carlo simulation based on a R_a-surface equation DOE experiment. The last research task foresees the optimization of the influence of process parameters on the surface roughness of the alloy. Furthermore, the process performance capability index (C_pk) shall be defined, which must show whether the optimal process parameters intervals give the roughness value for the given interval or the given default roughness value to achieve the given C_pk value.

For each of the three experiments (preliminary, DOE, and Monte Carlo), descriptive statistics results will be presented and analyzed in detail such as Mean, Standard Deviation, Variation Coefficient, Skewness, Kurtosis, Min, Max, Q₁, Median, Q₃, 95% CI-Mean, 95% CI-Median, AD Normality Test, Grubbs Outlier Test. In terms of graphical representations, the following diagrams will be shown: Frequency, Probability, and 95% Confidence Interval plots. For all three experiments, the stepwise regression procedure will be used to obtain surface equation regression models. The significance of the effect of process parameters on the roughness R_a for each surface equation model will be evaluated using analysis of variance (ANOVA) by comparing the F and p values. In addition, 3D surface plots, and R_amodel–R_a correlation diagrams will be shown. Values of R², R²(adj), and Test data R² will be calculated. Surface equation model result values will be compared using multiple sample comparison of roughness tests: ANOVA, Tukey’s pairwise comparisons, Tukey difference of means, Interval plot, Box plot, and Individual plot. In the end, the optimization of the process parameters will be carried out.

2.2. DOE Experiment

The DOE experiment was used in this paper to improve the regression model of R_aoLM obtained in the preliminary experiment where the values of the model’s input factors (P_o and W_o) were chosen by intuition. It is assumed that with the DOE experiment, better response variable surface roughness results (R_a) can be achieved as a result of a random choice of input variables (process parameters P and W) under the corresponding intervals.

DOE experiments are used in industry to systematically study the processes and factors (variables) that affect product quality. Once the process conditions and variables affecting the product quality are identified, the analysis of DOE experiment results can significantly improve the manufacturing process and make it more reliable, effective, and of better quality. In this paper, the full factorial DOE was used because of different levels of factors. The structure and the phases of DOE with default parameter values are shown below:

Step 1. Definition of the output variable (response variable).

Name	Analyze	Goal	Impact (1–5)	Sensitivity	Min	Max
Ra	Mean	Minimize	5.0	High	0.01	0.5

Step 2. Definition of input variables (input factors, process parameters) and their levels.

Name	Units	Type	Role	Levels
P	-	Categorical	Controllable	3; 5; 7
W	mW	Categorical	Controllable	52; 126; 186; 270

Step 3. Choice of the type of DOE experiment.

Name	DF	Randomized	Nr. of Replicates	Total Runs	Total Blocks
Full-Factorial	22	Yes	2	36	3

Step 4. Choice of the interaction model between factors.

Model	Interactions
2-factor interactions	P, W, PW

Step 5. Analyses of the experiment results. Some of the parameters from the experiment structure shall have the following meaning.

Parameters	Meaning
Analyze	Parameter of Interest
Goal	The goal of the experiment is Ra minimization
Impact	Relative importance of factors
Sensitivity	The importance of achieving the best-desired value
Min-Max	Value interval of the response variable
Type	Factor variable type
Role	The role of the factor in terms of adjustment
Levels	Factor levels

2.3. Monte Carlo Simulation Experiment

Monte Carlo simulation is used to generate a set of output (dependent) variables based on the distributions of input variable data, either known or estimated in advance, most of the time, to avoid performing a large number of experiments necessary to determine a mathematical model for assessing the output variable.

Monte Carlo simulation was conducted in this paper to improve the linear regression model of R_a. Due to many variations of input factors (P and W), it was assumed that the newly obtained MC regression model should be better at evaluating roughness R_a compared to the R_aLM DOE model.

The coupling of input variables (P, W) with the output variable (R_a) necessary for performing the simulation was made based on Formula (2) for evaluating roughness R_aLM using the obtained DOE model:

$$R_{\mathrm{aLM}} = \mathrm{f}(P_a, W_a)$$

(2)

The estimated roughness value $R_{\mathrm{amcLM}}$ of the MC model for input variables P and W was obtained based on the Uniform distribution by a random selection of the X (P or W) input variable value between the lower limit (a) and the upper (b) limit. The parameters (a) and (b) are the continuous uniform distribution parameters (Figure 4).

For any value of the X in the interval $a \le { x}_{ 1}{< x}_2 \le b$, the probability was calculated using the formula:

$$P\left(a \le X \le b\right) = \frac{x_{2 }-{ x}_1}{(b - a)}$$

(3)

The uniform distribution parameters (a) and (b) for the variables P and W were determined based on their minimum and maximum values as follows:

(a)_P = 3.0

(b)_P = 7.0

(a)_W = 52.0

(b)_W = 270.0

The estimated number of iterations in the Monte Carlo simulation was 10,000.

After performing the Monte Carlo simulation, the multiple regression method was used to generate a statistical model to determine the influence of two or more independent quantitative factors X for the dependent variable Y. The method included the option of performing a stepwise regression where a subset of variables X with the greatest impact on the change in the value of the dependent variable Y was selected. Moreover, the impact was also identified by the ANOVA test. The expected outcome was the response surface equation model for the surface roughness prediction model R_amcLM. To determine the efficiency of the model, the following parameters were determined: R², R²(adj), Test data R², and correlation coefficient R_amcLM vs. R_amc. The response surface (P_amc, W_mc) vs. R_amcLM. was represented graphically.

2.4. Multiple Sample Comparison Test of Roughness

This paper suggests the use of the Multiple Sample Comparison Test to compare the sample values of two or more independent variables (factors). In this case, factors were represented by the following roughness values: R_ao, R_aoLM, R_a, R_aLM, R_amc, and R_amcLM. The testing was carried out to establish the existence of significant differences between the means of populations from which the samples had been drawn. In addition, this paper intended to ascertain whether there was a significant difference between the linear models (R_aLM, R_amcLM) of roughness obtained from DOE and Monte Carlo experiments and the roughness model (R_aoLM) obtained as a result of the preliminary experiment for intuitively selected values of process parameters P and W. The testing of uniformity (homogeneity) and the comparison of samples were done by testing the mean values of factors using ANOVA under the hypothesis that all samples had the same mean value. All the equations necessary for the multiple sample comparison test and their meanings are given in Figure S1.

Furthermore, Tukey’s pairwise comparison tests (multiple range tests) were suggested to detect the means of those factors that differed significantly from others. To visually display the ANOVA results the following graphic representations were employed: Tukey’s Difference of Means, Interval plot, Individual plot, and Box-and-Whisker plot.

2.5. Optimization of Process Parameters

Surface texturing is an industrial process whereby certain roughness required by the client needs to be achieved in a given interval or with desired default values. Automatically, the manufacturer is required to identify and optimize the parameters of the production process to minimize unacceptable losses, that is, rejects. Given the request, the ultimate objective of this study was to determine optimal values of P and W parameters to achieve the desired surface roughness. As shown in Figure 1, two types of process parameters optimization method were foreseen:

• Determination of optimal P and W interval values based on the default R_a interval.
• Determination of optimal P and W interval values based on the desired default R_a value.

As the key performance indicator of optimization, the process capability ratio (C_pk) was suggested and calculated using the Formulas under (4)–(6) [38].

$$C_{\mathrm{pk}}=\min \left(C_{\mathrm{pu}}, C_{\mathrm{pl}}\right)$$

(4)

$$C_{\mathrm{pu}}=\frac{\mathit{USL} - \mu }{3\sigma }$$

(5)

$$C_{\mathrm{pl}}=\frac{\mu - \mathit{LSL}}{3\sigma }$$

(6)

where:

C_pk—Process capability ratio

C_pu—Process capability Upper specification limit (USL)

C_pl—Process capability Lower specification limit (LSL)

µ—Sample mean

σ—Sample standard deviation

The once-established C_pk represented a deviation from the optimized average surface roughness (R_a) and previously defined upper and lower specification limits (LSL/USL) for surface roughness. The benchmark was C_pk > 1.33. The greater the C_pk, the greater the likelihood that during the surface texturing of the sample, the values of R_a would remain within the default LSL/USL limits for parameters P and W.

3. Results and Discussion

3.1. Preliminary Experiment

The descriptive statistics results of the preliminary experiment for factors (P_o, W_o) and the response variable R_ao are shown in Figure 5, Figure 6 and Figure 7. A total of 27 experiments were performed to determine the surface roughness R_ao for the intuitively determined values of P_o and W_o factors. Considering that experiments were conducted at only three levels (3, 5, 7), the P_o factor had a high variation coefficient CoefVar = 33.28%, whereas the W_o factor had a mean coefficient of variation CoefVar = 19.54%. Based on Q₃, 75% of P_o data was below 7 passes, and the W_o was in 75% of cases lower than 128 mW. Based on the skewness indicator, P_o had a symmetric distribution, whereas W_o had a positive or right-skewed distribution. According to the AD Normality Test, neither factor had a normal distribution of data (p-Value < 0.005). Furthermore, according to the Grubbs’ outlier test, no outlier was found.

As shown in Figure 6, the response variable (R_ao) had an extremely high coefficient of variation (CoefVar = 91.27%), and this was probably due to a variation in factor values. The probability plot showed a huge waste of data, as almost all response variable values were outside the confidence interval lines and the mean line. It should be observed that most values (75%) of R_ao were less than 0.0820, whereas the maximum value was 0.2720. Based on the author’s experience, the lower roughness is the better to make the forming of wear-resistant layers more effective. It can be observed that, according to intuitively set factor values, the experiment generated significantly lower values of R_a concerning the maximum values obtained for particular combinations of factors. The response variable Ra did not have a normal distribution nor an outlier. The interaction plot in Figure 6 shows the effect of changing one factor and the response variable against the constant value of the other factor. Given that the interactive lines were not parallel, we concluded that there was an interaction between W_o and R_ao and a relationship between them and the value of P_o.

Given the interaction between the factors P_o and W_o and the R_ao, the stepwise regression was performed to obtain a linear model of R_aoLM = f (P_o, W_o) as well as the variance analysis to examine by means of p-values the relevance of the influence of P_o and W_o factors on the response variable value (R_ao). The obtained regression model is represented by Equation (7):

$$R_{\mathrm{aoLM}}= -0.3541 + 0.01480P_0 + 0.002041W_0$$

(7)

To obtain the regression model, 30% of data were omitted and the Forward Stepwise Selection procedure was applied by gradually introducing input variables in the model provided that they remained in the model if they had a significant effect on the response variable. From the analysis of variance, given that the p-values were relatively low compared to the standard value of 0.05, both variables, that is, factors, had a significant effect on the response variable (R_ao). The performance of the regression model R² = 40.53%, R²(adj) = 33.09%, test data R² = 60.10%, and correlation R_ao − R_aoLM = 69.60% was relatively low. As a result, the obtained R_aoLM model had a low possibility of prediction of R_ao depending upon the changes in R_o and W_o. According to the response surface plot, it can generally be concluded that higher values of R_ao were obtained with fewer passes P_o and lower power rates.

3.2. DOE Experiment

A 2-factor interaction DOE experiment with a total of 14 coefficients was performed. For the experiment, 1 response variable (R_a) and 2 experimental factors (P and W) were determined. Using the previously suggested experiment structure, a total of 36 experiments were performed. The results of the descriptive statistics of factors are shown in Figure 8. The values of the P and W factors were randomly selected from the interval [3…7] and the interval [52…270], respectively. According to the Q₁ parameter for P, 25% of the value was below 3 passes, whereas according to Q₃, 75% of the value was below 7 passes. For the W factor, according to the parameter Q₁, 25% of the value was below 101.50 mW, whereas according to Q₃, 75% of the value was below 207.0 mW. The AD normality test showed that none of the factors had a normal distribution of data (p-Value = 0), nor the outlier values according to the Grubbs’ outlier test (p-value = 1). Given the skewness values, both factors had a relatively symmetrical distribution of data (however, not a normal distribution considering their p-values). Factor P had a negatively skewed or left-skewed distribution, whereas factor W_mc had a positively skewed or right-skewed distribution of data. The Kurtosis values for both factors were between −2 and +2, as expected for data with a normal distribution. The figure also shows 95% confidence intervals for the mean and median values, and it also shows that the mean and median values are relatively close.

Figure 9 shows the descriptive statistics results for the R_a variable. The R_a had a negatively skewed or left-skewed distribution of data, and the distribution of data followed a normal distribution as shown by the AD normality test p-value of 0.3380, which is rather high. The normal distribution of data of the response variable R_a could be also inferred from the probability plot of R_a. The distribution of the response variable was such that it showed no outlier value, as illustrated by the Grubbs’ outlier test (p-Value = 1). Figure 9 also shows the interaction plot. It shows a change between one of the factors and the response variable compared to the constant value of the other factor. Given that the interactive lines are not parallel, and two lines (P = 5 and P = 7) intersect, it can be concluded that there is a relationship between the factor W and roughness, and that R_a depends on the value of the factor P.

Figure 10 shows the results of the multiple regression analysis which was performed to obtain the surface roughness model R_aLM = f (P, W) based on the results of the DOE experiment. To obtain the regression model, the fitting algorithm was used to perform the forward stepwise selection, starting from the model which included only the constant (Const.), by gradually adding one variable (factor) at a time, provided that after adding the variables they would be statistically significant (p-value = 0.00). The regression model R_aLM was obtained in two steps, provided that only in the second step (Step 2), all three steps (Const., P, and W) were significant for the response variable as confirmed by the p-value = 0 for factors P and W. The obtained regression model is represented by Formula (8):

$$R_{\mathrm{aLM} }= -0.0090 + 0.01034\mathit{P} + 0.000388\mathit{W}$$

(8)

The significant effect of all three factors on R_aLM was also confirmed by the variance analysis and by p-value = 0. In addition, the obtained values of R² = 89.31, R²(adj) = 88.13, and the correlation coefficient R_a − R_aLM = 92.40% are considerably higher than in the preliminary regression model R_aoLM. The test data value of R² = 70.40% was relevant which means that the R_aLM mode had a good roughness prediction even for those P and W values that were left out from the model-making process by a stepwise regression fitting algorithm. From the total number of datasets, 30% of data were omitted from the regression model creation procedure. The importance of the DOE assessment of the experimental values of R_a with the new model R_aLM is shown in the correlation R_a$-$R_aLM C and surface R_aLM plots.

Views of the surface roughness profile and SEM image following the laser texturing for three different laser powers (low, medium, and high laser power; 52 mW, 186 mW, and 270 mW) and passes (3, 5, and 7) are shown in Figure 11.

Figure 11 shows surface roughness profiles and SEM images of surfaces textured using three different laser powers (52 mW, 186 mW, and 270 mW) and irradiated using three different numbers of laser beam passes (3, 5, and 7). The representative profiles of each group of parameters (P and W) shown were randomly chosen because they approximately had the same profile, meaning that not the best, i.e., the most uniform profile within the given group had been chosen. The beams of 52 mW led to a decrease in surface roughness R_a (R_a = 0.016 ± 0.0031 μm) compared to the roughness of non-textured surfaces (R_a = 0.038 ± 0.0025 μm), as a lower laser beam power caused a smoother roughness profile. SEM images of the fs-laser textured samples using different laser beams powers show significant differences in the surface topography of samples processed using laser texturing with a laser power intensity of 52 mW compared to other intensities (186 mW and 270 mW). SEM images of samples subjected to texturing with an intensity of 52 mW show no changes in surface topography as compared to non-textured samples, as shown in Figure 9–52 mW. In other samples (subjected to laser texturing with laser beam intensities of 186 mW and 270 mW), ripple patterns in the ablated area as well as irregular laser deposition products such as microdroplets formed. The formation of loose microdroplets during the laser texturing process was the result of high peaks that arose on the roughness profiles at higher intensities and became more frequent as the number of laser beam passes increased. Micropores and microcracks contributed to the formation of deep grooves on the roughness profile. With the increase in power intensity of the laser beam (186 mW and 270 mW) and in the number of laser-beam passes, the intensity of the melting process increased (as a result, more microdroplets were generated). According to Pan et al. [39], the ablation that occurs at defined laser intensities occurs mainly due to spallation which is the result of the ejection of the surface layer of the material. As the melting and the formation of microdroplets on the sample surface at higher intensities increased (270 mW), it caused high surface roughness values R_a = 0.311 ± 0.0413 μm and expressed peaks on the roughness profile (Figure 11–250-7).

3.3. Monte Carlo Simulation Model

A Monte Carlo simulation was run with approximately 10,000 iterations. The assumed distribution of data for P and W was uniform within the interval of minimum and maximum factor values as indicated in the Section 2.3 Monte Carlo simulation experiment. The coupling of input factors (P_mc, W_mc) with the response variable (R_a) are connected by using Equation (8). The descriptive statistics results are shown in Figure 12 and Figure 13. The AD normality test results for factors P_mc and W_mc show that the distribution of data was not normal because the p-value was <0.005. The Grubbs’ outlier test shows that there were no outliers among the factors (p-value = 1.00). According to the Q₁ parameter for factor P_mc, 25% of the value was below 4 passes, and 75% of the value was above 4 passes. According to the Q₃ parameter, 75% of the value was below 6 passes, and 25% of the value was above 6 passes. Given the skewness values, both factors had a relatively symmetrical distribution of data (however, not a normal distribution considering their p-values). P_mc had a negatively skewed or left-skewed distribution, and W_mc had a positively skewed or right-skewed distribution. The Kurtosis values for both factors were between −2 and +2 as expected for data with a normal distribution. The figure also shows 95% confidence intervals for the mean and median values.

Figure 13 shows the descriptive statistics results of the Monte Carlo simulation for the response variable R_amc. Considering the foregoing considerations of the descriptive statistics parameters, R_amc had a positively or right-skewed distribution, which means that the distribution of data did not follow a normal distribution. The p-value of the AD normality test, which was <0.005, showed the same. In addition, the distribution of data for R_amc showed no outliers. For subsequent interval optimization of parameters of surface texturing, the consideration of the values of statistical parameters (Min, Max, Q_1, and Q₃) was of particular importance. The minimum value was 0.0139 mW, and the maximum value was 0.2028 mW. It is logical that the lower the value of roughness for the given P and W, the better, and according to Q₃, only 25% of data for R_amc had a value higher than 0.1302 mW. The intervals of input process parameters (P and W) established by optimization gave the lowest value of R_amc.

As stated above, following the Monte Carlo simulation, the multiple regression procedure was performed to obtain a surface roughness prediction model $R$_amcLM = f ($P$, $W$). As the model fitting procedure, forward stepwise selection was applied, starting with a model containing only a constant, and then adding variables (factors) to the model, provided that once they were added the variables became statistically significant ($p$-value = 0.00). The forward stepwise algorithm made it possible to remove variables from or add them to the analysis in one of the following steps if they were no longer statistically relevant (p-value > 0.05). The results of the Monte Carlo regression model are shown in Figure 14. For removing, that is, for adding the variables to the model, the threshold of $\alpha = 0.15$ was used. As much as 30% of the total number of datasets was used to test the model during the performance of the forward stepwise procedure. In this case, already in the first step (Step 1) for all the three factors (Const. P_mc, W_mc), the significant effect of the p-value was determined along with the surface equation R_amcLM (9).

$$R_{\mathrm{amcLM}}=-0.008080+0.010141{\mathit{P}}_{\mathrm{mc}}+0.000392{\mathit{W}}_{\mathrm{mc}}$$

(9)

The significant effect of all three factors on R_amcLM was confirmed by the analysis of variance (p-value = 0). Furthermore, it can be observed that the values of R² = 99.19 and R²(adj) = 99.19 were significantly higher than in the DOE regression model R_aLM, which justified the performance of the Monte Carlo simulation following the DOE experiment. The value of the test data R² = 99.18% is of particular importance, which means that the R_amcLM model had a good roughness prediction even for those P and W values which had been omitted from the process of model creation using the stepwise regression fitting algorithm. Typically, the value of the test data for R² was lower than R² obtained on training dataset, because the data from the test dataset were not included in the procedure of “learning” the R_amcLM model. The robustness of the R_amcLM model was also confirmed by a high correlation coefficient R_amcLM vs. R_amc = 99.60% as can be seen from the scatterplot in Figure 14.

3.4. Multiple Sample Comparison of Roughness

As can be seen from the multiple comparisons of roughness test in Figure 15, the ANOVA table shows the results of the comparison of means between the groups of factors (row factor) and within the factors themselves (row error). The F-value, which in this case amounted to 10.52, was the ratio obtained by making an intra- and an inter-group assessment. Given that the p-value of the F-test was lower than 0.05, we conclude that the null hypothesis was not true and that there was a significant difference between the means of the indicated factors at a significance level of 5%. As Tukey’s pairwise comparison test at 95% confidence level shows, there was no significant difference in means between the roughness obtained using the DOE experiment and the roughness obtained using the Monte Carlo simulation and the respective LM models (Group A). In addition, there was no significant difference between the preliminarily obtained roughness and the respective LM model (Group B). Tukey’s mean difference plot shows that by two-pairs comparison, all factors in Group A differ from those in Group B. Such a difference in the results of the comparison was expected considering that the values of P and W in the preliminary experiment were intuitively, not randomly selected.

3.5. Optimization of Process Parameters

The optimization intervals of input parameters P and W and the response variable R_a were determined based on the surface and contour plots of R_a as shown in Figure 16. The contour plot shows that higher roughness values were obtained in the R_a = [0.1, 0.15], P = [4,7] and W = [160,270] intervals, and lower in the R_a = [0.05, 0.100], P = [3,4] and W = [52,126] intervals.

In light of the above considerations, final optimization intervals were suggested: R_a = [0.06, 0.09], P = [3,4], W = [52,126]. Figure 17 shows the Monte Carlo simulation results (1000 iterations) for the Ra interval (R_a = 0.06, 0.09) and the non-optimized values of the process parameters (P and W). The C_pk is 0.4957, which is certainly lower than the acceptable value (1.33). In the case of a new experiment, for interval values P = [3,4], W = [52,126], we would have a waste of R_a outside the interval of 31.69%, which is unacceptable.

From the total of 10,000 samples, 1574 samples were below the LSL limit, and 3595 were above the USL limit. Following the Monte Carlo simulation of the non-optimized intervals of process parameters, the first type of parameter interval optimization was performed to obtain the highest possible C_pk and the lowest percentage of samples outside the intervals.

After consecutive optimization iterations where the values of P and W changed, and after conducting the sensitivity analysis where the standard deviations for P and W (−50% do + 50%) were changed, the results for optimized intervals of process parameters P and W for the predefined R_a interval (R_a = [0.06, 0.09]) are shown in Figure 18.

This optimization contributed to achieving a satisfactory value for C_pk (C_pk = 1.5 > 1.33), whereas for the samples out of limits, the value was 0.00%. In the case of a new experiment, the number of passes and intervals P = [3,4] as well as the intensity W = [105,132] would have to be selected. The application of P and W values from the indicated intervals during the experiment would in 95% of cases ensure R_a values within the LSL = 0.06 and USL = 0.09 intervals. The table in Figure 18 shows that the mean = 0.0749, min = 0.0646, and max = 0.0853 values practically correspond to the optimization area of R_a in Figure 16. The second type of optimization is even more important for practical application than the previous one because a certain surface roughness (R_a) of the sample is normally required in practice. Additionally in this case, the obtained process parameters of P and W guaranteed, with a statistical certainty of 95%, that the default R_a would be obtained. Given that a real production process is involved here, achieving the exact default R_a is not possible, but rather a value within the LSL and USL interval range. Figure 19 shows the simulation results for the default R_amcLM = 0.08. This latter was chosen because it is in the middle of the predefined interval for R_a as can be seen from the contour plot in Figure 16. In the case of a new experiment, the default value for R_amcLM = 0.08 would be achieved if we conducted the experiment with the following intervals of P = [4,4] and W = [111,138]. Moreover, the satisfactory value of C_pk = 1.36 > 1.33 was obtained, whereas the performance capability for the samples out of limits was 0.00%.

4. Conclusions

This study investigates the influence of process parameters (number of passes P and power W) on the average surface roughness values (R_a) of Ti alloys by applying the femtosecond laser texturing technology. The optimization of the P and W process parameters was carried out. To determine the effect of the P and W parameters on R_a, before optimizing the process parameters, three experiments were carried out: a preliminary experiment, a DOE experiment, and a Monte Carlo simulation experiment. Based on these experiments, dependency models of R_a = f (P, W) were generated, and the values of R_a generated by each model were compared. In addition, the preliminary experiment was used to measure the roughness value range R_a = 0.038 ± 0.0025 μm on the untreated sample and the surface was analyzed using SEM. Subsequently, 27 Ti-alloy samples were textured using a femtosecond laser with intuitively chosen process parameters (P_o = [3,5,7], W_o = [80,100,128]).

The texturing process was repeated three times for each value of the P_o interval with each value of the W_o interval. The descriptive analysis results show a significant variation of roughness (R_ao = 91.27%). The AD test shows that the distribution of R_ao is not normal. The ANOVA test shows that a change in the P_o and W_o parameters significantly affects the change of R_ao. The fitting performance of the generated regression model R_aoLM (7), R² = 40.53%, R²(adj) = 33.09%, test data R² = 60.10%, and correlation R_ao − R_aoLM = 69.60 is relatively low. The preliminary model of R_aoLM has a weak prediction performance of R_a. A 2-factor interactions DOE experiment showed a better result. A variation in roughness (R_a) of 39.11% is acceptable and the R_a values have a normal distribution as shown by the AD test (p = 0.3380). The fitting performance of the generated DOE regression model of R_aLM (8) is significantly better compared to the preliminary model of R_aoLM and is as follows: R² = 89.31%, R²(adj) = 88.13%, test data R² = 70.40%, and correlation R_a − R_aLM = 92.40%. The analysis of the roughness values obtained by the DOE shows that the samples textured at lower laser beam power values (e.g., 52 mW) show lower values of parameters (R_a = 0.016 ± 0.0031 μm), even in relation to the non-textured sample (R_a) = 0.038 ± 0.0025 μm.

SEM images of the surface topography of samples for other (higher) laser beam power values show changes in the Ti-alloy surface in the form of ripple patterns, and irregular laser deposition products such as microdroplets. Microdroplets cause changes in the roughness profile in the form of peaks which result in higher values of roughness (R_a = 0.311 ± 0.0413 μm for the sample 270-7).

The objective of the third experiment using MC simulation was to generate a more accurate roughness prediction model compared to the DOE R_aLM model. The objective was accomplished. The uniform distribution of data with a and b parameters (P_mc(3, 7), W_mc(52, 270)) was suggested. A total of 10,000 simulation iterations were performed and the coupling function was R_aLM. The fitting performance of the generated MC regression model of R_amcLM (9) was better compared to the DOE R_aLM model and, generally, the values were much higher: R² = 99.19%, R²(adj) = 99.19%, test data R² = 99.18%, and correlation R_a − R_aLM = 99.60%. Multiple sample comparison test of roughness (Figure 15) shows no difference in means between the DOE and MC models and displays uniformity of data values as indicated by Tukey’s pairwise comparison plot. It could therefore be concluded that the MC regression model will provide reliable predictions of roughness for any new combination of process parameters given a high value of Test data R². It can be concluded from Tukey’s difference of means plot (Figure 15) that there is a significant difference in means between the DOE R_aLM and the MC R_amcLM models in respect of the preliminary R_aoLM model. The optimization of process parameters was done based on the R_amcLM model and with regard to the method of defining the surface roughness default value. In the first method, the default surface roughness was given within the interval R_amcLM = [0.06, 0.09]. The optimized intervals of process parameters that contribute to achieving such an interval were: P = [3,4] and W = [105,132]. In the second method, the default surface roughness was R_amcLM = 0.08, and the optimized intervals of process parameters needed to achieve such an interval were P = [4,4] and W = [111,138].

For future research, we suggest optimizing the surface roughness of Ti alloys to achieve the highest wear resistance of the applied layer.