Prediction Method for Milling Stability of Thin-Walled Parts with Multiple Mass Points and Multiple Modes Considering Modal Coupling Effect

Abstract

The present invention provides a method for predicting the milling stability of thin-walled parts with multiple mass points and multiple modes considering the modal coupling effect, which includes the following steps: meshing the milling cutter and the thin-walled part to establish a multi-mass-point and multi-mode milling dynamics model; considering the modal coupling effect of each mass point and introducing a multi-order modal vibration mode matrix to correct the milling dynamics model; calculating the change of the instantaneous cutting-in angle and cutting-out angle according to the action of the milling cutter helix angle on the cutting surface of the cutter tooth and the thin-walled part micro-element; establishing a state space equation with periodic time delay according to the multi-mass-point and multi-mode milling dynamics model of the thin-walled part considering the modal coupling effect; constructing a reconstructed Newton-Cotes fitting method to calculate the transfer matrix according to the periodic time delay term of the discretized state space equation; outputting a three-dimensional stability lobe diagram according to the transfer matrix reduction and dimension reduction method and improving the output efficiency.

Description

Technical Field

[0001] The invention belongs to the technical field of thin-walled part milling, and in particular relates to a method for predicting the stability of multi-particle multi-mode milling of thin-walled parts taking into account modal coupling effects. Background Art

[0002] Thin-walled parts are characterized by light weight and large aspect ratio, and are widely used in aerospace, medical equipment, ships and other fields. However, thin-walled parts have weak rigidity and low strength, and are prone to chatter during milling. Chatter will affect the surface quality, reduce processing efficiency, and even cause equipment failure. Therefore, it is necessary to accurately model and predict the stability state of the milling process of thin-walled parts.

[0003] The dynamic system involved in milling thin-walled parts exhibits position dependence and aggregated modal characteristics. During the milling process, the low-order modes of the tool and the thin-walled part can be excited by the cutting force at the same time, resulting in a modal coupling effect. In addition, the dynamic response characteristics change with the cutting position. Therefore, the traditional single-point and single-mode dynamic models introduce significant approximation errors. To this end, the tool and the workpiece can be spatially discretized to establish a high-dimensional lumped parameter dynamic model, that is, a multi-point multi-modal lumped parameter model for milling dynamics and stability analysis of aluminum alloy thin-walled parts with integrated modal coupling effects, to describe the machining vibration of the low-rigidity milling system with high precision, and provide a reliable framework for predicting and improving the stability of thin-walled parts processing. Summary of the invention

[0004] In view of the deficiencies in the prior art, the present invention provides a method for predicting the stability of multi-particle multi-mode milling of thin-walled parts taking into account the modal coupling effect.

[0005] In the first aspect, the present invention provides a multi-particle multi-modal milling dynamics model of thin-walled parts considering the modal coupling effect, comprising: constructing an initial milling dynamics model including the excitation of the cutting force in the feed direction and the radial cutting force; discrete tools along the axial direction are N t Mass points, discrete thin-walled parts N along the feed direction wpmass points, obtain the dynamic parameters of each mass point, and establish a first milling dynamic model of the interface between the milling cutter and the thin-walled part at different processing positions in the feed direction according to the initial milling dynamic model to simulate the change trend of the stability state of the milling cutter cutting process after discreteness; convert the vibration displacement in the first milling dynamic model from Cartesian space to modal space, consider the coupling effect between the low-order mode and the high-order mode of each mass point, introduce a multi-order modal vibration matrix to correct the milling dynamic model, and obtain a second milling dynamic model considering the modal coupling effect; further, the construction of the initial milling dynamic model including the excitation of the cutting force in the feed direction and the radial cutting force includes: constructing the initial milling dynamic model expression including the excitation of the cutting force in the feed direction and the radial cutting force:

[0006] ,

[0007] Where M is the mass matrix in Cartesian space; C is the damping matrix in Cartesian space; K is the stiffness matrix in Cartesian space; is the vibration acceleration vector at time t; is the vibration velocity vector at time t; is the vibration displacement vector at time t; is the cutting force vector of the milling cutter at time t; is the axial cutting depth; is the number of milling cutter teeth; is the rotation matrix in Cartesian space; A window function to determine whether the cutter tooth is in the cutting state; is the cutting force coefficient matrix; is the rotation matrix of the jth tooth direction at the i-th node of the milling cutter at time t; T is the hysteresis period.

[0008] Furthermore, the axially discrete tool is N t Mass points, discrete thin-walled parts N along the feed direction wp mass points, obtain the dynamic parameters of each mass point, and establish the first milling dynamic model of the interface between the milling cutter and the thin-walled workpiece at different processing positions in the feed direction according to the initial milling dynamic model, so as to simulate the change trend of the stability state of the milling cutter cutting process after discrete cutting, including:

[0009] Using the finite element method, the tool is discretely divided into N along the axial direction. t The discrete thin-walled parts along the feed direction are N wp particles, forming N t ×N wp The dynamic parameters of each particle are obtained respectively, and the multi-particle multi-mode-milling dynamic equations are constructed based on the initial milling dynamics model:

[0010] ,

[0011] Where M tool is the tool quality matrix after discretization; C tool is the discretized tool damping matrix; K tool is the tool stiffness matrix after discretization; M wp is the mass matrix of the thin-walled part after discretization; C wp is the damping matrix of the discretized thin-walled component; K wp is the stiffness matrix of the thin-walled component after discretization; is the vibration acceleration vector of the milling cutter at time t after the discretization; is the vibration velocity vector of the milling cutter at time t after the discretization; is the vibration displacement vector of the milling cutter at time t after discretization; is the vibration acceleration vector of the thin-walled part at time t after discretization; is the vibration velocity vector of the thin-walled part at time t after discretization; is the vibration displacement vector of the thin-walled part at time t after discretization; is the cutting force vector of the milling cutter at time t after discretization. The mass matrix, damping matrix and stiffness matrix in the discretized milling dynamics system are related to the dynamic parameters of each particle. The milling cutter and thin-walled workpiece milling dynamics system are coupled into an overall dynamics system about the particle. The first milling dynamics model of the interface between the milling cutter and the thin-walled workpiece obtained by coupling is:

[0012] .

[0013] Furthermore, the vibration displacement in the first milling dynamics model at different processing positions in the feed direction is converted from Cartesian space to modal space, the coupling effect between the low-order mode and the high-order mode of each particle is considered, and the multi-order modal vibration matrix is ​​introduced to correct the milling dynamics model, so as to obtain the second milling dynamics model considering the modal coupling effect, including:

[0014] According to the modal superposition theory, the transformation relationship of the vibration displacement vector from Cartesian space to modal space is calculated as follows:

[0015] ,

[0016] In the formula, represents the vibration displacement expression in the modal space, and U represents the modal vibration matrix. The modal vibration matrix of the milling cutter and thin-walled parts is shown as follows:

[0017] ,

[0018] Where U tool Represents the milling cutter mode shape; U wp Represents the modal shape matrix of thin-walled parts; and Represent the j-order modal vibration vector in the x-direction and y-direction respectively. The multi-order modal vibration matrix of the low-order mode and high-order mode coupling effect of the milling cutter and the workpiece is shown as follows:

[0019] ,

[0020] ,

[0021] ,

[0022] ,

[0023] In the formula, represents the i-th vibration mode in the x direction; Represents the i-th order vibration mode in the y direction. Based on the conversion relationship between Cartesian space and modal space and the principle of modal coupling, the multi-particle multi-modal milling dynamics expression can be rewritten as the second milling dynamics model as follows:

[0024] ,

[0025] In the formula, represents the modal mass matrix; represents the modal damping matrix; represents the modal stiffness matrix; represents the vibration acceleration at time t in the modal space; represents the vibration velocity at time t in the modal space; represents the vibration displacement at time t in the modal space; a p represents the axial cutting depth; U represents the overall modal vibration matrix; U T represents the transpose of the overall modal vibration matrix; U wp Represents the modal shape matrix of thin-walled parts; represents the transpose of the modal vibration matrix of the thin-walled part; U tool represents the milling cutter mode shape; represents the transpose of the milling cutter mode shape matrix; Represents the milling force coefficient matrix of the jth infinitesimal element of the i-th tooth at time t.

[0026] In the second aspect, the present invention provides a method for predicting the stability of thin-walled workpiece milling taking into account the modal coupling effect, including: calculating the changes in the instantaneous cutting-in angle and cutting-out angle when the helix angle of the milling cutter acts on the cutter teeth and the micro-element cutting surface of the thin-walled workpiece; establishing a state-space equation with periodic time lag based on the second milling dynamics model; constructing a reconstructed Newton-Cotes fitting formula based on the periodic time lag term of the discretized state-space equation, and calculating the transfer matrix; reducing the dimension based on the transfer matrix, configuring a three-dimensional stability lobe diagram drawing method, and predicting the milling stability.

[0027] Further, the calculation of the instantaneous cutting angle and cutting-out angle changes when the milling cutter helix angle acts on the cutter tooth and the thin-walled workpiece micro-element cutting surface includes: constructing an instantaneous cutting angle model according to the helix angle on the cutter tooth and the thin-walled workpiece micro-element cutting surface as follows:

[0028] ,

[0029] In the formula, represents the instantaneous cutting angle at time t; Ω represents the spindle speed; i represents the i-th tooth; j represents the j-th discrete element along the axial direction of the milling cutter; m represents the axial discrete number of the milling cutter; N f Indicates the number of milling cutter teeth; a p It represents the axial cutting depth; β represents the helix angle of the milling cutter; and D represents the diameter of the milling cutter.

[0030] Furthermore, the state space equation with periodic time lag is established according to the second milling dynamics model, including:

[0031] ,

[0032] The second milling dynamics model considering the modal coupling effect can be transformed into a state space expression as follows:

[0033] ,

[0034] In the formula, Represents state variables; represents the vibration velocity at time t in the modal space; represents the vibration displacement at time t in the modal space; represents the first-order derivative of the state variable; A 0 represents the state matrix; B(t) represents the observation matrix; a p Indicates the axial cutting depth; represents the inverse of the modal mass; represents the modal stiffness; represents modal damping; I represents the unit matrix; U represents the overall modal vibration matrix; represents the transpose of the overall mode shape matrix; Represents the milling force coefficient matrix of the jth infinitesimal element of the i-th tooth at time t.

[0035] Further, the reconstructing the Newton-Cotes fitting formula and calculating the transfer matrix according to the periodic time lag term of the discretized state-space equation include: according to the numerical solution method of the second kind of Volterra equation, the state-space equation can be solved as the following expression:

[0036] ,

[0037] In the formula, Represents state variables; represents the i-th discretized state vector; t i represents the i-th discrete time delay; A 0 represents the state matrix; represents the observation matrix; a p represents the axial cutting depth; T represents the time delay; Represents the integral factor. According to the intermittent process of tooth cutting, the time lag is divided into a vibration-free stage and a self-excited vibration stage. The time interval of self-excited vibration is discretized to construct the transfer function. The reconstructed Newton-Cotes method is introduced into each discretized time interval to perform interpolation fitting on the integral part of the state space equation, and the following interpolation formula is obtained:

[0038] ,

[0039] Where, X i represents the state variable of the i-th integration interval; X i+1 represents the state variable of the i+1th integration interval; X i+2 represents the state variable of the i+2th integration interval; represents the i-th state vector; represents the i+1 / 4th state vector; represents the i+1 / 2th state vector; represents the i+3 / 4th state vector; represents the i+1th state vector; A 0 represents the state matrix; h represents the discretized time-delay interval; represents the kth barycentric Lagrange interpolation of the jth integration interval; a 0 、a 1 、a 2 、b 0 、b 1 、b 2 、c 0 、c 1 、c 2 represents the weight coefficient; t j and t k They represent the jth and kth discrete points in the center-of-gravity Lagrangian interpolation at time t. According to the transfer relationship of the state variables after discretization, the mapping relationship of all state variables in a single tooth-passing cycle can be expressed as follows:

[0040] ,

[0041] Where, I represents the identity matrix; C represents the first solution module matrix about the state matrix; D represents the second solution module matrix about the state matrix and the observation matrix; E represents the third solution module matrix about the free vibration stage; a p represents the axial cutting depth; Δt represents the discrete post-delay; Represents the first state vector; Represents the second state vector; Represents the n+1th state vector. Based on the mapping relationship of the state vector, the transfer matrix expression of a single tooth-through cycle can be obtained as follows:

[0042] ,

[0043] ,

[0044] ,

[0045] ,

[0046] In the formula, represents the transfer matrix of a single tooth-through cycle; I represents the identity matrix; C represents the first solution module matrix about the state matrix; D represents the second solution module matrix about the state matrix and the observation matrix; E represents the third solution module matrix about the free vibration stage; a p represents the axial cutting depth; Δt represents the time lag after discreteness; A 0 represents the state matrix; represents the observation matrix about the i-th discrete interval; a 0 、a 1 、a 2 、b 0 、b 1 、b 2 、c 0 、c 1 、c 2 represents the weight coefficient; t f Represents the free vibration time interval.

[0047] Furthermore, the method of reducing the dimension based on the transfer matrix and configuring a three-dimensional stability lobe diagram drawing method to predict the milling stability includes:

[0048] First, the transfer matrix expression of a single tooth pass cycle is Rewrite as , the dimensions of H matrix and R matrix are ;

[0049] Next, remove the 4i-1th to 4ith column vectors of the R matrix, where i=1,2,…,n / 2; Matrix, dimension is ;

[0050] Secondly, calculate Matrix and The matrix product gives ;

[0051] Secondly, remove the 4i-1th to 4ith column vectors of the L matrix, where i=1,2,…,n / 2; and obtain the reduced dimension Matrix, dimension is ;

[0052] Finally, the eigenvalues ​​of the transfer matrix after dimensionality reduction are calculated and the stability is evaluated; the stable domain and the unstable domain are divided, and a three-dimensional stability lobe diagram is output. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] In order to more clearly illustrate the technical solution of the present invention, the drawings required for use in the embodiments are briefly introduced below. Obviously, the drawings described below are only embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0054] Figure 1 Schematic diagram of milling processing and discrete particles of thin-walled parts provided by an embodiment of the present invention;

[0055] Figure 2 A simplified diagram of the modal vibration shapes of a milling cutter and a thin-walled part provided in an embodiment of the present invention;

[0056] Figure 3 A schematic diagram of mesh division of a thin-walled part provided in an embodiment of the present invention;

[0057] Figure 4 A real part data diagram of the particle frequency response function of a milling cutter and a thin-walled part provided in an embodiment of the present invention;

[0058] Figure 5 A schematic diagram of the instantaneous cutting process of a cross section of a tool provided in an embodiment of the present invention;

[0059] Figure 6 A flow chart of transfer matrix reduction and dimensionality reduction provided by an embodiment of the present invention;

[0060] Figure 7 A two-dimensional stability lobe diagram of a thin-walled component provided by an embodiment of the present invention;

[0061] Figure 8 A three-dimensional stability lobe diagram of a thin-walled component provided by an embodiment of the present invention;

[0062] Fig. 9A cross section of a three-dimensional stability lobe diagram of a thin-walled part provided by an embodiment of the present invention with a main shaft speed of 5000 rpm;

[0063] Fig.10 The actual measured acceleration and its spectrum data diagram at a cutting depth of 1 mm provided by an embodiment of the present invention;

[0064] Fig.11 The actual measured acceleration and its spectrum data diagram at a cutting depth of 1.3 mm provided by an embodiment of the present invention;

[0065] Fig.12 The actual measured acceleration and its spectrum data diagram at a cutting depth of 1.5 mm provided in an embodiment of the present invention. DETAILED DESCRIPTION

[0066] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0067] The present invention provides a method for predicting the stability of multi-particle and multi-modal milling of thin-walled parts considering the modal coupling effect, including: constructing an initial milling dynamics model including feed direction cutting force and radial cutting force excitation. During the milling process of thin-walled parts, periodic cutting forces are generated under the influence of cutting excitation. When the frequency of periodic changes in cutting force is close to the natural frequency of thin-walled parts, machining vibration, i.e., chatter, is very likely to occur. Assuming that the modal parameters and the spindle speed remain unchanged during the cutting process, a two-degree-of-freedom lumped parameter model including feed direction cutting force and radial cutting force excitation is constructed. The expression of the initial milling dynamics model including feed direction cutting force and radial cutting force excitation is:

[0068] ,

[0069] Where M is the mass matrix in Cartesian space; C is the damping matrix in Cartesian space; K is the stiffness matrix in Cartesian space; is the vibration acceleration vector at time t; is the vibration velocity vector at time t; is the vibration displacement vector at time t; is the cutting force vector of the milling cutter at time t; is the axial cutting depth; is the number of milling cutter teeth; is the rotation matrix in Cartesian space; A window function to determine whether the cutter tooth is in the cutting state; is the cutting force coefficient matrix; is the rotation matrix of the jth tooth direction at the i-th node of the milling cutter at time t; T is the hysteresis period.

[0070] The discrete tool along the axial direction is N t Mass points, discrete thin-walled parts N along the feed direction wp The dynamic parameters of each particle are obtained, and the first milling dynamic model of the interface between the milling cutter and the thin-walled workpiece at different processing positions is established according to the initial milling dynamic model to simulate the changing trend of the stability state of the milling cutter cutting process after discrete cutting.

[0071] Since the stiffness of thin-walled parts in the wall thickness direction is low and the milling cutter is subjected to the cutting force in the feed direction, chatter usually occurs in the wall thickness direction and the feed direction. During the clamping process, thin-walled parts are approximately considered to be a cantilever beam structure with strong stiffness near the clamping part and weak stiffness away from the clamping part. The milling process and discrete particles are shown in the figure. Figure 1 As shown in the figure, during the milling process, the dynamic parameters of the thin-walled part change with the machining position. A multi-particle multi-modal milling dynamics model is constructed to accurately describe the changes in the stability state of the machining process.

[0072] Using the finite element method, the overhanging tool is discretized along the axial direction. The discrete tool is N t The discrete thin-walled part along the milling cutter feed direction is N wp particles, divided into N on the thin-walled part t ×N wp The dynamic parameters of each particle are obtained respectively, and the multi-particle multi-mode milling dynamic equations are constructed based on the initial milling dynamics model:

[0073] ,

[0074] Where M tool is the tool quality matrix after discretization; C tool is the discretized tool damping matrix; K tool is the tool stiffness matrix after discretization; M wp is the mass matrix of the thin-walled part after discretization; C wp is the damping matrix of the discretized thin-walled component; K wp is the stiffness matrix of the thin-walled component after discretization; is the vibration acceleration vector of the milling cutter at time t after the discretization; is the vibration velocity vector of the milling cutter at time t after the discretization; is the vibration displacement vector of the milling cutter at time t after discretization; is the vibration acceleration vector of the thin-walled part at time t after discretization; is the vibration velocity vector of the thin-walled part at time t after discretization; is the vibration displacement vector of the thin-walled part at time t after discretization; is the cutting force vector of the milling cutter at the discretized moment t.

[0075] The mass matrix, damping matrix and stiffness matrix in the discretized milling dynamics system are related to the dynamic parameters of each mass point. The milling dynamics system of the milling cutter and the thin-walled workpiece is coupled into an overall dynamics system about the mass point. The first milling dynamics model of the interface between the milling cutter and the thin-walled workpiece obtained by coupling is:

[0076] .

[0077] Based on the modal superposition theory, the vibration displacement in the milling dynamics model is converted from Cartesian space to modal space. Considering the coupling effect between the low-order modes and high-order modes of each particle, the multi-order modal vibration matrix is ​​introduced to correct the milling dynamics equation, and the second milling dynamics model considering the modal coupling effect is obtained.

[0078] like Figure 2 As shown in the figure, the low-order modes and high-order modes of the milling cutter and the thin-walled workpiece far away from the clamping position of the tooling may produce modal coupling phenomenon, and the change of its modal vibration shape will affect the prediction accuracy of the stability state in the actual processing process. At the same time, the coupling effect at different particle positions is different, resulting in different contributions of low-order modes and high-order modes to processing stability. Therefore, it is necessary to calculate the modal coupling vibration shape through the modal superposition principle and construct a second milling dynamics model considering the modal coupling effect. According to the modal superposition theory, the conversion relationship of the vibration displacement vector from Cartesian space to modal space is calculated:

[0079] ,

[0080] In the formula, represents the vibration displacement expression in the modal space, and U represents the modal vibration matrix. The modal vibration matrix of the milling cutter and thin-walled parts is shown as follows:

[0081] ,

[0082] Where U tool Represents the milling cutter mode shape; U wp Represents the modal vibration matrix of thin-walled parts; and represent the j-order mode shape vectors in the x-direction and the y-direction respectively.

[0083] The multi-order modal vibration matrix of the low-order mode and high-order mode coupling effect of the milling cutter and the workpiece is shown as follows:

[0084] ,

[0085] ,

[0086] ,

[0087] ,

[0088] In the formula, represents the i-th vibration mode in the x direction; represents the i-th vibration mode in the y direction.

[0089] Based on the conversion relationship between Cartesian space and modal space and the principle of modal coupling, the dynamic expression of multi-particle multi-modal milling can be rewritten as:

[0090] ,

[0091] In the formula, represents the modal mass matrix; represents the modal damping matrix; represents the modal stiffness matrix; represents the vibration acceleration at time t in the modal space; represents the vibration velocity at time t in the modal space; represents the vibration displacement at time t in the modal space; a p represents the axial cutting depth; U represents the overall modal vibration matrix; U T represents the transpose of the overall modal vibration matrix; U wp Represents the modal vibration matrix of thin-walled parts; represents the transpose of the modal vibration matrix of the thin-walled part; U tool represents the milling cutter mode shape; represents the transpose of the milling cutter modal shape matrix; Represents the milling force coefficient matrix of the jth infinitesimal element of the i-th tooth at time t.

[0092] In order to illustrate the construction process of the multi-modal milling dynamic model of thin-walled parts with modal coupling effect, the milling of aluminum alloy 7075 thin-walled parts is taken as an example. The size of the selected thin-walled parts is 100 mm × 70 mm × 5 mm, the Young's modulus of the material is 68.9 GPa, and the Poisson's ratio is 0.33. The milling cutter is made of cemented carbide material, with 4 teeth, a diameter of 12 mm, a helix angle of 35°, and a tool overhang of 75 mm. The processing parameters are spindle speed of 5000 rpm, cutting depth of 0.2 mm, and feed rate of 0.04 mm / tooth.

[0093] First, the thin-walled part is divided into 10 equal parts along the machining feed direction. Based on the mesh division, the mass points are selected as follows: Figure 3 The first-order and second-order natural frequencies, stiffness and damping ratios of thin-walled parts in the wall thickness direction obtained through modal tests are shown in Table 1:

[0094]

[0095] Secondly, the two-dimensional natural frequency, stiffness and damping ratio of the milling cutter are obtained along the cutter axis direction as shown in Table 2:

[0096]

[0097] In order to determine whether there is modal coupling in the dynamic parameters of the tool and each mass point of the thin-walled workpiece, it is necessary to estimate the frequency response function by least squares. Taking the milling cutter x direction and the fourth mass point in the thin-walled workpiece as an example, the superposition result of the first-order frequency response and the second-order frequency response is calculated by the modal superposition principle, as shown in Figure 4 It can be seen that there is no coupling phenomenon between the real part curves of the low-order and high-order modal frequency response functions on the tool side, and the coupling phenomenon between the real part curves of the low-order and high-order modal frequency response functions of the thin-walled part is obvious near the high-order frequency.

[0098] In order to accurately construct the model, it is necessary to accurately identify the dynamic parameters of the milling system. The milling force coefficient is identified through slot milling cutting experiments. The cutting force coefficient of aluminum alloy 7075 is K t =159.8N / mm 2 , K n =66.5 N / mm 2 .

[0099] Exemplarily, the method for calculating the instantaneous cutting angle and cutting angle changes when the helix angle of the milling cutter acts on the cutter teeth and the micro-element cutting surface of the thin-walled workpiece includes:

[0100] like Figure 5 As described above, according to the helix angle between the cutter teeth and the micro-element cutting surface of the thin-walled workpiece, the instantaneous cutting angle model is constructed as follows:

[0101] ,

[0102] In the formula, represents the instantaneous cutting angle at time t; Ω represents the spindle speed; i represents the i-th tooth; j represents the j-th discrete element along the axial direction of the milling cutter; m represents the axial discrete number of the milling cutter; N f Indicates the number of milling cutter teeth; a p It represents the axial cutting depth; β represents the helix angle of the milling cutter; and D represents the diameter of the milling cutter.

[0103] Exemplarily, according to the second milling dynamics model, a method for establishing a state space equation with periodic time lag includes:

[0104] make , then the second milling dynamics model considering the modal coupling effect can be transformed into a state space expression as follows:

[0105] ,

[0106] In the formula, Represents state variables; represents the vibration velocity at time t in the modal space; represents the vibration displacement at time t in the modal space; represents the first-order derivative of the state variable; A 0 represents the state matrix; B(t) represents the observation matrix; a p Indicates the axial cutting depth; represents the inverse of the modal mass; represents the modal stiffness; represents modal damping; I represents the unit matrix; U represents the overall modal vibration matrix; represents the transpose of the overall mode shape matrix; Represents the milling force coefficient matrix of the jth infinitesimal element of the i-th tooth at time t.

[0107] Exemplarily, according to the periodic lag term of the discretized state-space equation, a method of constructing a reconstructed Newton-Cotes fitting formula and calculating the transfer matrix includes:

[0108] According to the numerical solution of the Volterra equation of the second kind, the state space equation can be solved as follows:

[0109] ,

[0110] In the formula, Represents state variables; represents the i-th discretized state vector; t i represents the i-th discrete time delay; A 0 represents the state matrix; represents the observation matrix; a p represents the axial cutting depth; T represents the time delay; Represents the integrating factor.

[0111] According to the intermittent process of tooth cutting, the time lag is divided into a vibration-free stage and a self-excited vibration stage. The time interval of self-excited vibration is discretized to construct the transfer function. The reconstructed Newton-Cotes method is introduced into each discretized time interval to perform interpolation fitting on the integral part of the state space equation, and the following interpolation formula is obtained:

[0112] ,

[0113] Where, X i represents the state variable of the i-th integration interval; X i+1 represents the state variable of the i+1th integration interval; X i+2 represents the state variable of the i+2th integration interval; represents the i-th state vector; represents the i+1 / 4th state vector; represents the i+1 / 2th state vector; represents the i+3 / 4th state vector; represents the i+1th state vector; A 0 represents the state matrix; h represents the discretized time-delay interval; represents the kth barycentric Lagrange interpolation of the jth integration interval; a 0 、a 1 、a 2 、b 0 、b 1 、b 2 、c 0 、c 1 、c 2 represents the weight coefficient; t j and t k They represent the j-th and k-th discrete points in the barycentric Lagrangian interpolation with respect to time t, respectively.

[0114] According to the transfer relationship of the discrete state variables, the mapping relationship of all state variables in a single tooth-passing cycle can be expressed as follows:

[0115] ,

[0116] Where, I represents the identity matrix; C represents the first solution module matrix about the state matrix; D represents the second solution module matrix about the state matrix and the observation matrix; E represents the third solution module matrix about the free vibration stage; a p represents the axial cutting depth; Δt represents the discrete post-delay; Represents the first state vector; Represents the second state vector; Represents the n+1th state vector.

[0117] Based on the mapping relationship of the state vector, the transfer matrix expression of a single tooth pass cycle can be obtained as follows:

[0118] ,

[0119] ,

[0120] ,

[0121] ,

[0122] In the formula, represents the transfer matrix of a single tooth-through cycle; I represents the identity matrix; C represents the first solution module matrix about the state matrix; D represents the second solution module matrix about the state matrix and the observation matrix; E represents the third solution module matrix about the free vibration stage; a p represents the axial cutting depth; Δt represents the time lag after discreteness; A 0 represents the state matrix; represents the observation matrix about the i-th discrete interval; a 0 、a 1 、a 2 、b 0 , b 1 、b 2 、c 0 、c 1 、c 2 represents the weight coefficient; t f Represents the free vibration time interval.

[0123] Exemplarily, the method for predicting milling stability by reducing the dimension based on the transfer matrix and configuring a three-dimensional stability lobe diagram drawing method includes:

[0124] like Figure 6 As shown, first, the transfer matrix expression of a single tooth pass cycle is Rewrite as , the dimensions of H matrix and R matrix are .

[0125] Next, remove the 4i-1th to 4ith column vectors of the R matrix, where i=1,2,…,n / 2; Matrix, dimension is .

[0126] Secondly, calculate Matrix and The matrix product gives ;

[0127] Secondly, remove the 4i-1th to 4ith column vectors of the L matrix, where i=1,2,…,n / 2; and obtain the reduced dimension Matrix, dimension is .

[0128] Finally, the eigenvalues ​​of the transfer matrix after dimensionality reduction are calculated to predict and evaluate the stability; the stable domain and the unstable domain are divided, and a three-dimensional stability lobe diagram is output.

[0129] In order to illustrate the implementation process of the stability prediction method for thin-walled workpiece milling considering the modal coupling effect, the stability is derived and analyzed based on the established multi-particle multi-modal milling dynamics model of thin-walled workpieces considering the modal coupling effect and the machining dynamics parameters of the milling cutter and thin-walled workpiece, and the two-dimensional stability lobe diagram of particle 4 is output, as shown in Figure 2. Figure 7 As shown. It can be seen that with the increase of spindle speed, its stability boundary generally shows a trend of first increasing and then decreasing. The solid line represents the processing stability boundary considering the modal coupling effect, and the dotted line represents the processing stability boundary considering only the first-order mode. Taking the spindle speed of 5000rpm as an example, the stability cutting depth considering the modal coupling effect of thin-walled parts is lower than the stability cutting depth considering only the single modal effect. This is because in the cutting process, the high-order natural frequency and the low-order natural frequency jointly play a dominant role, resulting in aggravated cutting instability.

[0130] To further verify the accuracy of the stability prediction method, the matrix reduction method is used to achieve rapid calculation and storage of eigenvalues, and output a three-dimensional stability lobe diagram, such as Figure 8 Select the section with a rotation speed of 5000rpm on the three-dimensional stability lobe diagram, as shown in Fig. 9 As shown in the figure, it can be seen that the stability boundary increases first and then decreases with the increase of cutting position. This is because the clamping method of thin-walled parts can be approximated as a symmetrical cantilever beam structure. As the cutting position approaches the symmetry center, its stiffness gradually increases, and its anti-vibration ability is also enhanced accordingly. To further verify the accuracy of the stability boundary, Figure 10-12 Three groups of axial cutting depths of 1.0mm, 1.3mm and 1.5mm were selected for milling experiments. The experiment was carried out on a vertical machining center, and the size and shape of the machining tool and thin-walled parts were the same as in the first example. During the experiment, an acceleration sensor was used to collect the acceleration signal of the thin-walled part, and the frequency spectrum of the acceleration was output through fast Fourier transform, as shown in Fig.10 As shown in the figure. At a cutting depth of 1mm, there is no chatter component in the acceleration spectrum, which belongs to the stable working condition. Fig.11 As shown in the figure, the acceleration spectrum at a cutting depth of 1.3 mm contains a small amount of chatter components, which belongs to a mild chatter condition. Fig.12 As shown in the figure, the acceleration spectrum at a cutting depth of 1.5 mm contains a large amount of chatter components, which belongs to a severe chatter condition. The three groups of experiments verified the machining stability state under these machining parameters, which is consistent with the stability lobe diagram results of the thin-walled part milling stability prediction considering the modal coupling effect proposed in the present invention, indicating the accuracy of the proposed stability prediction method.

[0131] In summary, compared with the prior art, the present invention proposes a method for predicting the stability of multi-particle multi-modal milling of thin-walled parts that takes into account the modal coupling effect. It comprehensively considers the characteristics of the change of dynamic parameters with the processing position during the milling of thin-walled parts and the influence of the coupling effect of low-order modes and high-order modes on the boundary of processing stability, and optimizes the modeling method of the milling dynamics of thin-walled parts. By reconstructing the Newton-Cotes formula and the reduction and dimensionality reduction method of the transfer matrix, the dynamic model containing the coupled modal vibration matrix can be quickly solved and the stability analysis can be achieved.

[0132] The present invention has been described in detail above in conjunction with specific implementations and exemplary examples, but these descriptions cannot be understood as limiting the present invention. Those skilled in the art understand that, without departing from the spirit and scope of the present invention, a variety of equivalent substitutions, modifications or improvements may be made to the technical solution of the present invention and its implementation methods, all of which fall within the scope of the present invention. The scope of protection of the present invention shall be subject to the attached claims.

Claims

1. A stability prediction method for multi-particle multi-modal milling of thin-walled parts considering modal coupling effect, the method comprising: Construct an initial milling dynamics model including cutting force in feed direction and radial cutting force excitation; The discrete tool along the axial direction is N t The discrete thin-walled parts along the feed direction are N wp mass points, obtain the dynamic parameters of each mass point, and establish the first milling dynamic model of the interface between the milling cutter and the thin-walled workpiece at different processing positions in the feed direction according to the initial milling dynamic model, so as to simulate the change trend of the stability state of the milling cutter cutting process after discrete cutting; The vibration displacement in the first milling dynamics model is converted from Cartesian space to modal space, the coupling effect between low-order modes and high-order modes of each particle is considered, and a multi-order modal vibration matrix is ​​introduced to correct the first milling dynamics model to obtain a second milling dynamics model considering the modal coupling effect; Calculate the instantaneous cutting angle and cutting angle changes when the milling cutter helix angle acts on the cutter teeth and the micro-element cutting surface of thin-walled parts; According to the second milling dynamics model, a state space equation with periodic time lag is established; According to the periodic time-delay term of the discretized state-space equation, the reconstructed Newton-Cotes fitting formula is constructed and the transfer matrix is ​​calculated; Based on the transfer matrix, dimensionality reduction is performed, a three-dimensional stability lobe diagram drawing method is configured, and milling stability is predicted.

2. According to the method for predicting the stability of multi-particle multi-modal milling of thin-walled parts considering modal coupling effect in claim 1, the step of constructing an initial milling dynamics model including feed direction cutting force and radial cutting force excitation specifically comprises: Construct the initial milling dynamics model expression including the cutting force in the feed direction and the radial cutting force excitation: ， Where M is the mass matrix in Cartesian space; C is the damping matrix in Cartesian space; K is the stiffness matrix in Cartesian space; is the vibration acceleration vector at time t; is the vibration velocity vector at time t; is the vibration displacement vector at time t; is the cutting force vector of the milling cutter at time t; is the axial cutting depth; is the number of milling cutter teeth; is the rotation matrix in Cartesian space; A window function to determine whether the cutter tooth is in the cutting state; is the cutting force coefficient matrix; is the rotation matrix of the jth tooth direction at the i-th node of the milling cutter at time t; T is the hysteresis period.

3. According to the method for predicting the stability of multi-particle multi-modal milling of thin-walled parts considering modal coupling effect in claim 2, the axially discrete tool is N t The discrete thin-walled parts along the feed direction are N wp The steps of obtaining the dynamic parameters of each particle point, and establishing the first milling dynamic model of the interface between the milling cutter and the thin-walled workpiece at different processing positions in the feed direction according to the initial milling dynamic model to simulate the change trend of the stability state of the milling cutter cutting process after discrete cutting specifically include: The discrete tool along the axial direction is N t The discrete thin-walled parts along the feed direction are N wp particles, forming N t ×N wp The mesh is used to obtain the dynamic parameters of each particle, and the milling dynamics equations are constructed based on the initial milling dynamics model: ， Where M tool is the tool quality matrix after discretization; C tool is the discretized tool damping matrix; K tool is the tool stiffness matrix after discretization; M wp is the mass matrix of the thin-walled part after discretization; C wp is the damping matrix of the discretized thin-walled component; K wp is the stiffness matrix of the thin-walled component after discretization; is the vibration acceleration vector of the milling cutter at time t after the discretization; is the vibration velocity vector of the milling cutter at time t after the discretization; is the vibration displacement vector of the milling cutter at time t after discretization; is the vibration acceleration vector of the thin-walled part at time t after discretization; is the vibration velocity vector of the thin-walled part at time t after discretization; is the vibration displacement vector of the thin-walled part at time t after discretization; is the cutting force vector of the milling cutter at the discrete moment t; where, The mass matrix, damping matrix and stiffness matrix in the discretized milling dynamics system are coupled with the dynamic parameters of each particle to obtain the overall dynamic model of the particle: 。 4. According to the method for predicting the stability of multi-particle multi-modal milling of thin-walled parts considering the modal coupling effect of claim 3, the steps of converting the vibration displacement in the first milling dynamics model from Cartesian space to modal space, considering the coupling effect between the low-order mode and the high-order mode of each particle, introducing a multi-order modal vibration matrix to correct the first milling dynamics equation, and obtaining the second milling dynamics model considering the modal coupling effect specifically include: According to the modal superposition theory, the transformation relationship of the vibration displacement vector from Cartesian space to modal space is calculated as follows: ， In the formula, represents the vibration displacement expression in the modal space, and U represents the modal vibration matrix; The modal vibration matrix of the milling cutter and thin-walled workpiece is expressed as: ， Where U tool Represents the milling cutter mode shape; U wp Represents the modal shape matrix of thin-walled parts; and represent the j-order mode shape vectors in the x-direction and the y-direction respectively; The multi-order modal vibration matrix of the low-order mode and high-order mode coupling effect of the milling cutter and the workpiece is shown as follows: ， ， ， ， In the formula, represents the i-th vibration mode in the x direction; represents the i-th vibration mode in the y direction; Based on the conversion relationship between Cartesian space and modal space and the modal coupling principle, the first milling dynamics model expression is rewritten into the second milling dynamics model as follows: ， In the formula, represents the modal mass matrix; represents the modal damping matrix; represents the modal stiffness matrix; represents the vibration acceleration at time t in the modal space; represents the vibration velocity at time t in the modal space; represents the vibration displacement at time t in the modal space; a p represents the axial cutting depth; U represents the overall modal vibration matrix; U T represents the transpose of the overall modal vibration matrix; U wp Represents the modal shape matrix of thin-walled parts; represents the transpose of the modal vibration matrix of the thin-walled part; U tool represents the mode shape of the milling cutter; represents the transpose of the milling cutter mode shape matrix; Represents the milling force coefficient matrix of the jth infinitesimal element of the i-th tooth at time t.

5. According to the method for predicting the stability of multi-particle multi-modal milling of thin-walled parts considering the modal coupling effect of claim 4, the step of calculating the instantaneous cutting angle and cutting angle changes when the helix angle of the milling cutter acts on the cutter teeth and the micro-element cutting surface of the thin-walled part is: According to the helix angle between the cutter teeth and the micro-element cutting surface of the thin-walled workpiece, the instantaneous cutting angle model is constructed as follows: ， In the formula, represents the instantaneous cutting angle at time t; Ω represents the spindle speed; i represents the i-th tooth; j represents the j-th discrete element along the axial direction of the milling cutter; m represents the axial discrete number of the milling cutter; N f Indicates the number of milling cutter teeth; a p It represents the axial cutting depth; β represents the helix angle of the milling cutter; and D represents the diameter of the milling cutter.

6. According to the method for predicting the stability of multi-particle multi-modal milling of thin-walled parts considering the modal coupling effect of claim 5, the step of establishing a state space equation with periodic time lag based on the multi-particle multi-modal milling dynamics model of thin-walled parts considering the modal coupling effect is: make , The second milling dynamics model considering the modal coupling effect is transformed into a state space expression as follows: ， In the formula, Represents state variables; represents the vibration velocity at time t in the modal space; represents the vibration displacement at time t in the modal space; represents the first-order derivative of the state variable; A0 represents the state matrix; B(t) represents the observation matrix; a p Indicates the axial cutting depth; represents the inverse of the modal mass; represents the modal stiffness; represents modal damping; I represents the unit matrix; U represents the overall modal vibration matrix; represents the transpose of the overall mode shape matrix; Represents the milling force coefficient matrix of the jth infinitesimal element of the i-th tooth at time t.

7. According to the method for predicting the stability of multi-particle multi-modal milling of thin-walled parts considering modal coupling effect in claim 6, the method for constructing and reconstructing the Newton-Cotes fitting formula and calculating the transfer matrix according to the periodic lag term of the discretized state space equation comprises: According to the numerical solution of the Volterra equation of the second kind, the state space equation is solved as follows: ， In the formula, Represents state variables; represents the i-th discretized state vector; t i represents the i-th discrete time delay; A0 represents the state matrix; represents the observation matrix; a p represents the axial cutting depth; T represents the time delay; represents the integrating factor; According to the intermittent process of tooth cutting, the time lag is divided into a vibration-free stage and a self-excited vibration stage. The time interval of self-excited vibration is discretized to construct a transfer function. The reconstructed Newton-Cotes method is introduced into each discretized time interval to perform interpolation fitting on the integral part of the state space equation, and the following interpolation formula is obtained: ， Where, X i represents the state variable of the i-th integration interval; X i+1 represents the state variable of the i+1th integration interval; X i+2 represents the state variable of the i+2th integration interval; represents the i-th state vector; represents the i+1 / 4th state vector; represents the i+1 / 2th state vector; represents the i+3 / 4th state vector; represents the i+1th state vector; A0 represents the state matrix; h represents the discretized time lag interval; represents the kth centroid Lagrange interpolation of the jth integral interval; a0, a1, a2, b0, b1, b2, c0, c1, c2 represent weight coefficients; t j and t k They represent the jth and kth discrete points about time t in the barycentric Lagrangian interpolation respectively; According to the transfer relationship of the discrete state variables, the mapping relationship of all state variables in a single tooth-through cycle is expressed as follows: ， Where, I represents the identity matrix; C represents the first solution module matrix about the state matrix; D represents the second solution module matrix about the state matrix and the observation matrix; E represents the third solution module matrix about the free vibration stage; a p represents the axial cutting depth; Δt represents the time lag after discreteness; Represents the first state vector; Represents the second state vector; represents the n+1th state vector; Based on the mapping relationship of the state vector, the transfer matrix expression of a single tooth-through cycle is obtained as follows: ， ， ， ， In the formula, represents the transfer matrix of a single tooth-through cycle; I represents the identity matrix; C represents the first solution module matrix about the state matrix; D represents the second solution module matrix about the state matrix and the observation matrix; E represents the third solution module matrix about the free vibration stage; a p represents the axial cutting depth; Δt represents the discrete post-delay; A0 represents the state matrix; represents the observation matrix for the ith discrete interval; a0, a1, a2, b0, b1, b2, c0, c1, c2 represent weight coefficients; t f Represents the free vibration time interval.

8. According to the method for predicting the stability of multi-particle multi-modal milling of thin-walled parts considering the modal coupling effect of claim 7, the step of reducing the dimension based on the transfer matrix and configuring a three-dimensional stability lobe diagram drawing method to predict the milling stability comprises: Step 1: Transform the transfer matrix expression of a single tooth cycle into Rewrite as , then the dimensions of the H matrix and the R matrix are ; Step 2, remove the 4i-1th to 4ith column vectors of the R matrix, where i=1,2,…,n / 2; get Matrix, dimension is ; Step 3, calculate Matrix and The matrix product gives ; Step 4: Remove the 4i-1th to 4ith column vectors of the L matrix, where i=1,2,…,n / 2; and obtain the reduced dimension Matrix, dimension is ; Step 5, calculate the eigenvalues ​​of the transfer matrix after dimensionality reduction, and predict and evaluate the milling stability.

Images