Robot grinding chatter prediction method considering discrete vibration speed

Abstract

The invention provides a robot grinding chatter prediction method considering discrete vibration speed. The method comprises the following steps: establishing a grinding dynamics differential equation of a composite five-degree-of-freedom operation robot; converting the grinding dynamics differential equation of the combined type five-degree-of-freedom operation robot into a time-delay grinding dynamics equation in a state space form; discretizing the delay time of the time-delay grinding kinetic equation, and obtaining a solution of the time-delay grinding kinetic equation in a corresponding time interval through direct integration; taking a time delay part and a state item in the expression of the solution as a whole, and processing by adopting third-order Hemet interpolation to obtain a discrete dynamic iterative formula; constructing a state transition matrix through a discrete dynamic iterative formula; and the flutter stability is judged through the Floquet theory. According to the method, the influence of the discrete vibration displacement and the discrete vibration speed is considered at the same time, the prediction error generated by independent discretization of the vibration state variable is avoided, and the prediction precision is high.

Description

Technical Field

[0001] The present invention relates to the technical field of chatter prediction in a grinding process, and in particular to a robot grinding chatter prediction method considering discrete vibration speeds. Background Art

[0002] Industrial robots are widely used in the machining of complex components due to their flexibility, wide workspace, and low cost. However, the relatively weak rigidity of the tandem structure of industrial robots during machining makes chatter more likely, leading to poor workpiece surface quality, reduced machining efficiency, and shortened tool life.

[0003] Predicting and avoiding grinding chatter not only allows industrial robots to operate in a stable machining state, achieving high-quality machined surfaces to meet the demands of precision manufacturing, but also reduces tool wear and extends the life of the equipment. Further research into chatter stability during grinding processes is of great significance for improving machining accuracy, stabilizing the operation of industrial robots, and extending the life of equipment.

[0004] The essence of time-domain discretization methods lies in using numerical methods to calculate the state transition matrix of a dynamic differential equation (DDE) with periodic coefficients over the entire cycle, thereby constructing an equivalent discrete dynamic system that approximates the original grinding dynamics over that cycle. Specifically, these methods involve converting the DDE into a state-space representation. The entire cycle is then divided into a series of small time intervals. Subsequently, the state terms, time delay terms, and periodic coefficient matrices in the integral form of the DDE solution are interpolated and partially or fully approximated using various numerical methods. The resulting recursive formula is transposed and combined to construct the state transition matrix for the entire cycle. According to Floquet theory, the spectral radius of the state transition matrix is ​​compared with 1 to evaluate the stability of the machining system.

[0005] In the time domain method, the state transfer matrix is ​​established and the Floquet theory is used to determine the chatter stability. However, only the discrete vibration displacement is usually considered, while the influence of the discrete vibration velocity is ignored. The present invention proposes a robot grinding chatter prediction method that considers the discrete vibration velocity, which avoids the prediction error caused by the separate discretization of the vibration state variables and has higher prediction accuracy. Summary of the Invention

[0006] To address the above issues, the present invention proposes a robot grinding chatter prediction method considering discrete vibration velocities, comprising the following steps:

[0007] Step S1: establishing the grinding dynamics differential equation of the composite five-degree-of-freedom working robot;

[0008] Step S2: transforming the grinding dynamics differential equation of the composite five-DOF working robot into a time-delay grinding dynamics equation in state space form;

[0009] Step S3: For the time-delay grinding dynamics equation in the state space form, first discretize the delay time, and then obtain the approximate solution of the equation in the corresponding time interval by direct integration;

[0010] Step S4: The time lag part and the state term in the solution of the time lag grinding dynamics equation are treated as a whole, and processed using third-order Hermite interpolation to obtain a discrete dynamic iterative formula;

[0011] Step S5: constructing a state transfer matrix through a discrete dynamic iterative formula;

[0012] Step S6: Determine the flutter stability using Floquet theory.

[0013] In one embodiment, step S1 includes:

[0014] Step S11: Establish the DH reference coordinate system and link parameter table of the composite five-DOF working robot, and obtain the following homogeneous transformation matrix of adjacent coordinate systems:

[0015]

[0016]

[0017]

[0018]

[0019]

[0020] Where, is the rotation angle of J1~J5 axis; α is the rotation angle of the turntable; θ i is the joint variable; d i is the joint offset; a i Connecting rod length; X1, X2 and X3 are the travels of the nut on the ball screws S1, S2 and S3, robot structural parameters e and tool compensation L P All are considered as 0.

[0021] Step S12: Based on the above homogeneous transformation matrix, the homogeneous transformation matrix of the coordinate system {O5} of the end effector of the composite five-degree-of-freedom working robot relative to the base coordinate system {O0} can be obtained:

[0022]

[0023] Where,

[0024]

[0025]

[0026]

[0027]

[0028] Step S13: Based on the obtained homogeneous transformation matrix, the Jacobian matrix can be directly obtained by using the vector product method without the need for derivation. The expression for solving the Jacobian matrix using the vector product method is:

[0029]

[0030] Where, J Pi 、J Oi is the i-th column of the Jacobian matrix of the robot's moving joints and rotating joints; P e is the position vector of the robot's end effector relative to the base coordinate system {O0}, which can be represented by the homogeneous transformation matrix The first three elements of the fourth column are obtained; P i-1 It can be represented by the homogeneous transformation matrix The first three elements of the fourth column are obtained; z i-1 Can be obtained by the rotation matrix The third column of is obtained.

[0031] Solving the J1 to J5 columns of the Jacobian matrix, we can obtain:

[0032] J1=[1 0 0 0 0 0] T

[0033] J2=[-(L3-L3cosα)L3sinα 0 00 1] T

[0034] J3=[0 0 1 0 0 0] T

[0035] J4=[0 0 0 0 0 1] T

[0036]

[0037] Therefore, the Jacobian matrix of the composite five-DOF operating robot is:

[0038] J=[J1 J2 J3 J4 J5]

[0039] Step S14: Establish the robot dynamics differential equation:

[0040]

[0041] Where M, C, and K are the 5×5 mass matrix, damping matrix, and stiffness matrix of the robot, respectively; x is the displacement vector; and F is the external force vector.

[0042] Where,

[0043] and,

[0044] Where c ij represents the element in row i and column j of the damping matrix C; Indicates the movement speed of the robot joint; coefficient c ijk It is called the Christoffel symbol of the first form and is given by:

[0045]

[0046] Where b ij represents the element in the i-th row and j-th column of the mass matrix, q k represents the joint angles of the robot, Represents the partial derivative of the elements in the mass matrix with respect to the joint angles of the robot.

[0047] Calculate the stiffness matrix of the composite five-DOF operating robot, and establish the static stiffness model of the robot based on the Jacobian matrix and Hooke's law:

[0048] F=J -T K θ J -1 X

[0049] Where F represents the external force vector of the robot during machining; K θ represents the stiffness matrix of the robot; X represents the deformation of the robot's end effector under the action of external force.

[0050] Expand the static stiffness model of the robot and get:

[0051]

[0052] Where, J ij is the i-th row and j-th column of the Jacobian matrix, F i is the i-th element of the external force vector acting on the robot during machining.

[0053] Through the joint stiffness identification experiment of the robot, the external forces applied to the robot during the processing at different postures and the deformation of the robot end caused by the external forces are measured, and the stiffness matrix of the robot can be obtained.

[0054] Step S15: Establish the robot grinding dynamics differential equation:

[0055]

[0056] Where M D 、C D , K D and F D They respectively represent the mass matrix, damping matrix, stiffness matrix and grinding force vector of the Cartesian space at the end of the composite five-DOF working robot.

[0057] Where,

[0058]

[0059] F(t)=k m hba(t)

[0060] a(t)=a0-[X(t)-X(tT)]

[0061] T=60 / Ω

[0062] Where F(t) is the dynamic grinding force, k m represents the grinding force coefficient of the workpiece, h represents the grinding depth of the workpiece, b represents the contact width of the grinding, a(t) represents the surface vibration of the workpiece, a0 represents the initial surface vibration of the workpiece, X(t) represents the surface amplitude of the workpiece at time t, T is the rotation period of the spindle, and Ω represents the spindle speed.

[0063] In one embodiment, step S2 includes converting the differential equation of grinding dynamics of the composite five-DOF working robot into a time-delay grinding dynamics equation in state space form:

[0064]

[0065] Where,

[0066] and,

[0067] In one embodiment, step S3 includes discretizing the delay time of the time-delay grinding dynamics equation in the state space form, that is, dividing the delay time τ0 into n intervals with a width of Δt, τ0 = nΔt, and t t =iΔt. And by direct integration, the solution of the time-delay grinding dynamics equation in the approximate state space form in the corresponding time interval is obtained:

[0068]

[0069] Let τ=iΔt+Δt-ξ, so the above formula becomes:

[0070]

[0071] In one embodiment, step S4 includes:

[0072] Step S41: The time lag part and the state term in the solution of the time lag grinding dynamics equation are taken as a whole, that is, X(iΔt+Δt-ξ-τ0) and X(iΔt+Δt-ξ) are taken as a whole, and the third-order Hermite interpolation polynomial is used to approximate the term X(kh+h-ξ-T)-X(kh+h-ξ). The third-order Hermite interpolation polynomial is defined as follows:

[0073]

[0074]

[0075]

[0076]

[0077] Where H1 and H2 are adjacent time t i and t i+1 The relative displacement coefficients at the adjacent time t i and t i+1 The relative velocity coefficient at .

[0078] Time interval [t i , t i+1 ] can be expressed as:

[0079]

[0080] Substituting the interpolated result into X(iΔt+Δt), we get:

[0081]

[0082] Expressed in matrix form:

[0083]

[0084] Where, 0 is represented by a 6×6 zero matrix, X i+1-n Expressed as X(iΔt+Δt-τ0), X i+1 Expressed as X(iΔt+Δt), X i-n Expressed as X(iΔt-τ0), X i Expressed as X(iΔt).

[0085] Step S42: Substitute the equation obtained after the above interpolation into X(iΔt+Δt) to obtain the discrete dynamic iteration formula:

[0086] X i+1 =P0X i +R i (a0+X i-n -X i )+R i+1 (a0+X i+1-n -X i+1 )

[0087] Where,

[0088] P0=e AΔt

[0089]

[0090]

[0091]

[0092]

[0093]

[0094]

[0095]

[0096] Where C i，k Indicated as H k The coefficient preceding the i-th power of ξ.

[0097] The above discrete dynamic iteration formula is sorted out as follows:

[0098] X i+1 =S i X i +T i X i-n +L i X i+1-n +a0(T i +L i )

[0099] Where,

[0100] S i =(I+R i+1 ) -1 (P0-R i )

[0101] T i =(I+R i+1 )-1 R i

[0102] L i =(I+R i+1 ) -1 R i+1

[0103] In one embodiment, step S5 includes:

[0104] Step S51: The following matrix sequence is derived through the sorted discrete dynamic iterative formula:

[0105]

[0106] Step S52: Through the matrix sequence derived above, the state transfer matrix can be obtained as follows:

[0107] Φ=D m-1 D m-2 …D1D0

[0108] In one embodiment, step S6 includes solving the modulus of the state transfer matrix according to Floquet theory to determine the vibration stability. When the modulus of the eigenvalue of the state transfer matrix Φ is greater than 1, the composite five-degree-of-freedom working robot grinding processing system is in a vibration state; when the modulus of the eigenvalue of the state transfer matrix Φ is less than 1, the composite five-degree-of-freedom working robot grinding processing system is in a stable state.

[0109] Compared with the prior art, the present invention avoids the prediction error caused by the separate discretization of vibration state variables and has higher prediction accuracy and reliability. BRIEF DESCRIPTION OF THE DRAWINGS

[0110] The present invention will be described in more detail below based on embodiments and with reference to the accompanying drawings.

[0111] Figure 1 represents the DH reference coordinate system of the 3T2R configuration robot;

[0112] FIG2 shows a schematic diagram of the mechanism of the 3T2R configuration robot;

[0113] Figure 3 Indicates the external cylindrical cutting grinding model;

[0114] Figure 4 Represents a flow chart for drawing stability lobe diagrams;

[0115] Figure 5 represents the stability lobe diagram;

[0116] Figure 6 It shows the comparison between the stability lobe diagram drawn by this method and EFDM;

[0117] Figure 7 It is a robot experimental platform; DETAILED DESCRIPTION

[0118] The following describes the implementation of the present invention through specific embodiments and in conjunction with the accompanying drawings.

[0119] The present invention proposes a robot grinding chatter prediction method considering discrete vibration speeds, comprising the following steps:

[0120] Step S1: Establishing the grinding dynamics differential equation of the composite five-degree-of-freedom working robot, including the following steps:

[0121] Step S11: Create Figure 1 The DH reference coordinate system and link parameter table of the composite five-DOF operating robot are shown, and the homogeneous transformation matrix of the adjacent coordinate systems is obtained as follows:

[0122]

[0123]

[0124]

[0125]

[0126]

[0127] Where, is the rotation angle of J1~J5 axis; α is the rotation angle of the turntable; θ i is the joint variable; d i is the joint offset; a i Connecting rod length; X1, X2 and X3 are the travels of the nut on the ball screws S1, S2 and S3, robot structural parameters e and tool compensation L P All are considered as 0.

[0128] Step S12: Based on the above homogeneous transformation matrix, the homogeneous transformation matrix of the coordinate system {O5} of the end effector of the composite five-degree-of-freedom working robot relative to the base coordinate system {O0} can be obtained:

[0129]

[0130] Where,

[0131]

[0132]

[0133]

[0134]

[0135] Step S13: Based on the obtained homogeneous transformation matrix, the Jacobian matrix can be directly obtained by using the vector product method without the need for derivation. The expression for solving the Jacobian matrix using the vector product method is:

[0136]

[0137] Where, J Pi 、J Oi is the i-th column of the Jacobian matrix of the robot's moving joints and rotating joints; P e is the position vector of the robot's end effector relative to the base coordinate system {O0}, which can be represented by the homogeneous transformation matrix The first three elements of the fourth column are obtained; P i-1 It can be represented by the homogeneous transformation matrix The first three elements of the fourth column are obtained; z i-1 Can be obtained by the rotation matrix The third column of is obtained.

[0138] according to Figure 1 It can be seen from the DH reference coordinate system shown in Figure 2 and the schematic diagram of the composite five-degree-of-freedom working robot mechanism shown in Figure 2 that joint 1 and joint 2 work together to make the connecting rod 1 move horizontally along the direction of X0 and rotate around the direction of Z0. Therefore, the translation and rotation are decoupled, and joint 1 is equivalently treated as a moving joint and joint 2 is treated as a rotating joint.

[0139] Solving the J1 to J5 columns of the Jacobian matrix, we can obtain:

[0140] J1=[1 0 0 0 0 0] T

[0141] J2=[-(L3-L3cosα)L3sinα 0 0 0 1] T

[0142] J3=[0 0 1 0 0 0] T

[0143] J4=[0 0 0 0 0 1] T

[0144]

[0145] Therefore, the Jacobian matrix of the composite five-DOF operating robot is:

[0146] J=[J1 J2 J3 J4 J5]

[0147] Step S14: Establish the robot dynamics differential equation:

[0148]

[0149] Where MC and K are the 5×5 mass matrix, damping matrix, and stiffness matrix of the robot, respectively; x is the displacement vector; and F is the external force vector.

[0150] Where,

[0151] and,

[0152] Where c ij represents the element in row i and column j of the damping matrix C; Indicates the movement speed of the robot joint; coefficient c ijk It is called the Christoffel symbol of the first form and is given by:

[0153]

[0154] Where b ij represents the element in the i-th row and j-th column of the mass matrix, q k represents the joint angles of the robot, Represents the partial derivative of the elements in the mass matrix with respect to the joint angles of the robot.

[0155] Calculate the stiffness matrix of the composite five-DOF operating robot, and establish the static stiffness model of the robot based on the Jacobian matrix and Hooke's law:

[0156] F=J -T K θ J -1 X

[0157] Where F represents the external force vector of the robot during machining; K θ represents the stiffness matrix of the robot; X represents the deformation of the robot's end effector under the action of external force.

[0158] Expand the static stiffness model of the robot and get:

[0159]

[0160] Where, J ij is the i-th row and j-th column of the Jacobian matrix, F i is the i-th element of the external force vector acting on the robot during machining.

[0161] Through the joint stiffness identification experiment of the robot, the external forces applied to the robot during the processing at different postures and the deformation of the robot end caused by the external forces are measured, and the stiffness matrix of the robot can be obtained.

[0162] Step S15: Establish the robot grinding dynamics differential equation:

[0163]

[0164] Where M D 、C D , K D and F D They respectively represent the mass matrix, damping matrix, stiffness matrix and grinding force vector of the Cartesian space at the end of the composite five-DOF working robot.

[0165] Where,

[0166]

[0167] according to Figure 3 As shown in the figure, the grinding force during the grinding process is mainly radial force;

[0168] F(t)=k m hba(t)

[0169] a(t)=a0-[X(t)-X(tT)]

[0170] T=60 / Ω

[0171] Where F(t) is the dynamic grinding force, k m represents the grinding force coefficient of the workpiece, h represents the grinding depth of the workpiece, b represents the contact width of the grinding, a(t) represents the surface vibration of the workpiece, a0 represents the initial surface vibration of the workpiece, X(t) represents the surface amplitude of the workpiece at time t, T is the rotation period of the spindle, and Ω represents the spindle speed.

[0172] Step S2: The grinding dynamics differential equation of the composite five-DOF robot is transformed into a time-delay grinding dynamics equation in state space form:

[0173]

[0174] Where,

[0175]

[0176]

[0177] Step S3 includes: discretizing the delay time of the time-delay grinding dynamics equation in the state space form, that is, dividing the delay time τ0 into n intervals with a width of Δt, τ0 = nΔt, and t i =iΔt. And by direct integration, the solution of the time-delay grinding dynamics equation in the approximate state space form in the corresponding time interval is obtained:

[0178]

[0179] Let τ=iΔt+Δt-ξ, so the above formula becomes:

[0180]

[0181] Step S4: constructing a discrete dynamic iterative formula, including the following steps:

[0182] Step S41: The time lag part and the state term in the solution of the time lag grinding dynamics equation are taken as a whole, that is, X(iΔt+Δt-ξ-τ0) and X(iΔt+Δt-ξ) are taken as a whole, and the third-order Hermite interpolation polynomial is used to approximate the term X(kh+h-ξ-T)-X(kh+h-ξ). The third-order Hermite interpolation polynomial is defined as follows:

[0183]

[0184]

[0185]

[0186]

[0187] Where H1 and H2 are adjacent time t i and t i+1 The relative displacement coefficients at the adjacent time t i and t i+1 The relative velocity coefficient at .

[0188] Time interval [t i , t i+1 ] can be expressed as:

[0189]

[0190] Substituting the interpolated result into X(iΔt+Δt), we get:

[0191]

[0192] Write the above formula in matrix form:

[0193]

[0194] Where, 0 is represented by a 6×6 zero matrix, X i+1-n Expressed as X(iΔt+Δt-τ0), X i+1 Expressed as X(iΔt+Δt), X i-n Expressed as X(iΔt-τ0), X i Expressed as X(iΔt).

[0195] Step S42: Substitute the above interpolated equation into X(iΔt+Δt) to obtain:

[0196] X i+1 =P0X i +R i (a0+X i-n -X i )+R i+1 (a0+X i+1-n -X i+1 )

[0197] Where,

[0198] P0=e AΔt

[0199]

[0200]

[0201]

[0202]

[0203]

[0204]

[0205]

[0206] Where C i，k Indicated as H k The coefficient preceding the i-th power of ξ.

[0207] The above discrete dynamic iteration formula is sorted out as follows:

[0208] X i+1 =S i X i +T i X i-n +L i X i+1-n +a0(T i +L i )

[0209] Where,

[0210] S i =(I+R i+1 ) -1 (P0-R i )

[0211] T i =(I+R i+1 ) -1 R j

[0212] L i =(I+R i+1 ) -1 R i+1

[0213] Step S5: constructing a state transfer matrix, including the following steps:

[0214] Step S51: The following matrix sequence is derived through the sorted discrete dynamic iterative formula:

[0215]

[0216] Step S52: Through the matrix sequence derived above, the state transfer matrix can be obtained as follows:

[0217] Φ=D m-1 D m-2 …D1D0

[0218] Step S6: According to Figure 4 The flowchart is used and the stability lobe diagram is drawn on MATLAB according to the Floquet theory to judge the flutter stability.

[0219] Figure 5 The stability lobe diagram drawn for this method shows the chatter region, which is the area above the stability lobe diagram. Within this region, chatter will occur at any grinding depth and any spindle speed. The stable region is the area below the stability lobe diagram. Within this region, chatter will not occur at any grinding depth and any spindle speed. Therefore, in actual grinding, it is best to select machining parameters within the stable region to avoid chatter.

[0220] Figure 6The stability lobe diagrams generated by this method were compared with those generated by the EFDM method. Because the EFDM chatter prediction method does not consider the influence of discrete vibration velocities and only considers discrete vibration displacements, it results in abrupt changes in the critical grinding depth corresponding to the spindle speed range of 1800 to 2400 r / min, thereby reducing its reliability and accuracy. The overall discretization method proposed in this paper, which considers discrete vibration velocities, does not exhibit abrupt changes in the critical grinding depth corresponding to the spindle speed range of 1800 to 2400 r / min. Therefore, it has higher reliability and accuracy than the EFDM chatter prediction method.

[0221] Figure 7 A robot experimental platform was used to carry out grinding experiments. The collected acceleration signals were subjected to fast Fourier transform to determine whether chatter occurred, which verified the effectiveness of the method.

[0222] These embodiments are merely exemplary and serve only as an illustrative example. On this basis, various replacements and improvements may be made to the present invention, all of which fall within the scope of protection of the present invention.

Claims

1. A robot grinding chatter prediction method considering discrete vibration speeds, characterized in that: The following steps are involved: Step S1: establishing the grinding dynamics differential equation of the composite five-degree-of-freedom working robot; Step S2: transforming the grinding dynamics differential equation of the composite five-DOF working robot into a time-delay grinding dynamics equation in state space form; Step S3: For the time-delay grinding dynamics equation in the state space form, first discretize the delay time, and then obtain the approximate solution of the equation in the corresponding time interval by direct integration; Step S4: The time lag part and the state term in the solution of the time lag grinding dynamics equation are treated as a whole, and processed using third-order Hermite interpolation to obtain a discrete dynamic iterative formula; Step S5: constructing a state transfer matrix through a discrete dynamic iterative formula; Step S6: Determine the flutter stability using Floquet theory.

2. The robot grinding chatter prediction method considering discrete vibration speeds according to claim 1, characterized in that: The step S1 comprises: Step S11: Establish the DH reference coordinate system and link parameter table of the composite five-DOF working robot, and obtain the following homogeneous transformation matrix of adjacent coordinate systems: Where, is the rotation angle of J1~J5 axis; ɑ is the rotation angle of the turntable; θ i is the joint variable; d i is the joint offset; a i Connecting rod length; X1, X2 and X3 are the travels of the nut on the ball screws S1, S2 and S3, robot structural parameters e and tool compensation L P All are considered as 0; Step S12: Based on the above homogeneous transformation matrix, the homogeneous transformation matrix of the coordinate system {O5} of the end effector of the composite five-degree-of-freedom working robot relative to the base coordinate system {O0} can be obtained: Where, Step S13: Based on the obtained homogeneous transformation matrix, the Jacobian matrix is ​​directly obtained by using the vector product method without the need for derivation. The expression for solving the Jacobian matrix by the vector product method is: Where, J Pi 、J Oi is the i-th column of the Jacobian matrix of the robot's moving joints and rotating joints; P e is the position vector of the robot's end effector relative to the base coordinate system {O0}, which can be represented by the homogeneous transformation matrix The first three elements of the fourth column are obtained; P i-1 It can be represented by the homogeneous transformation matrix The first three elements of the fourth column are obtained; z i-1 Can be obtained by the rotation matrix The third column of is obtained; Solving the J1 to J5 columns of the Jacobian matrix, we can obtain: J1＝[1 0 0 0 0 0] T <h2 style=";text-align:left;direction:ltr">J2=[-(L3-L3cosα) L3sinα 0 0 0 1]<h2 style=";text-align:left;direction:ltr"> T J3＝[0 0 1 0 0 0] T J4＝[0 0 0 0 0 1] T Therefore, the Jacobian matrix of the composite five-DOF operating robot is: J=[J1 J2 J3 J4 J5] Step S14: Establish the robot dynamics differential equation: Where M, C, and K are the 5×5 mass matrix, damping matrix, and stiffness matrix of the robot, respectively; x is the displacement vector; and F is the external force vector. Where c ij represents the element in row i and column j of the damping matrix C; Indicates the movement speed of the robot joint; coefficient c ijk It is called the Christoffel symbol of the first form and is given by: Where b ij represents the element in the i-th row and j-th column of the mass matrix, q k represents the joint angles of the robot, Represents the partial derivative of the elements in the mass matrix with respect to the joint angles of the robot; Calculate the stiffness matrix of the composite five-DOF operating robot, and establish the static stiffness model of the robot based on the Jacobian matrix and Hooke's law: F=J -T K θ J -1 X Where F represents the external force vector of the robot during machining, K θ represents the stiffness matrix of the robot, and X represents the deformation of the robot's end effector under the action of external force. Expanding the above formula, we get: Where, J ij is the i-th row and j-th column of the Jacobian matrix, F i is the i-th element of the external force vector acting on the robot during machining; Through the robot joint stiffness identification experiment, the external forces on the robot during the processing in different postures and the deformation of the robot end caused by the external forces are measured, so as to obtain the stiffness matrix of the robot; Step S15: Establish the robot grinding dynamics differential equation: Where M D 、C D , K D and F D They represent the mass matrix, damping matrix, stiffness matrix, and grinding force vector of the composite five-DOF working robot end in Cartesian space respectively; F(t)=k m hba(t) a(t)=a0-[X(t)-X(tT)] T=60 / Ω Where F(t) is the dynamic grinding force, k m represents the grinding force coefficient of the workpiece, h represents the grinding depth of the workpiece, b represents the contact width of the grinding, a(t) represents the surface vibration of the workpiece, a0 represents the initial surface vibration of the workpiece, X(t) represents the surface amplitude of the workpiece at time t, T is the rotation period of the spindle, and Ω represents the spindle speed.

3. The robot grinding chatter prediction method considering discrete vibration speeds according to claim 1, characterized in that: The step S2 includes converting the differential equation of grinding dynamics of the composite five-DOF working robot into a time-delay grinding dynamics equation in state space form: Where, 4. The robot grinding chatter prediction method considering discrete vibration speeds according to claim 1, characterized in that: The step S3 includes: discretizing the delay time of the time-delay grinding dynamics equation in the state space form, that is, dividing the delay time τ0 into n intervals with a width of Δt, τ0 = nΔt, and t i =iΔt, and the solution of the time-delay grinding dynamics equation in the approximate state space form in the corresponding time interval is obtained by direct integration: Let τ=iΔt+Δt-ξ, so the above formula becomes:

5. The robot grinding chatter prediction method considering discrete vibration speeds according to claim 1, characterized in that: The step S4 comprises: Step S41: The time lag part and the state term in the solution of the time lag grinding dynamics equation are taken as a whole, that is, X(iΔt+Δt-ξ-τ0) and X(iΔt+Δt-ξ) are taken as a whole, and the third-order Hermite interpolation polynomial is used to approximate the term X(kh+h-ξ-T)-X(kh+h-ξ). The third-order Hermite interpolation polynomial is defined as follows: Where H1 and H2 are adjacent time t i and t i+1 The relative displacement coefficients at the adjacent time t i and t i+1 Relative velocity coefficient at ; Time interval [t i , t i+1 ] can be expressed as: Substituting the interpolated result into X(iΔt+Δt), we get: B[X(iΔt+Δt-ξ-τ0)-X(iΔt+Δt-ξ)] =H1(ξ)B[X(iΔt+Δt-τ0)-X(iΔt+Δt)] +H2(ξ)B[X(iΔt-τ0)-X(iΔt)] +H3(ξ)DI(iΔt+Δt)+H4(ξ)DI(iΔt) so, Where, 0 is represented by a 6×6 zero matrix, X i+1-n Expressed as X(iΔt+Δt-τ0), X i+1 Expressed as X(iΔt+Δt), X i-n Expressed as X(iΔt-τ0), X i Expressed as X(iΔt); Step S42: Substitute the above interpolated equation into X(iΔt+Δt) to obtain: X i+1 =P0X i +R i (a0+X i-n -X i )+R i+1 (a0+X i+1-n -X i+1 ) Where, P0=e AΔt Where C i，k Indicated as H k The coefficient preceding the i-th power of ξ in the corresponding column; The above discrete dynamic iteration formula is sorted out as follows: X i+1 =S i X i +T i X i-n +L i X i+1-n +a0(T i +L i ); Where, S i =(I+R i+1 ) -1 (P0-R i ) T i =(I+R i+1 ) -1 R i L i =(I+R i+1 ) -1 R i+1 6. The robot grinding chatter prediction method considering discrete vibration speeds according to claim 1, characterized in that: The step S5 comprises: Step S51: The following matrix sequence is derived through the sorted discrete dynamic iterative formula: Step S52: Through the matrix sequence derived above, the state transfer matrix can be obtained as follows: Φ=D m-1 D m-2 …D1D0 7. The robot grinding chatter prediction method considering discrete vibration speeds according to claim 1, characterized in that: The step S6 includes solving the modulus of the state transfer matrix according to the Floquet theory to determine the vibration stability. When the modulus of the eigenvalue of the state transfer matrix Φ is greater than 1, the composite five-degree-of-freedom operation robot grinding processing system is in a vibration state; when the modulus of the eigenvalue of the state transfer matrix Φ is less than 1, the composite five-degree-of-freedom operation robot grinding processing system is in a stable state.

Images