Accuracy of the subsurface damage parameters calculated by the finite difference algorithm

An important approach to characterize the full three-dimensional information of subsurface damage is to simulate the etching process of a sample reversely. The simulation starts from the morphology of the sample after the subsurface damage micro cracks being opened totally. During the etching experiment, it is possible for us to get the surface morphology at any moment. This paper presents a ﬁnite difference algorithm to simulate the morphology evolution during the etching process and then the surface’s morphology of the sample at a speciﬁc time can be obtained. Comparison between the simulated morphology and the measured one provides the clue of improving the precision of the ﬁnite difference algorithm. This method is kind of the fast calculation. In addition, the accuracy of this calculation of the corrosion model needs to be ensured. In order to improve the precision of calculation, the time interval should be set as the appropriate value by comparison and analysis. In this paper, the accuracy can be calculated through comparing the simulated result with the experimental result, and the maximum error of this method can be gained.


INTRODUCTION
Subsurface damage (SSD) refers to the residual fractured and deformed material locating under the surface region of the brittle optical materials during the shaping, grinding and lapping processes.However, the existence of SSD has a significant impact on the various properties for high-precision optical components, such as reducing the transmission performance and anti-laser damage threshold, affecting the component's strength directly, reducing the component's life and long-term stability, and so on.Thus, the measurement and characterization of SSD should be made for offering quantitative and full evaluation of the optical glasses' SSD.
The SSD layer usually locates about 1-100 µm below the surface and is composed of the crack layer and the deformation layer [1].Up to now, many methods have been utilized to measure the SSD in optical glasses, and they are generally classified into destructive and nondestructive evaluation methods.The chemical etching method, one of the typical destructive evaluation methods, is a simple, low cost and intuitive way to measure SSD depth.However, it has many drawbacks in utilization, such as difficulties of controlling etching process and partial information of the SSD.For the nondestructive evaluation methods, the main feature is that a sample will not be damaged.But with these methods, the detailed analysis can-not be given [2,3].To characterize SSD by indirect measurement of surface morphology, P.E.Miller had suggested the empirical and semi-empirical correlations to allow one to estimate the depth of SSD [4].Lambropoulos's study shows that abrasive processes has an important effect on the depth of SSD and gives us an empirical formula [5].Wesley B, et al, proposed a method of using quantum dots to tag SSD in lapped and polished glass samples, which is a quick way to detect SSD, and their work takes a new insight into how material is removed during lapping and polishing processes [6].Marcus Trost and Tobias Herffurth evaluated SSD by the light scattering techniques, and their research shows that the light scattering measurements has an advantage in fast, nondestructive, accurate and detailed measurement [7].
Many efforts have been made to characterize the SSD precisely and fully.Tayyab Suratwala had provided several new rules of thumb, such as the relationships between length and lap properties, between width and depth, between width and rogue particle size [8].A.Esmaeilzare provided the relationship between surface roughness of the ground surface which has maximum un-deformed chip thickness and the distribution of micro-cracks in the ground surface.He gave us the parameters of SSD depth, SR (surface roughness) values (pv), length of all cracks, and dispersion of damages in various depths using angle polishing technique [9].However, due to the micro scale and position feature of SSD, it is still difficult to get a detailed and precise characterization of the SSD of high optical components.
One method has been proposed to obtain the 3D (threedimensional) information of the micro cracks of SSD with FD (finite difference) algorithm in our previous work [10].In this method, accurate description of 3D evolution of the etching process is significant for the calculation.In this paper, we calculate the etching surface morphology at a certain time by means of the five-point FD algorithm model .After analyzing the 3D morphology of micro cracks at the specific time, parameters of the SSD including the width, length, can be obtained.Through comparing calculated SSD parameters with direct measured ones, some useful conclusions can be obtained.

EXPERIMENT AND METHODOLOGY
To get a complete characterization of the SSD, the five-point finite different algorithm model is programmed, which starts from the morphology of a sample with totally opened SSD cracks gained by etching process.The FD algorithm can be used to deduce the morphology of SSD at a certain time by simulating the etching process reversely.When the morphology at a given time is obtained by the five-point finite different algorithm model, we can stipulate some rules of extracting parameters to characterize it.The parameters not only refer to the traditional depth of SSD, but also include the width, length, distribution and more features of the SSD.To verify the validity of this model, the observed result during the experiment was used to show the real morphology of the sample at different etched time.The following part will give a description of the experiment and the FD algorithm.

Experiment
The K9 glass samplesK9: SiO 2 =69.13%B 2 O 3 =10.75%BaO=3.07%Na 2 O=10.40%K 2 O=6.29%As 2 O 3 =0.36%with size of φ10×2 mm and marked with a circle φ5 mm by lasing was made 9×9 humanoid cracks with the same lengthdepth artificially by micro indentation with standard Vickers tip at load of 0.0981 N, on the micro hardness tester (HXD-1000TMC).The Vickers hardness indenter is triangular pyramid diamond indenter with tip nominal radius of curvature of 50 nm.Then, the samples were etched by BOE solution (buffered oxide etching; 40% NH 4 F, 49% HF ratio of 12:1).During the etching process, the morphology of the sample's surface at different time were measured by the confocal laser scanning microscope (OLS4000).The measured results provided data for the corrosion model.A FD algorithm, according to reference [11], was programmed to simulate the etching process.

Methodology
The finite difference method uses the corresponding variable discrete values to replace the continuous values as independent variables in the differential equations.The surface function which was used to describe the morphology is shown in Eq. ( 1). S = S(x, y, t) Where x, y are the coordinate of the plane, S is the surface coordinate at time t and at coordinate (x, y) and is a function of time t, the coordinate x and y.
The morphology at a certain time can be gained through the five-point finite different algorithm as long as the morphology data at the previous moment is known.The formula of computation was shown as Eq. ( 2) Where S(x i , y j , t) is the surface height at time t and at coordinate (x i , y j ), S(x i , y j , t + t) is the surface height at time t + t and at coordinate (x i , y j ), r b is the corrosion rate, n z represents the direction of growth at the point (x i , y j , t).
In order to discretize values of the independent variable x, y, the entire surface is divided into equal steps 1024×1024 grids network in parallel to the direction of x and y axis, according to the image taken by the laser scanning confocal microscope (OLS4000).With the 3D morphology, we can examine cracks of surface and find out important parameters to characterize the cracks such as depth, length, width, density and so on.If the morphology data is fully gained, the SSD of the sample can also be calculated.
As shown in Figure 1, the procedure of deducing morphology of surface after a required etching time from the measured raw data was performed using the following steps.
Step 1: The data of original measurement is loaded into MAT-LAB and drawn into the three-dimensional graphics in order to be observed clearly.
Step 2: The effective data is chosen.The etching rate rb and etching time ∆t should be set.According to the differential corrosion of sample, r b is about 0.088 µm/min.
Step 3: The five-point difference algorithm is established to get the coordinate values of every point during the entire etching process.In order to set a five-point difference algorithm, a none-boundary point should be selected as the center point at first.Then the nearest four points to the center point are found.The five points will form a rectangular pyramid and the four cam faces normal vector of the pyramid can be obtained and be summed up to calculate the growth direction.After all of those, the growth direction, the etching rate r b and etching time ∆t are used to calculate the position at a certain time.
Step 4: Considering the enhancement of etching rate on the convex surface in our model, the correction vector Q(i, j) is added.
Step 5: In order to turn these data back to the original uniform system, the five points should be kept as a unit to calculate the nearest plane to the original point.Then the coordinate in Z direction on this plane is got.
Step 6: According to the internal time during the etching process and the etching step which are set by us, the iterations can be determined.If cycle times is less than or equal to the number of iteration, the program will turn back to step 3. Else the data will be output and drawn into 3D graphics.FIG. 2 The observed morphology and the calculated morphology at different moment.

Experiment results and alignment work
After etching and observation by the confocal laser scanning microscope, the data of the morphology of the sample's surface at different time could be gained.The data was then put into MATLAB and the morphology at specific moment would be calculated by the five-point finite different algorithm.The observed morphology and the calculated ones were put together, as shown in Figure 2 The alignment work should be done first before making comparison of the experiment and calculated morphology in the vertical and horizontal directions.The least square method was used to find the four experimental cracks corresponding to those in the calculated cracks.
Where i, j are the coordinate of X and Y, m and n are from 0 to a and b, S 0 (i + m, j + n) is the data of experiment, S 1 (i, j) is the data of calculation.The matched cracks were found when the min c was the minimum value.We accounted that the error of calculated crack was the smallest at this position.The distances of corrosion evolution morphology moving along the X − Y direction are m and n.
After finding the coincident four cracks, the morphology of experimental and calculated at the time of 10 min, 20 min, 30 min were put into the same coordinate system respectively, as shown in Figure 3.
The calculated morphology was lower than experimental morphology at the same moment, because the entire surface was also etched during the etching process.To solve the problem, the datum plane of experimental and calculated morphology, namely the highest plane of none cracks, should be calculated.The averaged value is the datum value.
Where h is the height value of the datum, S(i, j) is the coordinate values of topography, g is the number of points.The difference between the datum planes (∆h) should be taken into consideration Where h 1 is the value of experimental datum plane, h 2 is the value of calculated datum plane.
After the alignment of the image in horizontal and vertical directions, as shown in Figure 4 (the aligned picture at 20 min was taken as an example), the next step is to compare the experimental with the calculated data to evaluate the length, depth, and the sum of Residual Square of the subsurface damage.

Comparison and analysis
In this paper, the etching process of a sample with SSD had been simulated reversely.At first, the cracks were opened quickly.Then, the cracks fused together due to the little space between each other until the etched surface became flat.In order to improve the precision of calculation, the time interval should be compared to gain the appropriate value.The calculated images at 20 min with different time intervals were shown as Figure 5.The smoothness of the image and the computational difficulty were two of the factors to decide the appropriate value of time interval.
In the process of calculation, each condition with different corrosion time interval was calculated respectively.The results were gained and then were compared with the original data to get the value of error.The values of error with different time intervals were compared to find the appropriate time interval to ensure the error precision in the future calculation.
Results at different moment were compared with the original image, taking the variation curve error value into consideration.In the assessment of the error value, the location of the length-calculation is in the center of crack and parallel with x axis and y axis.The position of the depth-calculation is consistent with the calculation of length.If the five-point finite different algorithm model is used properly, the information of SSD at every moment can be inferred accurately.
In the model to evaluate the parameters' accuracy, the parameters will be calculated including depth, length and the residual sum of square about SSD.
Calculation for the error of depth: First, in the comparison process, a value, p, can be defined as a desired error.The value of p can be set as 10% of the value of depth of crack D or a specific numerical value 0.1 µm to get the different error index value.∆D(i, j) = S i (i, j) − S 0 (i, j) Where S 0 (i, j) presents the value of crack depth at calculated condition, S 1 (i, j) means the value of crack depth at experimental condition, ∆D is the absolute difference between the crack depth .
In the calculation process, if the value ∆D > p, the point is put into points of error, then s is added by 1 (s = s + 1).and the value of ∆D is recorded and accumulated.The sum of absolute difference between the crack depth ∆D is the depth of error value about the experiment data and is defined as D 0 .
Where q is the number of error point.
Then calculate the error of the length until the end of the cycle.
Where L e is crack length in corrosion evolution model, L c is the crack length under experimental conditions, ∆L is the absolute length error of crack.Calculation of the residual sum of squares (RSS): In the calculation, the residual sum of squares of data can reflect the fitting degree of the calculated and experiment data.
The smaller the RSS is, the better the two sets of data match.
As shown in Figure 6, the value of error increases with the increase of time interval.In the 10 and 30 min, the error value of length ranges from 0.25 µm to 1.49 µm.The residual sum of squares varies in a wide range as it has the max value of 3.51 µm 2 at 10 min and 38.46 µm 2 at 30 min, which is shown in Figure 7.The Figure 8 and Figure 9 show that the percentages of error about length and depth also increase with the increase of the time interval.In 10 min, the minimum and maximum percentages of the error of length are 6.07%and 15.18%.
Through the comparison, the results of calculated morphology were consistent with those of experimental morphology basically, and the errors between them were so little that they would not have big influence on the accuracy of the simulation.
There are several factors which will influence the accuracy in this study: FIG. 8 The percentage error of length at 10 and 30 min.
FIG. 9 The percentage error of depth at 10 and 30 min.
1.The measurement accuracy of the confocal laser scanning microscope is only 0.12 µm in X direction in the measurement.
2. The round-off error in the process of iterations will induce the error in the calculated model.3.There are some conditions that cannot be controlled easily in the corrosion, like the impurity of solution, the hardness of controlling the very time span of corrosion, and so on.

CONCLUSIONS
The FD algorithm is a useful tool to simulate the etching process.The measurement results of SSD are consistent with the preconceived values of this model.With the 3D morphology simulated by the model, the cracks of surface and their important parameters can be obtained to characterize SSD, such as depth, length, width, density, and so on.Based on the application of the finite different algorithm and the experimental research presented in this paper, the following conclusions can be drawn: 1.The accuracy of the model proposed in this paper can be got by comparing experimental with calculated results.Through analyzing the result of this research, the error is found to be very little.
2. The corrosion model based on the finite difference five point algorithm was established.The model takes into account the influence of time interval on the accuracy of the algorithm.By comparing the value of error at different time interval, the appropriate time interval can be got to make the calculation results closer to the real measurement.The value of error between the measured and calculated data increased with the time interval, so in order to keep the percentage of the error of length and depth approximately stable at 10%, the time interval of 0.2 min can be chosen in the process of calculation.
In future work, more research should be taken to decrease the error of the 3D finite differential model.The parameters calculated by this model can provided a theoretical basis to optimize the mechanical process.Besides, this work lays a foundation for the future reverse calculation.

FIG. 1
FIG.1The procedures of calculating model by using the five-point finite different algorithm.

FIG. 3
FIG.3The morphology evolution of sample calculated by FD algorithm at different etching time.

FIG. 5
FIG.5The images of different time intervals at 20 min.

FIG. 6
FIG.6The error of length and depth at 10 and 30 min.