Three-dimensional surface residual stress prediction method based on end milling thickness differential calculation and moving heat source simulation

By combining differential calculation of end milling thickness with simulation of moving heat source, the problem of predicting three-dimensional residual stress during end milling was solved, and accurate prediction of three-dimensional residual stress was achieved. In particular, in the processing of thin-walled quasi-shaped parts and titanium alloy materials, the prediction accuracy and temperature prediction accuracy were improved.

CN121503106APending Publication Date: 2026-02-10SHANGHAI JIAOTONG UNIV
View PDF 0 Cites 0 Cited by

Patent Information

Application Number
CN202411089115.2
Authority / Receiving Office
CN · China
Patent Type
Applications(China)
Current Assignee / Owner
Filing Date
2024-08-09
Publication Date
2026-02-10

AI Technical Summary

Technical Problem

Existing technologies struggle to accurately predict three-dimensional residual stress during end milling, especially in the machining of thin-walled quasi-shaped parts and the cutting of titanium alloy materials. The complexity of the cutting method and thermal stress distribution makes analytical calculation of residual stress difficult.

Method used

A method based on differential calculation of end milling thickness and simulation of moving heat source is adopted. By calculating the cutting width and mechanical stress distribution by differentiation and combining it with thermal stress simulation, the three-dimensional residual stress is coupled and predicted.

Benefits of technology

It improves the prediction accuracy of residual stress during end milling, and can accurately predict the residual stress distribution of the machined surface in three-dimensional space. In particular, it takes into account the influence of the end milling side on the bottom residual stress, and improves the accuracy of temperature prediction.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN121503106A_ABST
    Figure CN121503106A_ABST
Patent Text Reader

Abstract

The invention relates to a three-dimensional surface residual stress prediction method based on end milling thickness differential calculation and moving heat source simulation. The method comprises the following steps: calculating a tool path in an end milling process; differentiating the machined surface generated in one rotation period, and calculating the cutting width of any differential unit; the shear force and the ploughing force of any infinitesimal are calculated through the constitutive relation and coordinate conversion; calculating the mechanical stress distribution of any infinitesimal in a rotation period; setting a moving heat source to simulate and calculate the thermal stress distribution of any infinitesimal; mechanical stress and thermal stress are coupled, the residual stress value of each infinitesimal element is calculated, and stress elastic-plastic loading and unloading calculation in a rotation period is carried out; inputting the residual stress distribution into a side wall influence function, and calculating to obtain final residual stress distribution of any infinitesimal in a rotation period; and cyclic calculation is completed by using a superposition method, and the residual stress distribution of the fully-machined surface is obtained. Compared with the prior art, the method has the advantages of improving the precision of residual stress calculation and the like.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention relates to the field of machining, and in particular to a method for predicting three-dimensional surface residual stress based on end milling thickness differential calculation and moving heat source simulation. Background Technology

[0002] Predicting residual stress in parts is a crucial evaluation indicator for designing machining processes and parameters, as the distribution of residual stress significantly impacts part machining quality and subsequent service performance. Liang et al., in their paper "Residual Stress Modeling in Orthogonal Machining. (CIRP Annals 2007. Volume 56. Issue 1.65-68.)", proposed a model for predicting residual stress in orthogonal cutting. This model uses process conditions as input to predict the distribution of surface and subsurface residual stress generated during machining. Xiaoliang Liang et al., in their paper "Prediction of residual stress with multi-physics model for orthogonal cutting Ti-6Al-4V under various tool wear morphologies. (Journal of Materials Processing Tech. 288(2021)116908)", considered tool geometry changes caused by wear and predicted the residual stress distribution of Ti-6Al-4V during orthogonal turning. In their paper "Theoretical prediction of machining-induced residual stresses in three-dimensional oblique milling processes. (International Journal of Mechanical Sciences 133(2017)426–437.)", Min Wan et al. achieved the prediction of three-dimensional oblique milling residual stress by introducing a helix angle, thus enabling the prediction of the depth direction of residual stress during side milling. Although the cutting model can achieve three-dimensional calculation, the output result is still a one-dimensional residual stress distribution, which cannot solve the problem of predicting three-dimensional residual stress.

[0003] End milling is an important finishing milling process, especially for the final pass of thin-walled, near-circular parts. The residual stress generated during end milling is a crucial indicator for subsequent workpiece processing and quality assessment. The cutting method of the end mill during end milling differs significantly from orthogonal turning or side milling. End milling involves dual cutting of the bottom and sides, and because the tool path exhibits rotary feed, the cutting thickness within the same rotational cycle is not constant but continuously varying. Therefore, the distribution of residual stress within the same rotational cycle is also uneven, making analytical calculation of residual stress quite difficult. For end milling materials with strong tool affinity and significant thermal adhesion effects, the intense and uneven temperature distribution during machining, coupled with the influence of chips, makes accurate calculation of thermal stress analysis models difficult. Therefore, other auxiliary methods are needed to apply thermal stress. How to analytically predict residual stress in three-dimensional space under the complex cutting conditions of end milling is of considerable significance for evaluating the machining process of thin-walled end milling, especially for titanium alloys. Summary of the Invention

[0004] The purpose of this invention is to provide a method for predicting three-dimensional surface residual stress based on end milling thickness differential calculation and moving heat source simulation.

[0005] The objective of this invention can be achieved through the following technical solutions:

[0006] A method for predicting residual stress on a three-dimensional surface based on differential calculation of end milling thickness and simulation of a moving heat source includes the following steps:

[0007] Step S1: Obtain the end milling process parameters and the geometry of the end mill bottom edge, and calculate the tool path during the end milling process;

[0008] Step S2: Differentiate the machined surface produced in one rotation cycle and calculate the cutting width of any differential unit;

[0009] Step S3: Based on the cutting width calculated in step S2, calculate the shear force and plowing force of any infinitesimal element using constitutive relations and coordinate transformation;

[0010] Step S4: Calculate the mechanical stress distribution of any micro-element within one rotation cycle based on shear force and plowing force;

[0011] Step S5: Set up a moving heat source to simulate and calculate the thermal stress distribution of any micro-element;

[0012] Step S6: Couple mechanical stress with thermal stress, calculate the residual stress value of each micro-element, and perform stress elastoplastic loading and unloading calculations within one rotation cycle;

[0013] Step S7: Input the residual stress distribution obtained in step S6 into the sidewall influence function to calculate the final residual stress distribution of any micro-element within one rotation cycle;

[0014] Step S8: Use the superposition method to complete the iterative calculation and obtain the residual stress distribution on the fully machined surface.

[0015] Step S2 specifically includes the following steps:

[0016] Step S21: Discretize one rotation period into N s There are n time increments, where i is the number of time increments for any infinitesimal element, i∈[1,N]. s Discretize the depth of any infinitesimal element from the surface H into N. H There are n time increments, where j is the depth increment of any infinitesimal element, j∈[1,N]. H ], k represents the k-th tooth of a bottom plane, calculate its immersion angle. :

[0017] φ p =2π / N t

[0018]

[0019] Where n is the rotational speed, t i For the current moment accumulated from i time increments, N t This represents the total number of teeth on the cutting tool.

[0020] Step S22: Complete the differentiation of the bottom machined surface, and calculate the static cutting width w(i,k) of any infinitesimal element P(i,j,k) of each tooth in one rotation cycle:

[0021]

[0022] Among them, f t is the feed per tooth, m is a dimensionless integer, and r is the radius of each bottom tooth, i.e., the end mill radius.

[0023] Specifically, step S3 involves: based on the constitutive relation of the material, using the constitutive model to solve the shear flow stress of the material deformation, calculating the shear force for any point P(i,j,k) within a rotational cycle path, and calculating the plowing force of each micro-element using the slip line field method.

[0024] The shear force F s The calculation method is as follows:

[0025]

[0026] Where, τ s L is the shear yield strength.AB Let φ be the shearing zone length, φ be the shearing angle, w(i,k) be the cutting width, and r be the shearing angle. e For end milling thickness, shear force F s The distribution is uneven on the shear plane.

[0027] The method for calculating the plowing force is as follows:

[0028]

[0029] Where φ is the shear angle during end milling, i.e., the angle at which the material actually undergoes plastic deformation; ζ is the tool rake angle, i.e., the angle between the tool edge and the cutting plane; θ is the angle between the chip thickness and the cutting speed direction; η is the angle between the chip thickness and the normal to the worktable surface; and r... eq τ is the equivalent radius of the transition edge of the plowing bottom edge during end milling. s Let w(i,k) be the shear yield strength, and w(i,k) be the cutting width.

[0030] Specifically, step S4 involves: calculating the distribution of shearing force and plowing force on the tool trajectory based on step S3, transforming the mechanical stress generated by shearing in the local coordinate system to the global coordinate system, and solving for the stress distribution at any point P(i,j,k) on the machined surface within one rotation cycle.

[0031] Step S5 specifically includes the following steps:

[0032] Step S51: Analyze and calculate the heat source during the milling process, and calculate the temperature rise acting on the machined surface;

[0033] Step S52: Input the temperature rise history of each micro-element into the simulation software, set the movement trajectory of the heat source and milling and the material properties, use the heat transfer analysis step, and set the mesh type to obtain the temperature distribution of the entire single-pass milling process through simulation.

[0034] Step S53: Import the temperature distribution over time as the initial state into the next simulation step, set the static general analysis step, constrain the bottom surface of the workpiece to be fixed, and obtain the thermal stress distribution of the workpiece due to temperature during the milling process.

[0035] Step S54: Export the thermal stress loading dataset for any infinitesimal element P(i,j,k) from the running results within one rotation cycle.

[0036] Step S6 specifically involves: sequentially coupling the mechanical stress set and the thermal stress set, and calculating the residual stress value of each micro-element P(i,j,k); for the stress within one rotation cycle, the process is divided into loading and unloading. During loading, the state of P(i,j,k) determines whether it is elastic loading or plastic loading. When it is elastic loading, the stress increment is calculated according to P(i,j,k)-P(i-1,j,k). When it is plastic loading, the stress increment is recalculated according to the yield criterion; during unloading, the Ducker postulate determines whether to enter the unloading process, and stress unloading is completed when the boundary conditions are met.

[0037] Step S7 specifically involves: forming a two-dimensional distribution of residual stress P(i) and a depth distribution P(j) on the end-milled surface, thereby forming a three-dimensional residual stress distribution on the machined surface. Considering the influence of sidewalls and other factors on the residual stress distribution, a sidewall function H(i,j) is set and substituted into the calculation to obtain the latest three-dimensional distribution of residual stress.

[0038] Step S8 specifically involves the following steps: After the residual stress of any P(i,j,k) within a rotation cycle is calculated, in the subsequent milling process, a cyclic calculation is performed along the width and length directions. The calculated residual stress data within the rotation cycle is combined into the first cycle. Subsequent cycles are obtained by repeatedly superimposing the stress from the first cycle. If the results of the cyclic calculation are superimposed, the latter is taken. The superposition method is used to obtain the final three-dimensional distribution cloud map of the residual stress on the entire surface.

[0039] Compared with the prior art, the present invention has the following beneficial effects:

[0040] (1) End milling process is different from traditional cutting calculation. In the in-plane cutting process of micro-element, the thickness changes continuously. The method of this invention expands the dimension of residual stress calculation.

[0041] (2) End milling involves double cutting of the bottom and side surfaces during the bottom milling process. This invention considers the influence of the end milling side surface on the residual stress of the machined bottom surface, thus improving the accuracy of residual stress prediction.

[0042] (3) This invention introduces the simulation of moving heat source in the field of welding, and innovatively combines analytical mechanical stress calculation with moving heat source simulation fitting of thermal stress, optimizes the residual stress prediction logic, and improves the accuracy of temperature prediction. Attached Figure Description

[0043] Figure 1 This is a flowchart of the method of the present invention;

[0044] Figure 2 This is a tool trajectory diagram for end milling of the entire surface in one embodiment;

[0045] Figure 3 This is a schematic diagram illustrating the division of the machined surface into micro-elements within one rotation cycle of end milling in one embodiment;

[0046] Figure 4 This is a schematic diagram of shear force during end milling in one embodiment;

[0047] Figure 5 This is a schematic diagram of the plowing force during end milling in one embodiment. Detailed Implementation

[0048] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments. These embodiments are based on the technical solution of the present invention and provide detailed implementation methods and specific operating procedures. However, the scope of protection of the present invention is not limited to the following embodiments.

[0049] This invention relates to a method for predicting three-dimensional surface residual stress based on differential calculation of end-milling thickness and simulation fitting of a moving heat source. The method first obtains the end-milling process parameters and end-milling cutter geometry, and calculates the tool trajectory across the entire surface. Then, it differentiates the cutting thickness within one rotation cycle, and for any differential unit, calculates shear force and plowing force to complete the loading of mechanical stress during end-milling. Finally, it combines moving heat source simulation to realize thermal stress loading, and uses elastoplastic unloading to achieve three-dimensional residual stress prediction based on dynamic thickness changes within one rotation cycle. Compared with existing technologies, this invention, on the one hand, utilizes an innovative differential calculation method to calculate three-dimensional residual stress during end-milling, predicting both the in-plane distribution of residual stress on the machined surface and the distribution of residual stress in the depth direction; on the other hand, it improves the accuracy of thermal stress prediction through simulation calculation with a moving heat source, and further improves the accuracy of residual stress calculation through sequential coupling.

[0050] Specifically, such as Figure 1 As shown, the method includes the following steps:

[0051] Step S1: Obtain end milling process parameters (end milling depth is r) e Feed rate V, spindle speed n, depth of feed z) and end mill bottom edge geometry (tool radius R, total number of teeth N) t Using the center of the end mill's bottom surface as the origin of the coordinate system, the trajectory of any infinitesimal element in a single pass is obtained. Using the tool spindle as the coordinate center, a global coordinate system is established, yielding the global trajectory of the tool movement, as shown below. Figure 2 .

[0052] Step S2: Differentiate the machined surface generated in one rotation cycle. Since the cutting thickness is consistent and the stress is uniformly distributed in the radial direction, the change in radial coordinate is ignored. Represent each infinitesimal element P on the path as P(i,j,k), and calculate the cutting width of any infinitesimal element.

[0053] Step S2 specifically includes the following steps:

[0054] Step S21: Discretize one rotation period into N s There are n time increments, where i is the number of time increments for any infinitesimal element, i∈[1,N]. s ], which also represents the change in position of the bottom cutting edge of the k-th tooth within one rotation cycle; the depth of any infinitesimal element at surface H is discretized into N. H There are n time increments, where j is the depth increment of any infinitesimal element, j∈[1,N]. H ], k represents the k-th tooth on a bottom plane (generally a four-tooth end mill is selected), and the immersion angle of the k-th tooth in one rotation cycle. for:

[0055] φ p =2π / N t

[0056]

[0057] Where n is the rotational speed, t i For the current moment accumulated from i time increments, N t This represents the total number of teeth on the cutting tool.

[0058] Step S22: Complete the differentiation of the bottom machined surface. The differentiation state of the bottom machined surface is as follows: Figure 3 Calculate the static cutting width of any infinitesimal element P(i,j,k) for each tooth within one rotation cycle. Substituting equation (1) into the calculation, the static cutting width of any infinitesimal element for each tooth within one rotation cycle is:

[0059]

[0060] Among them, f t is the feed per tooth, m is a dimensionless integer, and r is the radius of each bottom tooth, i.e., the end mill radius.

[0061] Step S3: Based on the cutting width calculated in step S2, calculate the shear force and plowing force of any infinitesimal element using constitutive relations and coordinate transformation.

[0062] Specifically, based on the constitutive relation of material TC4, the shear stress and plowing force are determined using the JC constitutive model. For any point P(i,j,k) within a rotational cycle path, the stress mainly consists of shear force and plowing force, with the shear force F... s The milling direction and the normal direction can be represented as follows, and their distribution on the shear plane is not uniform:

[0063]

[0064] Where, τs L represents the shear yield strength of TC4. AB r is the length of the shearing region. e For end milling thickness, φ is the shear angle, such as Figure 4 As shown.

[0065] The plowing force at any point P(i,j,k) is calculated using Waldorf's slip-line field method:

[0066]

[0067] Among them, such as Figure 5 As shown, φ is the shear angle during end milling, i.e., the angle at which the material actually undergoes plastic deformation; ζ is the tool rake angle, i.e., the angle between the tool edge and the cutting plane; θ is the angle between the chip thickness and the cutting speed direction; η is the angle between the chip thickness and the normal to the worktable surface; and r... eq This is the equivalent radius of the transition edge of the plowing bottom edge during end milling.

[0068] Step S4: Calculate the mechanical stress distribution of any micro-element within one rotation cycle based on shear force and plowing force.

[0069] Specifically, based on step S3, the distribution of shearing force and plowing force on the tool trajectory is calculated. The mechanical stress generated by shearing in the local coordinate system is transformed into the global coordinate system, and the stress distribution at any point P(i,j,k) on the machined surface within one rotation cycle is obtained.

[0070] Let x′ be the global coordinates after coordinate transformation during the shearing calculation. The stress distribution at any point P on the machined surface can be obtained by solving equations (3) and (4).

[0071] For time i, when the tool is located at point P(i,j,k), solve for the stress distribution at time i+m on the trajectory of the same rotation period, where the x-direction distance from the position at time i+m is:

[0072]

[0073] Where r is the radius of the bottom cutting edge of the milling cutter.

[0074] Step S5: Set up a moving heat source to simulate and calculate the thermal stress distribution of any micro-element.

[0075] Step S5 specifically includes the following steps:

[0076] Step S51: Analyze and calculate the heat source during the milling process, and calculate the temperature rise acting on the machined surface.

[0077] Step S52: The DFLUX subroutine is used to set up the moving heat source. The temperature rise history of each micro-element is input into the simulation software. The movement trajectory of the heat source and the milling operation, as well as the TC4 material properties, are set. The Heat Transfer analysis step is used, and the mesh type is set to DC3D8. The temperature distribution of the entire single-pass milling process is obtained through simulation. In this embodiment, the workpiece size is set to 5×40×20mm. The Fortran file is written as the input file based on the coordinates, the Gaussian heat source equation, and the trajectory changes over time.

[0078] Step S53: Import the temperature distribution over time as the initial state (Predefined Field) into the next simulation step, set the static general analysis step (Static, General), constrain the bottom surface of the workpiece to be fixed and use the mesh type C3D8R to obtain the thermal stress distribution of the workpiece due to temperature during the milling process.

[0079] Step S54: Export the thermal stress loading dataset for any infinitesimal element P(i,j,k) from the running results within one rotation cycle.

[0080] Step S6: Couple mechanical stress with thermal stress, calculate the residual stress value of each micro-element, and perform stress elastoplastic loading and unloading calculations within one rotation cycle.

[0081] Specifically, the mechanical stress set and the thermal stress set are sequentially coupled to calculate the residual stress value of each micro-element P(i,j,k). The stress within one rotation cycle is divided into loading and unloading processes. During loading, the state of P(i,j,k) determines whether it is elastic or plastic loading. When it is elastic loading, the stress increment is calculated according to P(i,j,k) - P(i-1,j,k). When it is plastic loading, the stress increment is recalculated according to the yield criterion. During unloading, the Ducker postulate determines whether to enter the unloading process; stress unloading is completed when the boundary conditions are met.

[0082] Step S7: Input the residual stress distribution obtained in step S6 into the sidewall influence function to calculate the final residual stress distribution of any micro-element within one rotation cycle.

[0083] Specifically, the residual stress distribution of the end-milled surface is formed in two dimensions, P(i) and depth, P(j), thus forming the residual stress distribution of the three-dimensional machined surface. Considering the influence of the sidewall and other factors on the residual stress distribution, the sidewall function H(i,j) is set and substituted into the calculation to obtain the latest three-dimensional distribution of residual stress. In particular, attention should be paid to the transformation between different coordinate systems of P(i,j,k).

[0084] Step S8: Use the superposition method to complete the iterative calculation and obtain the residual stress distribution on the fully machined surface.

[0085] Specifically, once the residual stress of any P(i,j,k) within a rotation cycle is calculated, in subsequent milling processes, such as... Figure 2 The prop trajectory shown has consistent material and process parameters in each rotation cycle, so it is possible to perform cyclic calculations along the width and length directions. The residual stress data in the calculated rotation cycles are combined into the first cycle, and subsequent cycles are obtained by repeatedly superimposing the stress in the first cycle. If the results of the cyclic calculations are superimposed, the latter is taken. The superposition method is used to obtain the final three-dimensional distribution cloud map of the residual stress on the entire surface.

[0086] The preferred embodiments of the present invention have been described in detail above. It should be understood that those skilled in the art can make numerous modifications and variations based on the concept of the present invention without creative effort. Therefore, all technical solutions that can be obtained by those skilled in the art based on the concept of the present invention through logical analysis, reasoning, or limited experimentation on the basis of existing technology should be within the scope of protection defined by the claims.

Claims

1. A method for predicting residual stress on a three-dimensional surface based on differential calculation of end-milling thickness and simulation of a moving heat source, characterized in that, Includes the following steps: Step S1: Obtain the end milling process parameters and the geometry of the end mill bottom edge, and calculate the tool path during the end milling process; Step S2: Differentiate the machined surface produced in one rotation cycle and calculate the cutting width of any differential unit; Step S3: Based on the cutting width calculated in step S2, calculate the shear force and plowing force of any infinitesimal element using constitutive relations and coordinate transformation; Step S4: Calculate the mechanical stress distribution of any micro-element within one rotation cycle based on shear force and plowing force; Step S5: Set up a moving heat source to simulate and calculate the thermal stress distribution of any micro-element; Step S6: Couple mechanical stress with thermal stress, calculate the residual stress value of each micro-element, and perform stress elastoplastic loading and unloading calculations within one rotation cycle; Step S7: Input the residual stress distribution obtained in step S6 into the sidewall influence function to calculate the final residual stress distribution of any micro-element within one rotation cycle; Step S8: Use the superposition method to complete the iterative calculation and obtain the residual stress distribution on the fully machined surface.

2. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 1, characterized in that, Step S2 specifically includes the following steps: Step S21: Discretize one rotation period into Ns time increments, where i is the time increment number of any infinitesimal element, i∈[1,N]. s Discretize the depth of any infinitesimal element from the surface H into N. H Let j be the depth increment of any infinitesimal element, j∈[1, NH], and k represent the k-th tooth of a bottom plane. Calculate its immersion angle. φp=2π / Nt Where n is the rotational speed, t i Let Nt be the current moment accumulated over i time increments, and Nt be the total number of teeth on the tool. Step S22: Differentiate the bottom machined surface and calculate the static cutting width w(i,k) of any infinitesimal element P(i,j,k) for each tooth in one rotation cycle: Among them, f t is the feed per tooth, m is a dimensionless integer, and r is the radius of each bottom tooth, i.e., the end mill radius.

3. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 1, characterized in that, Step S3 specifically involves: based on the constitutive relation of the material, using the constitutive model to solve the shear flow stress of the material deformation, calculating the shear force for any point P(i, j, k) within a rotational cycle path, and calculating the plowing force of each micro-element according to the slip line field method.

4. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 3, characterized in that, The shear force F s The calculation method is as follows: Where, τ s L is the shear yield strength. AB Let φ be the shearing zone length, φ be the shearing angle, w(i,k) be the cutting width, and r be the shearing angle. e For end milling thickness, shear force F s The distribution is uneven on the shear plane.

5. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 3, characterized in that, The method for calculating the plowing force is as follows: Where φ is the shear angle during end milling, i.e., the angle at which the material actually undergoes plastic deformation; ζ is the tool rake angle, i.e., the angle between the tool edge and the cutting plane; θ is the angle between the chip thickness and the cutting speed direction; η is the angle between the chip thickness and the normal to the worktable surface; and r... eq τ is the equivalent radius of the transition edge of the plowing bottom edge during end milling. s Let w(i, k) be the shear yield strength, and w(i, k) be the cutting width.

6. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 1, characterized in that, Specifically, step S4 involves: calculating the distribution of shearing force and plowing force on the tool trajectory based on step S3, transforming the mechanical stress generated by shearing in the local coordinate system to the global coordinate system, and solving for the stress distribution at any point P(i, j, k) on the machined surface within one rotation cycle.

7. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 1, characterized in that, Step S5 specifically includes the following steps: Step S51: Analyze and calculate the heat source during the milling process, and calculate the temperature rise acting on the machined surface; Step S52: Input the temperature rise history of each micro-element into the simulation software, set the movement trajectory of the heat source and milling and the material properties, use the heat transfer analysis step, and set the mesh type to obtain the temperature distribution of the entire single-pass milling process through simulation. Step S53: Import the temperature distribution over time as the initial state into the next simulation step, set the static general analysis step, constrain the bottom surface of the workpiece to be fixed, and obtain the thermal stress distribution of the workpiece due to temperature during the milling process. Step S54: Export the thermal stress loading dataset for any infinitesimal element P(i,j,k) from the running results within one rotation cycle.

8. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 1, characterized in that, Step S6 specifically involves: sequentially coupling the mechanical stress set and the thermal stress set, and calculating the residual stress value of each micro-element P(i,j,k); for the stress within one rotation cycle, it is divided into loading and unloading processes. During the loading process, the state of P(i,j,k) is used to determine whether it is elastic loading or plastic loading. When it is elastic loading, the stress increment is calculated according to the stress increment of P(i,j,k)-P(i-1,j,k). When it is plastic loading, the stress increment is recalculated according to the yield criterion. During the unloading process, the Ducker postulate is used to determine whether to enter the unloading process. Stress unloading is completed when the boundary conditions are met.

9. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 1, characterized in that, Specifically, step S7 involves: forming a two-dimensional distribution of residual stress P(i) and a depth distribution P(j) on the end-milled surface, thereby forming a three-dimensional residual stress distribution on the machined surface. Considering the influence of sidewalls and other factors on the residual stress distribution, a sidewall function H(i,j) is set and substituted into the calculation to obtain the latest three-dimensional distribution of residual stress.

10. The method for predicting three-dimensional surface residual stress based on end-milling thickness differential calculation and moving heat source simulation according to claim 1, characterized in that, Step S8 specifically involves the following steps: After the residual stress of any P(i, j, k) within a rotation cycle is calculated, in the subsequent milling process, a cyclic calculation is performed along the width and length directions. The calculated residual stress data within the rotation cycle is combined into the first cycle. Subsequent cycles are obtained by repeatedly superimposing the stress of the first cycle. If the results of the cyclic calculation are superimposed, the latter is taken. The superposition method is used to obtain the final three-dimensional distribution cloud map of the residual stress on the entire surface.