A geometric prediction model of surface morphology in micro-EDM considering stochastic characteristics of discharge crater size

Accurate prediction of surface morphology is closely related to the application range of micro-EDM. Most of the existing prediction models do not consider the stochastic characteristics and uncertainty of the discharge crater, which can produce great differences with the experimental results. Aiming at these problems, the stochastic characteristic of the discharge crater size was studied, and a simulation model was established to predict the surface topography of EDM considering the stochastic characteristics of crater size. Firstly, the stochastic characteristics of the crater size were investigated through the finite-successive pulses discharge method and characterized from the perspective of probability theory. The stochastic characteristics were explained by analyzing discharge waveforms. Secondly, the response surface method (RSM) was used to map the corresponding relationship between the input machining parameters and the distribution characteristics of crater size. It was observed that the RSM model can accurately predict the distribution characteristics of crater size in the range of parameters selected in this study. Finally, models based on stochastic distribution of crater size (MSDS) and fixed crater size (MFS) were established respectively. By analyzing the characteristics of the simulated surfaces and the machined surface, it was found that the surface features of MSDS are closer to the real surface features. Compared with MFS, the roughness prediction error and the average of overall error of MSDS are reduced by 1.01% and 18.97%, respectively. The results of this work are helpful to understand the randomness and realize the controllability of EDM.


Introduction
Micro-electrical discharge machining (micro-EDM) is a noncontact machining technology which uses electro-thermal erosion to remove materials [1]. The machining mechanism makes it very advantageous in difficult-to-cut materials such as titanium alloy and cemented carbide [2]. It has been widely used in aerospace and medical fields due to its special advantages [3][4][5]. The EDM surface consists of overlapping discharge craters, creating a special microstructure with pits and protrusions. The microstructure can be used as a good basis for fabricating superhydrophobic surface, which is beneficial to improve the water repellent of the surface. Another application of EDM technology is to modify the surface of artificial joint and other implants, because the special morphology of EDM surface is helpful to promote cell growth. However, whether superhydrophobic surface or biocompatible surface, there are certain size requirements for surface morphology. Therefore, precise control of surface morphology is of great significance for the application of EDM technology in the manufacture of superhydrophobic surfaces, biocompatible surfaces, and other functional surfaces.
Simulation model of micro-EDM is an effective method to predict and control the morphology of micro-EDM surface. Generally, the reports on simulation model can be classified into two categories. One is based on single spark. Kiran and Joshi [6] developed an evaluation model of surface roughness based on single spark crater, and their model generated quite remarkable prediction of the surface roughness in the experiment. Salonitis et al. [7] developed a thermal based model to predicate the material removal rate and surface roughness and deduced the average surface roughness when two or more discharge craters overlapped based on the single discharge crater obtained through the model. Zhang et al. [8] observed through microscope that the discharge crater depth is about 3 times that of the average surface roughness, and their roughness model based on individual craters is in good agreement with the actual surface roughness to a certain extent. The other is based on multi-spark discharges. Tan and Yeo [9] used a symmetric distribution of discharge craters to estimate the upper and lower bounds of R max and found that the prediction of the model is more accurate in the case of larger discharge energy. Kurnia et al. [10] proposed a geometric model to estimate the surface roughness based on the uniformly distributed craters, with the average error of R a and R t being 6.5% and 4.5% separately compared with the experimental results. However, in the actual EDM process, the location distribution of craters is stochastic, which must be considered to improve the prediction accuracy of the model. Wang et al. [11] developed a multi-spark model to study the machining process of ultrasonic vibration-magnetic field assisted WEDM, the discharge located at the point where the gap width is shortest between the wire electrode and workpiece. Izquierdo et al. [12] have taken the random distribution of discharge positions into account, and the adjusted simulation results show that the errors are below 10% compared with the measured results. Liu and Guo [13] developed a multi-spark model for die-sinking EDM considering the discharge positions are randomly distributed, to simulate the evolutionary process of EDM surface morphology. However, it is not enough to consider only the stochastic nature of discharge positions; in another work, Jithin et al. [14,15] assumed that the multiple sparks are randomly generated with respect to discharge location and discharge energy, and the distribution of surface roughness shows agreement with that of experimental results under low pulse width, while the error increases under larger pulse width. However, the crater size has great randomness, and the existing models still have some gaps for improvement due to the lack of the assumptions of crater size which is closer to the real situation; furthermore, these models tend to focus on the morphology characteristics perpendicular to the surface (such as surface roughness) but pay little attention to the morphology characteristics parallel to the surface.
In EDM, change of machining parameters can change the discharge crater size, hence the surface morphology. To better control the surface morphology, scholars have done many researches on the factors affecting the discharge crater size. Li and Bai found that the discharge crater became smaller as the gap width increases in the case of the same discharge energy [16]. Mohammadreza et al. [17] compared the crater sizes obtained by thermal simulations and experiments to analyze the influence of the machining parameters on the flushing efficiency of plasma tunnel; they found that flushing efficiency of plasma tunnel increases with the increase of discharge current and decreases with the increase of pulse width. Carlos and Patricio [18] found that the crater size increases with the increase of discharge current, and they take the diameter of the crater as a good index to characterize the surface morphology of parts in die-sinking EDM. Liu et al. [19] investigated the effect of pulse width on the discharge crater diameter. It can be concluded from the above literature that the gap width, pulse width, and current have a significant impact on the size of the discharge crater, so the experiment will be designed with these three parameters as factors in this paper; furthermore, the average value of discharge crater size is taken as the statistical result, but this statistical method cannot reflect the stochastic characteristics of EDM process, because even under the same machining parameters, the crater size also has great differences [18,20,21]. Therefore, a more effective statistical method is needed to describe the size discrepancy of discharge craters.
In this paper, firstly, the stochastic distribution characteristics of crater size are calculated and described by the probability density function. Secondly, the influence of peak current, pulse width, and gap width on crater size distribution is studied by response surface design (RSM), and the regression equation for predicting crater size is established. Finally, the geometric model of EDM surface topography is established, and the crater size described by probability density function is applied to the model. The effectiveness of the proposed simulation model is verified by comparing with ordinary simulation model and the experimental results.

Experimental equipment
The experiments were conducted with a self-made micro-EDM machine tool as shown in Fig. 1a. The resolution of the machine tool is 0.1μm in X, Y, and Z directions. The maximum rotational velocity of the spindle is 30 000rpm, of which the radial runout is less than 1 μm. The tool electrode manufacturing module and vision inspection module based on the CCD camera are integrated into the machine tool, enabling the tool electrode to be fabricated and measured easily. In the experiment, a self-made impulse generator controlled by onechip computer was used, the power supply could generate any number of pulses, the circuit diagram for the experiment was illustrated in Fig. 1b, and the finite-successive pulses discharge method was used.

Electrode and workpiece
The diameter of the cylindrical tool electrode was 500 μm with tungsten as material, and the material for the workpiece was Ti-6Al-4V. The tip of the tool electrode was fine ground by wire electrical discharge grinding (WEDG) method, more details on the fabrication process can be found in our previous research [16]. The surface of the workpiece was ground with 2000 grit emery paper to ensure the original surface roughness had no influence on the experimental results.

Experimental procedures
The aim of this section is to investigate the stochastic characteristics of discharge crater size. The three machining parameters, including peak current (i p ), pulse duration (t on ), and discharge gap width (d), were taken as factors in the design of the experiments. The RSM was used to analyze the influence of the three factors on crater size. The discharge parameters could be adjusted easily; however the determination of the gap width between the electrode and the workpiece is a bit more difficult; hence, a measurement method based on ohmmeter was used. Firstly, the one ends of the tool electrode and workpiece were connected by metal wire, and the other ends were connected by the two pins of an ohmmeter to measure the reference resistance R g as shown in Fig. 2a; then the metal wire was removed and the tool electrode moved down 1 μm every step as shown in Fig. 2b; when the resistance value displayed on the ohmmeter was equal to R g , the gap width was 0 μm as shown in Fig. 2c; the tool electrode moved up to a certain distance d afterwards, and the target gap width was obtained and the experiment can be carried out. The top view of the discharge craters was observed with an ultra-depth 3D microscopy system. The diameter of each crater was measured manually for 3 times with the ruler of the microscopy system, and the reported diameter was the average value. The 3D morphology of the crater was obtained by probe scanning of an AFM automatically, and then the depth information of discharge crater was obtained. Fig. 2 e shows the top view of ultra-depth microscope image of discharge craters; it can be seen that the discharge crater is almost circular or elliptical; hence it is characterized by diameter (for circular crater) and equivalent diameter (the average value of major axis and minor axis for elliptical crater).
In this study, the finite-successive pulses discharge method was used, and the number of pulses is set to 9. One hundred discharge craters were obtained under a parameter combination through the repeated experiments, of which the equivalent diameters were selected to characterize the size distribution. The experimental parameters are shown in Table 1, and the variable parameters and their levels are listed in Table 2. The experiments were designed based on a central composite design method, and the experimental sequence is shown in Table 3.

Prediction model of discharge crater size
The response surface method (RSM) is a statistical test method to optimize stochastic processes. The objective is to find the quantitative law between the test index and each factor and find the best combination at each factor level. When the machining parameters are determined, the size characteristics of discharge craters can be predicted by fitting mathematical models as follows: where Y(o) is the target response; x i is coded values of the i th machining process parameters; a, b i , c i , and d ij are the regression coefficients; e is the corresponding experimental error of observation; and n is the number of variable parameters.

Characterization of the stochastic characteristic of discharge crater size
The main types of discharge that produce discharge craters are shown in Fig. 3. Fig. 3 a shows the normal spark discharge; the discharge process continues until the end of the pulse. Differently, there are two kinds of spark discharges as shown in Fig. 3 b and c: one is the spark with interruption, i.e. the discharge process is terminated before the end of the pulse; the other is transient spark, i.e. the discharge channel extinguishes almost immediately after it is formed. Fig. 3 d shows the discharge with a very short ignition delay and without ignition delay, namely transition arcing and arcing respectively. It can be seen from Fig. 3 that the processing parameters have significant effects on the discharge waveform. To clarify the influence mechanism of these parameters on the discharge waveform, we analyzed the discharge maintenance process from the perspective of the discharge channel. Based on Townsend discharge theory, the discharge maintenance can be described as below [16]: where γ is the secondary electron emission probability caused by the collision of positive ions with the cathode; α is the electron collision ionization coefficient; and d is the discharge gap width. γ and α increase with the increase of electric field strength. (e αd -1) is the number of positive ions formed by the α process during the emission of an electron from the cathode to disappearance at the anode. γ (e αd -1) is the number of electrons emitted by the γ process. The spark process as shown in Fig.3a can be explained as follows: after the discharge channel is formed, the maintenance of the discharge channel depends on the secondary electron emission until the end of the pulse. It is worth noting that in the discharge circuit as shown in Fig. 1b, the current limiting resistor (slide rheostat), DC power supply, and discharge gap are connected in series, and the open circuit voltage (is constant) is provided by DC power supply, and the peak current is controlled by adjusting the resistance of the current limiting resistor. The output voltage of DC power supply is distributed to both ends of current limiting resistor and discharge gap. Therefore, the current limiting resistor has a negative feedback effect on the discharge gap voltage, and the larger the resistance value, the more significant the negative feedback effect [22]. After the gap is broken down, the current in the discharge circuit increases suddenly, which leads to the increase of the voltage at both ends of the current limiting resistor. Accordingly, the voltage at both ends of the gap decreases. So, when the peak current drops to 0.5A as shown in Fig. 3b, the negative feedback effect of current limiting resistor is stronger, which leads to the excessive drop of gap voltage, the excessive decrease of gap electric field strength, and the weakness of γ process. Hence the secondary electron emission fails, and the discharge channel is easier to be extinguished, and that is the reason why so many transient sparks exist. In contrast, with the parameters shown in Fig. 3a, the negative feedback effect of current limiting resistor is weak, and the discharge channel can be maintained. Compared with Fig. 3a, only the gap width increases in Fig. 3c, which leads to the decrease of gap electric field strength, and the processes of α and γ are weakened; as a result, the maintenance capability of discharge channel is decreased, and it is easy to be extinguished on this condition; hence the discharge duration is very short. As shown in Fig. 3d, the electric field intensity will increase when the gap width decreases, resulting in the enhancement of α and γ processes; on the other hand, the decrease of gap width will also decrease the insulation strength of the gap, which makes the gap easier to be broken down. Therefore, the difficulty of gap breakdown is reduced, and the ignition delay of spark is reduced or even disappeared, and the arcing is formed. Although the peak current in Fig. 3d is smaller than that in Fig. 3a, there is no interruption in the discharge process, which indicates that the negative feedback effect of the current limiting resistor is not enough to weaken the maintenance capability of the discharge channel in such circumstances. Fig. 4 a shows the discharge craters, and the diameters of the craters were analyzed with Minitab software. The histogram of diameter size distribution for experiment #2 is shown in Fig. 4b, it can be found that the histogram is inverted bell shaped and roughly symmetrical, hence we intend to use the normal probability density function to fit the distribution of the histogram, and the goodness-of-fit test for the normal distribution was carried out. The probability map of the fit goodness test is shown in Fig. 4c, the median line is the expected percentile of the distribution based on maximum likelihood parameter estimation, and the upper and lower side lines represent the lower and upper limits of the confidence interval for each percentile; the probability map shows that the data points are approximately a straight line and within the confidence interval. The data of group #2 conform to the normal distribution, because P-value satisfies Eq. (3): where P-value (probability value) is the least significant level at which the original hypothesis can be rejected, which is obtained from the sample observations of the test statistic, and α is the minimum significance level for rejecting the original hypothesis, and the data can be considered accord with the normal distribution when P>α.  Table 3), and except for the 8th group of experimental data, the P-values of the other experimental groups are all greater than 0.05, which indicates that the size of the discharge crater obeys normal distribution except for exp. no. #8. In other words, when the normal probability density function is used to fit the size distribution of the discharge crater, only 5% of the experimental results have large errors, while the remaining 95% can be well described using the normal distribution. Therefore, under the condition of accepting a 5% error, the normal distribution function can be used to characterize the size distribution of discharge craters.
The general equation of normal distribution curve is described as [23]: where μ is the mean value of the sample; σ is the standard deviation of the sample; x is the random variable, i.e., the diameter of discharge crater in the study; and f(x) is the probability density of x. It can be obtained from Eq. (4) that the shape of normal distribution curve is controlled by μ and σ.

Prediction model for size distribution of discharge crater
Through the analysis of the previous section, it was found that the stochastic characteristics of discharge crater size could be described by the normal probability density function; hence the mean value and standard deviation of the crater size were used to quantitatively describe the stochastic distribution characteristics of crater size and were taken as the response of RSM. Fig. 6 depicts the output response of the RSM. Across the results, it was obtained that the mean value and standard deviation of crater size are smaller under 3#, 10#, 16#, 18#, and 20# parameter combinations, in which the current levels are low or gap width levels are high, hence discharge current and gap width have a significant influence on the crater size, and the influence of pulse width on crater size is relatively small.
A second order polynomial is used to clarify the relationship between machining parameters and the size distribution of discharge craters, as shown in Eqs. (5) and (6). where Y 1 is the mean value; Y 2 is standard deviation; and x 1 , x 2 , and x 3 are pulse width, peak current, and discharge gap width, respectively.
In order to verify the validity of the prediction model given in Eqs. (5) and (6), the test experiment was conducted as shown in Table 4. The results of the test experiment and the prediction are shown in Fig. 7, where MV1 and MV2 are the mean values of crater diameter obtained from the prediction model and the test experiment, respectively, and SD1 and SD2 are the standard deviation of crater diameter obtained from the prediction model and the test experiment, respectively. The maximum errors for the predictive values of mean value and standard deviation are 5.69% and 6.71%, respectively, and the average errors for the two are 4.29% and 5.49%, respectively. Thus, the prediction errors are acceptable, indicating that the prediction model is effective.

Effects of machining parameters on the discharge crater size distribution
The crater size is related to the discharge energy distributed to the surface of the workpiece, which can be calculated as [11]: where E W is the discharge energy distributed to the surface of the workpiece, η is the energy distribution coefficient, E is the discharge energy, and u, i, and t are average discharge voltage, average discharge current, and discharge duration, respectively. Generally, u is regarded as a constant, and i increases with the increase of peak current during discharge process. The larger the E W is, the larger the crater size is. Therefore, the discharge energy can be described by the crater size. The effects of pulse width, peak current, and discharge gap width on the response, i.e., mean value and standard deviation of discharge crater diameter, have been investigated with contour map obtained by RSM as follows: (1) Effects of parameters on the mean value of discharge crater diameter Fig. 8 shows the interaction of pulse width and peak current on the mean value. It is obtained that the mean value increases with the increase of pulse width and peak current and that the sensitivity of mean value to peak current increases under the high level of pulse width, that is because the discharge energy in a single pulse increases with the increase of peak current, which results in an increase of the discharge crater size. However, the pulse width has little effect on the average size of the discharge crater at low peak current level, this is due to the maintenance ability of the discharge channel decreases when the peak current is low, and the discharge channel is easily extinguished while the pulse is still on, as shown in Fig. 3b. Under this circumstance, the discharge duration t and discharge energy will not increase with the increase of pulse width. Under the low level of pulse width, the mean value increases slightly with the increase of peak current,   Standard deviation (μm) Fig. 7 Experimental and predictive results for crater size and the sensitivity of mean value to peak current is very small, that is because the small pulse width limits the increase of discharge energy E in a single pulse. Fig. 9 shows the effect of the interaction of pulse width and gap width on the mean value. It has been observed that the mean value is minimum with the combination of high-level gap width and low-level pulse width, and the mean value is maximum with the combination of low-level gap width and high-level pulse width. The influence of gap width on mean value can be considered from two aspects: firstly, with the increase of the gap width, the energy distribution coefficient η decreases [24], and the discharge channel expands and the heat flux decreases [25], which will cause the decrease of E W ; secondly, the increase of the gap width will lead to the instability of the discharge channel and the weakening of the discharge maintenance ability; according to the comparison between Fig. 3a and Fig. 3b, increasing the gap width will make the discharge channel easier to extinguish within the pulse duration when other conditions remain unchanged.
(2) Effects of parameters on the standard deviation of discharge crater diameter Standard deviation is a measure of the dispersion for a group of data. A large standard deviation represents a large difference between most discharge crater sizes and their mean value, and a smaller standard deviation indicates that these sizes are closer to the mean value. Fig. 10 shows the influences of pulse width and peak current on standard deviation; the sensitivity of standard deviation to pulse width is low, but it is high to peak current. The standard deviation increases with the increase of peak current, that is because the discharge power increases with the increase of peak current, and the discharge duration t is random even if the processing parameters are the same; therefore, when i increases, a small change of t will lead to a significant change of E as shown in Eq. (7), which leads to the increase of the size difference of discharge craters. On the other hand, according to Fig. 3, the relationship between discharge duration t and pulse width t on can be obtained as follows: where t = 0 means that it is an open circuit pulse and t = t on means that there is no ignition delay and the discharge channel does not distinguish during the pulse duration. When t on increases, the value range of t also increases, which leads to the further increase of the dispersion of discharge energy E. Therefore, the maximum standard deviation appears in the combination of high-level current and high-level pulse width. It can be seen from Fig. 11 that the influence of gap width on standard deviation increases with the increase of pulse width; this is due to the discharge duration t tends to be more consistent under the low-level pulse width and gap width, and with the increase of gap width, the maintenance ability of the discharge channel decreases, and it is more likely to be  Fig. 10 Contour plot for standard deviation vs i p and t p extinguished, which will increase the dispersion of the crater size. At the circumstance of high-level pulse width, the influence of gap width on standard deviation increases, which is also related to the maintenance ability of the discharge channel. The maintenance ability of the discharge channel decreases with the increase of discharge gap width, which leads to a greater dispersion of discharge duration t and the discharge energy E in a single pulse.

Simulation model for micro-EDM
EDM is an important method for manufacturing superhydrophobic surfaces, surface enhanced Raman scattering, and other functional surfaces. Surface morphology plays an important role in these functions. In order to accurately manufacture these functional surfaces, it is necessary to predict the surface morphology; hence a simulation model of the EDM process was developed. In this model, the tool and workpiece can be divided into two-dimensional meshes [26,27]: the electrode and workpiece are divided into squares of the same size as shown in Fig. 12.
The model proposed in this study is based on the following assumptions: (1) There is only one discharge that occurs in every time step, and the time step is described as Eq. (9).
where Δt is the time step and f is the pulse frequency.
(2) Only spark discharge is considered, and the spark discharge occurs in the shortest gap width between the two points of the electrode and the workpiece.
The stochastic distribution of the crater size is considered, and size distribution characteristics are obtained according to the experiment. (4) The diameter of the plasma channel is not considered, and the discharge energy is delivered point to point. (5) The shape of the discharge crater is cone-shaped, and the discharge debris is ignored.
The model includes horizontal feed and vertical feed, and the flow chart of the simulation process is shown in Fig. 13.
The electrode and workpiece are represented by black square pixels, the gap is represented by white pixels, and each pixel has three indexes, i.e., abscissa, ordinate, and gray values. The positions of the electrode and workpiece are determined by the abscissa and ordinate, and the contour of the electrode and workpiece is judged by the gray value. The workpiece and the electrode are in the same set of absolute coordinates. In Zhang's model, a fixed size of the discharge crater is used to simulate the EDM process so as to simplify the calculation. This method can obtain certain details such as material removal rate, electrode wear rate, and electrode shape, but it is difficult to simulate the surface morphology of EDM because in actual machining, the size of the discharge craters varies greatly even if using the same parameter combination, as mentioned above. Therefore, assumption (3) can in theory make the model closer to the actual EDM surface morphology in theory.

Searching for discharge position
In order to simplify the calculation, the lower boundary of the electrode and the upper boundary of the workpiece are

Feed direction Electrode
Workpiece Grid Fig. 12 The 2D simulation model of EDM milling process determined first, and the point pair from the two boundaries separately with the minimum distance is searched then. When there are many point pairs that meet the requirements, one of which is randomly selected as the discharge position. The determination method of discharge position is calculated according to Eq. (10) [26]: where d is the distance between the two points of the electrode and workpiece, (x i,e , y i,e ) is the coordinates of any point on the electrode contour, (x j,w , y j,w ) is the coordinates of any point on the workpiece contour, and D is the set discharge gap width.

Material removal
The actual crater (shown in Fig. 14a) is shaped like a spherical cap as shown in Fig. 14b, and it can be simplified to a geometric model as shown in Fig. 14c, where R is the diameter of the bulge and h is the depth between the top and bottom of the crater. Ten randomly selected craters were measured by AFM, and the average value of R/h is about 10.
Due to the discharge crater on the electrode being not easy to measure, the discharge crater on the electrode can be obtained according to the Eq. (11) [26]: where r and r′ are the crater diameters on electrode and workpiece separately and θ is the electrode wear rate, and it is obtained by experiment. In this model, the size of the crater diameter is normally distributed, and the normal probability density function is divided into 11 rectangles according to the Riemann integral in order to generate the discharge crater with normal distribution probability. The ordinate of each rectangle represents the probability of its occurrence, and the abscissa represents the crater diameter under this probability. The flow chart for judging discharge crater size is shown in Fig. 15.

The forming of EDM surface morphology using simulation
In order to demonstrate the formation of EDM surface, the simulation research is carried out, and the simulation parameters are shown in Table 5.
In the existing geometry simulation, the discharge crater is set to a fixed size; in order to investigate the influence of normal distributed crater size on the simulation results, two models, the model based on stochastic distribution of crater size (MSDS) and the model based on fixed crater size (MFS),  were developed. The simulation processes of the two models are shown in Fig. 16 and Fig. 17 separately. In the simulation, the first discharge location is stochastic, and the first overlap occurs at the 9th spark in MSDS, and the first overlap occurs at the 6th spark in MFS; with the increase of the number of sparks, the number of overlapping craters gradually increases until the surface of the workpiece is completely covered by craters as shown in Fig. 16c and Fig. 17c. Because the size of the electrode and workpiece is much larger than the crater size, some details cannot be shown in Fig. 16 and Fig. 17, hence the EDM surface morphology data for Fig. 16c and Fig. 17c were extracted, and the broken line diagram was drawn in Fig. 16d and Fig. 17d separately. It is found that the simulated surface of the two models is significantly different by comparing Fig. 16d and Fig. 17d, the degree of the surface undulation of MSDS is larger compared with that of MFS, and the simulated surface of MFS is regularly wavy; this is because the crater size in MSDS shows considerable differences, and the difference in depth of crater causes the fluctuation of surface height. Mohammad et al. [28] thought that the traditional roughness (Ra) could not accurately describe the surface morphology characteristics of EDM. They machined two groups of surfaces with the same Ra, but with different characteristics of discharge craters. This is because Ra primarily focuses on the characteristics of the morphology in the direction perpendicular to the surface, but does not pay attention to the morphological characteristics parallel to the surface. Currently, The mean value and standard deviation of the diameter of the discharge pit.

Determine the normal probability function y = f (x)
The normal function f (x) region is divided into 11 parts  Therefore, in order to describe the characteristics of EDM surface morphology more accurately and compare the applicability of the two surface morphology prediction models, the surface height and adjacent peaks interval distribution are used to describe the surface morphology features of the two surface morphology prediction models. As shown in Fig. 18, in the schematic diagram of cross section morphology, the height of the lowest point A is defined as 0, and L i and H i are the interval between adjacent peaks and the height of the data point, respectively. It should be noted that using this calculation method, the mean value of surface height is equal to the roughness Ra. Fig. 19 shows the surface height and peak interval of the simulated surfaces. As shown in Fig. 19 a and b, it can be seen that the surface height and the distance of adjacent peaks are normally distributed on the simulated surface of MSDS. As for MFS, the surface height and the distance of adjacent peaks are irregularly distributed as shown in Fig. 19 c and d. From  Fig. 19c, it can be observed that the surface height of MFS is relatively concentrated, which indicates that the surface is relatively flat. Fig. 19 d shows that the surface peak interval is about 25 μm, compared with Fig. 19b; it can be seen that the simulated surface of MFS is more regular that of MSDS; this feature can also be obtained by comparing Fig. 16d and Fig.  17d.

Validation experiment for the model
In order to verify the effectiveness of the surface morphology prediction model, a verification experiment was carried out, and the machining parameter combination adopted is #2. The experimental result is shown in Fig. 20a, and the data of surface morphology were obtained by laser confocal microscope. Because the EDM surface is isotropic, the surface profile data of y = 300 μ m section were selected for analysis as shown in Fig. 20b. The surface height and adjacent peaks' distance characteristics of the machined surface are shown in Fig. 21, and the surface height and peak interval are normal distributions. By comparing Fig. 19 and Fig. 21, it can be seen that the simulated surface of MSDS has some characteristics similar to that of the machined surface, while the simulated surface of MFS does not have these characteristics. It is observed from Fig. 20b that the fluctuation of the machined surface is obvious, and compared with Fig. 16d and Fig. 17d, we can conclude that the simulated result of MSDS is closer to the experimental result than that of MFS. H 1 H i Fig. 18 The schematic diagram of cross section morphology As shown in Fig. 22a, MV′1, MV′2, and MV′3 are the mean values of surface height for the results of MSDS, MFS, and experiment separately; SD′1, SD′2, and SD′3 are the standard deviation of surface height for the results of MSDS, MFS, and experiment separately. It could be obtained that the simulated surface of MSDS is closer to the machined surface; SD′2 is less than SD′1 and SD′3, which indicates that the surface height of MFS is more consistent; this can also be seen in Fig. 19c. Fig. 22 b shows the characteristic of adjacent peak intervals for simulated surfaces and machined surface, where MV″1, MV″2, and MV″3 are the mean values of peak interval for the results of MSDS, MFS, and experiment separately and SD″1, SD″2, and SD″3 are the standard deviation of peak interval for the results of MSDS, MFS, and experiment separately. It can be seen that the disparity between MV″1, MV″2, and MV″3 is less distinctive, while SD ″2 is lower than SD″1 and SD″3, which indicates that the peak interval for the simulated surface of MFS is more consistent. The peak interval of the machined surface depends on the size of the discharge crater. Therefore, we compared the average peak interval with the average diameter of the craters and obtained that the average diameter of the discharge crater is about 1.6 times the average peak interval in the calculation. The average prediction errors of MSDS and MFS are listed in Table 6, and it can be seen that the prediction errors of MV′ (Ra) of the two models are almost the same, while the other prediction errors are quite different; this is because the average surface height of the two models is almost the same, although the surface morphology characteristics are different. As shown in Table 6, the overall average prediction error is defined as: where e MV' and e SD' are average prediction error for mean value (MV′) and standard deviation (SD′) of surface height respectively, e MV" and e SD" are average prediction error for mean value and standard deviation of peak intervals respectively, and e o is the average of overall prediction error. As shown in Table 6, it can be seen that MSDS has higher prediction accuracy.

Conclusions
In this study, the stochastic distribution characteristics of the discharge crater size were investigated, and a geometric prediction model of surface morphology for micro-EDM was developed. The main conclusions are drawn as follows: 1. When the machining parameters is unchanged, the crater size distribution can be characterized by normal probability density function, and the main reason for the different crater sizes is the different discharge durations. 2. The influence of pulse width on the average crater size decreases with the decrease of peak current. Under small peak current and large gap width, the maintenance capability of the discharge channel is weak, which leads to the discharge process is easily interrupted. The standard deviation of the crater size increases with the increase of the pulse width and decrease of the maintenance capability of the discharge channel. 3. The surface height and adjacent peak interval of the machined surface are normally distributed, and these characteristics could be simulated by MSDS (the model based on stochastic distribution of crater size), while MFS (the model based on fixed crater size) could not simulate the distribution characteristics of surface height and peak interval. The peak interval is influenced by crater size, and the average diameter of crater is about 1.6 times the average peak interval under the experimental conditions of this study.  . The prediction errors for the surface roughness of MSDS and MFS are 7.50% and 8.51%, respectively, while the overall average errors for the two models are 7.65% and 26.62%, respectively, which indicates that MSDS has higher simulation accuracy.