A robot milling posture optimization method based on principal stiffness difference
By discretizing the milling path and utilizing the redundant angle characteristics to optimize the robot milling posture, the main stiffness difference and vibration stability index are calculated, and the problem of modal coupling vibration in robot milling processing is solved, achieving improved stability and stiffness of the processing process.
Patent Information
- Application Number
- CN202411851195.0
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2024-12-16
- Publication Date
- 2025-09-23
- Estimated Expiration
- 2044-12-16
AI Technical Summary
In existing robotic milling techniques, modal coupling chatter is suppressed by changing the workpiece clamping position or adjusting the feed direction, which limits machining flexibility and fails to fully utilize the advantages of robotic machining.
The milling machining path is discretized into multiple trajectory points. The redundant angle characteristics of the robot are utilized to optimize the robot milling posture by calculating the principal stiffness difference (PSD) and the chatter stability index (CSI) to maximize the overall stiffness performance and avoid modal coupling chatter.
The robot achieves improved stability during milling without changing the feed direction or adjusting the workpiece clamping, giving full play to the robot's flexibility, improving machining stability and rigidity, and avoiding modal coupling chatter.
Smart Images

Figure CN119458406B_ABST
Abstract
Description
Technical Field
[0001] The invention belongs to the field of robot milling processing, and in particular relates to a robot milling posture optimization method based on main stiffness difference. Background Art
[0002] Industrial robots offer advantages such as a large workspace, excellent flexibility, high production efficiency, and excellent safety performance. They can also be combined with mobile guides for expanded applicability and offer cost advantages over large, specialized machine tools. However, due to their inherent serial structure, industrial robots typically exhibit weak rigidity. Generally speaking, the rigidity of a 6-DOF industrial robot is only 1 / 50 that of a traditional CNC machine tool. When a robot performs metal milling, the large and constantly changing milling forces acting on the tool tip at the end of the robot can cause modal coupling chatter. Modal coupling chatter causes the robot to experience violent low-frequency vibrations, forming an elliptical vibration trajectory with continuously increasing major and minor axes at the tool tip, severely affecting the machining accuracy of the workpiece and even damaging the robot's main structure.
[0003] A common approach to addressing modal-coupled chatter in robotic milling is to establish a two-degree-of-freedom kinematic model for the tool-workpiece system at the end. Based on the resulting dynamic equations, the mechanism of modal-coupled chatter is analyzed. Based on the derived stability criteria for modal-coupled chatter in robotic milling, modal-coupled chatter in robotic milling is suppressed by modifying the workpiece clamping, optimizing milling process parameters, and optimizing the robot's spatial posture during milling.
[0004] Currently, published methods for optimizing the workpiece clamping position in robotic milling improve machining stability by changing the workpiece clamping position based on a calculated chatter stability lobe diagram. Published methods for suppressing modal coupling chatter based on the robot's stiffness characteristics suppress modal coupling chatter during robotic milling by adjusting the feed direction based on the stability criterion of the milling dynamics equation. However, these two methods, by employing methods such as changing the workpiece clamping position and feed direction, impose certain limitations on robotic machining and fail to fully utilize the flexibility of robotic machining. Summary of the Invention
[0005] Technical issues to be solved:
[0006] To overcome the shortcomings of existing technologies, the present invention provides a robot milling posture optimization method based on principal stiffness differences. This method discretizes the milling path into a number of machining trajectory points. At each trajectory point, the robot's redundant angle characteristics are utilized to optimize the robot's posture. The optimization goal is to maximize the sum of principal stiffness differences at each point along the entire machining trajectory. This allows the robot to perform milling with optimal stiffness performance. This ensures a smooth milling process and avoids modal coupling chatter.
[0007] The technical solution of the present invention is: a robot milling posture optimization method based on main stiffness difference, the specific steps are as follows:
[0008] Divide the robot milling machining trajectory into equal distances and number each trajectory point;
[0009] At each machining trajectory point, the tool tip stiffness ellipsoid of the robot is calculated under the posture corresponding to the specific redundant angle;
[0010] Project the tool tip point stiffness ellipsoid onto the milling plane to obtain the robot's principal stiffness difference PSD;
[0011] Taking the maximum sum of the principal stiffness differences (PSDs) corresponding to all machining trajectory points as the optimization goal, under the constraints of physical interference, flexible space, and joint speed, the optimal redundant angles of all trajectory points are traversed to solve the optimal robot milling posture.
[0012] A further technical solution of the present invention is: the milling trajectory is obtained according to the actual milling task; the trajectory points are numbered as [P1, P2, ···, P i ,···,P N-1 ,P N ] T .
[0013] A further technical solution of the present invention is: the process of solving the tool tip stiffness ellipsoid of the robot in the posture corresponding to the specific redundant angle is as follows:
[0014] According to the processing trajectory points, the corresponding redundant angle is determined;
[0015] Determine the robot's spatial posture through coordinate transformation;
[0016] According to the inverse kinematics of the robot, the joint angle configuration of the redundant angle corresponding to the machining trajectory point is solved;
[0017] According to the joint angle configuration, solving the end tool tip point stiffness matrix of the robot in the corresponding posture; the end tool tip point stiffness matrix includes a force / translational displacement stiffness submatrix, a force / rotational displacement stiffness submatrix, a moment / translational displacement stiffness submatrix, and a moment / rotational displacement stiffness submatrix;
[0018] According to the eigenvalues and eigenvectors of the force / translation displacement submatrix, the stiffness ellipsoid of the tool tip of the robot at the corresponding posture is solved.
[0019] A further technical solution of the present invention is: the coordinate conversion process is to define the machining coordinate system in the robot milling process as the coordinate system {M}, where x M Axis along the tool feed direction, z M Axis along the direction of tool axis vector; the machining coordinate system {M} is rotated around z M The new coordinate system obtained after the axis rotates the redundant angle is defined as the tool coordinate system {C}. The tool coordinate system {C} is then used to reflect the orientation information of the tool axis vector and determine the spatial posture of the robot. The homogeneous transformation matrix of the tool coordinate system {C} relative to the robot base coordinate system {B} is recorded as
[0020] A further technical solution of the present invention is: the calculation formula for the joint angle configuration of the redundant angle corresponding to the processing trajectory point is:
[0021] ξ i =f IK (δ R,i )
[0022] Among them, f IK (δ R,i ) is based on the trajectory point P i The redundancy angle δ R,i Joint angle configuration ξ i the generalized function to be solved;
[0023] According to the joint angle configuration, the stiffness matrix of the robot's end tool tip in the corresponding posture is solved, and the formula is as follows:
[0024]
[0025] Among them, K c K is the tool tip stiffness matrix of the robot; FD is the force / translation displacement stiffness submatrix; K Fθ is the force / rotational displacement stiffness submatrix; K MD is the moment / translation displacement stiffness submatrix; K Mθ is the moment / rotational displacement stiffness submatrix; in the subscripts of the internal submatrix, F and M represent the force and moment applied to the tool tip, and D and θ represent the translational displacement and rotational displacement generated by the tool tip; J c is the Jacobian matrix at the robot tool tip, K ξ is the joint stiffness matrix of the robot;
[0026] According to the force / translation displacement submatrix K FDThe eigenvalue and eigenvector of the robot are used to obtain the stiffness ellipsoid equation of the tool tip at the end of the robot under the corresponding posture. The calculation formula is:
[0027]
[0028] Among them, σ1, σ2, σ3 are the matrix K FD The three eigenvalues sorted from large to small correspond to the lengths of the three semi-axes of the tool tip stiffness ellipsoid; the corresponding unit eigenvectors are r1, r2, and r3, corresponding to the directions of the three semi-axes of the tool tip stiffness ellipsoid; define the ellipsoid coordinate system {E}, whose x E 、y E 、z E The axis direction is the direction of the three semi-axes, and the ellipsoid equation has a standard form in the ellipsoidal coordinate system {E}.
[0029] A further technical solution of the present invention is: the process of projecting the tool tip point stiffness ellipsoid onto the milling processing plane to obtain the corresponding ellipse intersection equation is as follows:
[0030] The rotation matrix of the ellipsoidal coordinate system {E} relative to the base coordinate system {B} Denoted as:
[0031] The rotation matrix of the machining coordinate system {M} relative to the base coordinate system {B} Expressed as:
[0032] The relationship between the ellipsoid coordinate system {E} and the machining coordinate system {M} is obtained as follows:
[0033] according to The tool tip stiffness ellipsoid equation is converted from the ellipsoid coordinate system {E} to the machining coordinate system {M}. The calculation formula is: in, E p=[x E ,y E ,z E ] T 、 M p=[x M ,y M ,z M ] T are the coordinate vectors in the ellipsoid coordinate system {E} and the machining coordinate system {M} respectively;
[0034]
[0035] Among them, k ij is the rotation matrix Each element in , i = 1, 2, 3; j = 1, 2, 3);
[0036] The relationship between the coordinate vectors in the ellipsoid coordinate system {E} and the machining coordinate system {M} is expressed as:
[0037]
[0038] The standard equation of the tool tip stiffness ellipsoid is converted from the ellipsoid coordinate system {E} to the machining coordinate system {M} for expression:
[0039]
[0040] Let z M =0, the elliptical intersection equation of the rigidity ellipsoid in the milling plane is obtained. The specific formula is:
[0041]
[0042] A further technical solution of the present invention is: the calculation process of the robot's main stiffness difference PSD is as follows:
[0043] By combining like terms and making equivalent substitutions, the ellipse intersection equation can be written in matrix form as follows:
[0044]
[0045] make:
[0046] Transform
[0047] For the real symmetric matrix Perform feature decomposition and the calculation formula is:
[0048]
[0049] Among them, λ1 and λ2 are matrices The corresponding eigenvalue, K is the matrix composed of its unit eigenvector;
[0050] Through the matrix K, the coordinate vector [x M ,y M ] T By performing a rotation transformation, the ellipse intersection equation can be converted from the machining coordinate system {M} to the principal stiffness coordinate system {P}; the specific formula is as follows:
[0051] P p=K M p
[0052]
[0053] After expansion:
[0054]
[0055] in, M p=[x M ,y M ] T 、 P p=[x P ,y P ] T are the coordinate vectors in the machining coordinate system {M} and the principal stiffness coordinate system {P} respectively;
[0056] In the principal stiffness coordinate system {P}, the intersection of the tool tip stiffness ellipsoid and the machining plane is written as a standard form of the ellipse equation. The principal stiffness difference PSD of the robot on the machining plane is quantified by the difference between the major and minor axes. The formula is as follows:
[0057]
[0058] in,
[0059] A further technical solution of the present invention is: the sum of the principal stiffness differences PSD corresponding to all the processing trajectory points is defined as a flutter stability index CSI, and the maximum flutter stability index CSI is used as an optimization target.
[0060] A further technical solution of the present invention is: the optimization model of the optimal robot milling posture is:
[0061]
[0062] Among them, [-90°, 90°] is the limit range of the robot's rotation around the tool axis vector under the interference of physical space after the milling tool is installed at the end of the robot; ξ Flexible It is a flexible workspace where the robot is far away from the singularity; is the rated load speed of the jth joint; μ is the set processing safety factor; ΔT is the processing time between adjacent milling trajectory points; ΔS is the interval distance between adjacent milling trajectory points; f is the set milling feed speed.
[0063] A robot milling posture optimization system based on principal stiffness difference comprises at least one processor and a memory communicatively connected to the at least one processor; wherein the memory stores a computer program executable by the at least one processor, and the computer program is executed by the at least one processor so that the at least one processor can execute the robot milling posture optimization method based on principal stiffness difference.
[0064] Beneficial effects
[0065] The beneficial effects of the present invention are as follows: This invention addresses the problem of modal coupling chatter in milling, defines the robot's primary stiffness difference index (PSD) and chatter stability index (CSI) for the milling process, and further proposes a milling posture optimization algorithm based on the robot's redundant angle characteristics, thereby improving the machining process stiffness and achieving chatter avoidance during the robot milling process. Specific advantages are as follows:
[0066] 1. This invention uses an ellipsoid to reflect the robot's tool tip's ability to resist deformation in various directions under a specific spatial posture. The points on the ellipsoid's surface constitute a set of forces required to produce unit deformation of the robot's tool tip in each direction. A greater vector length from the origin to the ellipsoid's surface indicates a greater force required to produce unit deformation of the tool tip in that direction, i.e., a greater stiffness of the robot's tool tip in that direction. Conversely, a smaller force required to produce unit deformation of the tool tip in that direction indicates a lower stiffness of the robot's tool tip in that direction.
[0067] 2. The present invention does not need to limit the milling feed direction of the robot, nor does it need to adjust the workpiece clamping or correct the process parameters, thus giving full play to the flexibility advantage of the robot compared to the machine tool.
[0068] 3. This invention discretely divides the machining trajectory and sets constraints on the corresponding joint angle changes between adjacent machining trajectory points. This approach can be applied to robotic milling of large structural parts. While maintaining overall smooth motion, it improves the horizontal stiffness of the robot's end-tool tip, achieving enhanced milling stability during continuous machining. BRIEF DESCRIPTION OF THE DRAWINGS
[0069] Figure 1 This is an overall flow chart of a robot milling posture optimization method based on principal stiffness difference of the present invention;
[0070] Figure 2 Schematic diagram of milling machining trajectory points in an embodiment of the present invention;
[0071] Figure 3 Schematic diagram of the tool tip point stiffness ellipsoid in an embodiment of the present invention;
[0072] Figure 4 1 is a schematic diagram of the projection of the tool tip point stiffness ellipsoid processing plane in an embodiment of the present invention;
[0073] Figure 5 Schematic diagram of the main stiffness of the robot processing plane in an embodiment of the present invention;
[0074] Figure 6 This is a flow chart of the robot milling posture optimization algorithm in an embodiment of the present invention;
[0075] Figure 7It is the two-degree-of-freedom machining dynamics model of "tool-workpiece" in robot milling in an embodiment of the present invention;
[0076] Figure 8 This is a diagram of the experimental environment for verifying the robot milling posture optimization based on the main stiffness difference implemented in the present invention;
[0077] Figure 9 is the principal stiffness difference PSD at the trajectory point in the embodiment of the present invention i and the redundancy angle δ R,i Correspondence diagram;
[0078] Figure 10 1 is a time-domain and frequency-domain signal diagram of milling acceleration when milling posture optimization is not performed in an embodiment of the present invention;
[0079] Figure 11 It is the time and frequency domain signal diagram of milling acceleration after milling posture optimization;
[0080] Explanation of the accompanying symbols: 1. Acceleration sensor, 2. Electric spindle, 3. Milling cutter, 4. Vise, 5. Workpiece, 6. Industrial robot, 7. Data acquisition equipment, 8. Computer. DETAILED DESCRIPTION
[0081] The embodiments described below with reference to the accompanying drawings are exemplary and are intended to explain the present invention, but should not be construed as limiting the present invention.
[0082] Based on the problems that the existing technology adopts methods such as changing the workpiece clamping and changing the feed direction to suppress the modal coupling chatter of robot milling, which imposes certain restrictions on robot processing and fails to fully utilize the flexibility advantages of robot processing, the present invention provides a robot milling posture optimization method based on the main stiffness difference, and the specific steps are as follows:
[0083] Step 1: According to the actual milling task, obtain its machining trajectory. Divide it into equidistant points to obtain several machining trajectory points, which are numbered in sequence and recorded as [P1, P2, ···, P i ,···,P N-1 ,P N ] T , which is convenient for subsequent optimization calculations.
[0084] Step 2: Based on the robot's redundant angle characteristics, calculate the tool tip stiffness ellipsoid at the end of the robot at a specific machining trajectory point and in the spatial posture corresponding to the specific redundant angle. The specific solution process is:
[0085] For the i-th processing trajectory point P i , and the corresponding redundancy angle is δ R,i The coordinate system {M} is defined as the machining coordinate system in the robot milling process, where xM Axis along the tool feed direction, z M The axis is along the direction of the tool axis vector. M Axis rotation redundancy angle δ R,i The new coordinate system obtained is defined as the tool coordinate system {C}. The tool coordinate system {C} can be used to reflect the orientation information of the tool axis vector and determine the spatial posture of the robot. The homogeneous transformation matrix of the tool coordinate system {C} relative to the base coordinate system {B} is According to the robot inverse kinematics, the corresponding joint angle configuration ξ can be solved i , the calculation formula is:
[0086] ξ i =f IK (δ R,i )
[0087] Among them, f IK (δ R,i ) is based on the trajectory point P i The redundancy angle δ R,i Joint angle configuration ξ i The generalized function to be solved.
[0088] According to the joint angle combination ξ obtained by the solution i , we can further obtain the stiffness matrix of the robot's end tool tip in the corresponding posture, and the calculation formula is:
[0089]
[0090] Among them, K c K is the tool tip stiffness matrix of the robot. FD is the force / translation displacement stiffness submatrix; K Fθ is the force / rotational displacement stiffness submatrix; K MD is the moment / translation displacement stiffness submatrix; K Mθ is the moment / rotational displacement stiffness submatrix; in the internal submatrix subscripts, F and M represent the force and moment applied to the tool tip, and D and θ represent the translational displacement and rotational displacement generated by the tool tip. c is the Jacobian matrix at the robot tool tip, K ξ is the joint stiffness matrix of the robot.
[0091] According to the force / translation displacement submatrix K FD The eigenvalue and eigenvector of can be used to obtain the stiffness ellipsoid equation of the robot tool tip in the corresponding posture. The calculation formula is:
[0092]
[0093] Among them, σ1, σ2, σ3 are the matrix KFD The three eigenvalues sorted from largest to smallest correspond to the lengths of the three semi-axes of the tool tip stiffness ellipsoid. The corresponding unit eigenvectors are r1, r2, and r3, corresponding to the directions of the three semi-axes of the tool tip stiffness ellipsoid. Define the ellipsoid coordinate system {E}, where x E 、y E 、z E The axis directions correspond to the directions of the three semi-axes. The ellipsoid equation has a standard form in the ellipsoidal coordinate system {E}. The points on the ellipsoidal surface constitute the set of forces required to cause unit deformation of the robot's tool tip in each direction. A larger vector length from the origin to the ellipsoidal surface indicates a greater force required to cause unit deformation of the tool tip in that direction, i.e., a greater stiffness of the robot's tool tip in that direction. Conversely, a smaller force required to cause unit deformation of the tool tip in that direction indicates a lower stiffness of the robot's tool tip in that direction. Clearly, this ellipsoid reflects the robot's ability to resist deformation in all directions at its tool tip in a specific spatial posture.
[0094] Step 3: Project the tool tip point stiffness ellipsoid onto the robot milling processing plane to obtain the corresponding ellipse intersection equation. According to its major and minor axis characteristics, calculate the corresponding principal stiffness difference index PSD. i The specific calculation process is as follows:
[0095] The rotation matrix of the ellipsoidal coordinate system {E} relative to the base coordinate system {B} It can be written as:
[0096]
[0097] From this, the relationship between the ellipsoid coordinate system {E} and the machining coordinate system {M} can be obtained, which can be expressed as:
[0098]
[0099] according to The ellipsoid equation can be converted from the ellipsoid coordinate system {E} to the machining coordinate system {M} using the following formula:
[0100]
[0101]
[0102] in, E p=[x E ,y E ,z E ] T 、 M p=[x M ,y M ,z M ] Tare the coordinate vectors in the ellipsoid coordinate system {E} and the machining coordinate system {M} respectively.
[0103] Let z M = 0, we can get the equation of the elliptical intersection line of the rigidity ellipsoid in the milling plane. The specific formula is:
[0104]
[0105] By combining similar terms and making equivalent substitutions, the ellipse intersection equation can be written in matrix form. The calculation formula is:
[0106]
[0107] For the real symmetric matrix Perform feature decomposition and the calculation formula is:
[0108]
[0109] Among them, λ1 and λ2 are matrices The corresponding eigenvalues, K is the matrix composed of its unit eigenvectors.
[0110] The ellipse intersection equation can be converted from the machining coordinate system {M} to the principal stiffness coordinate system {P} through the matrix K. The specific formula is as follows:
[0111] P p=L M p
[0112]
[0113] in, M p=[c M ,y M ] T 、 P p=[x P ,y P ] T are the coordinate vectors in the ellipsoid coordinate system {E} and the machining coordinate system {M} respectively.
[0114] According to the characteristics of the robot's serial structure, the end effector has a relatively long overhang distance, so its rigidity is relatively poor. Therefore, modal coupling chatter is prone to occur during the milling process, and it occurs near the robot's lowest natural frequency. In order to mathematically analyze the robot's chatter stability problem, the present invention establishes the following Figure 7 The two-degree-of-freedom milling dynamics model of "tool-workpiece" is shown. In the machining coordinate system {M}, the dynamic equations of the "tool-workpiece" milling system can be established:
[0115]
[0116] Where: M, C, K are the mass, damping, and stiffness matrices of the robot tool tip; Δ M The dynamic displacement vector of the tool tip in the machining coordinate system {M} is denoted as Δ M =[Δx M ,Δy M ] T ; ΔF is the dynamic milling force vector acting on the tool tip of the robot.
[0117] To simplify the analysis process, the following assumptions are made for the milling system: the structural damping of the system will always improve the stability of the milling process. This embodiment analyzes under the idealized condition of ignoring the structural damping. As for the magnitude of the dynamic milling force ΔF acting on the tool tip, this invention only considers its relationship with the tool tip at y. M Dynamic displacement Δy in the direction M There is a direct proportional relationship between them, that is:
[0118] |ΔF|=K cut Δy M
[0119] Among them: K cut The rigidity of the cutting process depends on factors such as the workpiece and tool material, spindle speed, feed rate, and axial depth of cut.
[0120] In this embodiment, the mass matrix M is simplified into a diagonal matrix with equal diagonal elements.
[0121] Based on the above assumptions, the kinetic equation can be rewritten as:
[0122]
[0123] Among them: K cut_M is the cutting process stiffness matrix in the machining coordinate system {M}.
[0124] Obviously K cut_M The form is:
[0125]
[0126] Therefore, the dynamic equation can be expanded as:
[0127]
[0128] It can be seen that the stiffness matrix K is M and y MThere is coupling in the direction. To further simplify the subsequent flutter stability analysis of the system, it is necessary to decouple the stiffness matrix K. Since the stiffness matrix K is a 2×2 real symmetric matrix, it can be similarly diagonalized using an orthogonal matrix. Reflecting it in two-dimensional space is to perform a rotation transformation on it.
[0129] like Figure 7 As shown, the regulation x P Axis along the tool tip point in the machining plane x M O M y M The direction of maximum internal stiffness, y P The axis is along the direction where the tool tip has the smallest stiffness in the machining plane, z P Axis and machining coordinate system {M} z M The axes coincide. Note x P Axis and x M The angle between the axes is β, and the milling force F is related to x P The angle between the axes is γ. Clearly, the principal stiffness coordinate system {P} can be viewed as a rotational transformation of the machining coordinate system {M}. The robot's stiffness matrix K is a diagonal matrix within the principal stiffness coordinate system {P}. Therefore, by transforming the dynamic equations from the machining coordinate system {M} to the principal stiffness coordinate system {P} for analysis, the stiffness matrix K can be decoupled.
[0130] Assume that the dynamic displacement vector of the tool tip in the machining coordinate system {M} and the main stiffness coordinate system {P} satisfies the conversion relationship:
[0131] Δ M =VΔ P
[0132] Where: Δ P The dynamic displacement vector of the tool tip in the main stiffness coordinate system {P} is denoted as Δ P =[Δx P ,Δy P ] T ;
[0133] V is the rotation transformation matrix between the machining coordinate system {M} and the main stiffness coordinate system {P}, according to Figure 8 The spatial geometric relationship shown in the figure, the rotation transformation matrix V can be written as:
[0134]
[0135] Simplifying the dynamic equations, we can get:
[0136]
[0137] Among them: K cut_P—The cutting process stiffness matrix in the main stiffness coordinate system {P} is obviously in the form of:
[0138]
[0139] Similarly, the dynamic equation can be expanded as:
[0140]
[0141] Where: k Pmax is the maximum principal stiffness of the robot in the processing plane; k Pmin is the minimum principal stiffness of the robot in the machining plane.
[0142] Thus, the milling dynamic equation in the main stiffness coordinate system {P} is obtained:
[0143]
[0144] make:
[0145]
[0146] Therefore, the milling dynamics equation in the main stiffness coordinate system {P} can be written as
[0147]
[0148] The dynamic displacement of the tool tip in the main stiffness coordinate system {P} is solved as Δ P Denoted as:
[0149]
[0150] Substitute it into In the above example, we can get:
[0151]
[0152] Written in the form of a system of equations:
[0153]
[0154] For the convenience of calculation, let:
[0155]
[0156] The system of equations can be written as:
[0157]
[0158] Right now:
[0159] (H-λ 2 I) A=0
[0160] Among them: A=[A x ,A y ] T —Dynamic displacement solution of tool tip Δ P The amplitude vector of . Obviously, the condition for the matrix equation shown to have a non-zero solution is: |H-λ 2 I|=0
[0161] Right now:
[0162]
[0163] Expanding it gives the characteristic equation of the milling dynamics model:
[0164]
[0165] Obviously, its characteristic roots satisfy the equation:
[0166]
[0167] λ 2 There are two possible values (λ 2 )1 and (λ 2 )2, we can further obtain the four characteristic roots λ 11 ,λ 12 ,λ 21 ,λ 22 .
[0168] Generally speaking, the robot structural stiffness k Pmax 、k Pmin Much larger than the cutting process stiffness K cut . k′ can be Pmax , k′ Pmin As a positive real number, then λ 2 It must be a complex number with a negative real part. Therefore, the values of the four characteristic roots depend on the signs inside the radical term, which can be discussed in the following two cases:
[0169] when When , the milling process will be in a stable state. When |k′ Pmax -k′ Pmin The larger the | is, the more likely it is that the value inside the square root term will be greater than 0.
[0170] |k′ Pmax -k′ Pmin |=|k Pmax -k Pmin +K cut sin(γ-β)|≈|k Pmax -k Pmin |
[0171] The difference in the robot's main stiffness on the machining plane |k Pmax -k Pmin When | is larger, the value inside the square root term is more likely to be non-negative, and the probability of modal coupling chatter in the robot milling system is lower. In other words, the difference in the main stiffness of the robot's machining surface |k Pmax -k Pmin | is positively correlated with the flutter stability. This can be achieved by increasing the main stiffness difference of the robot processing plane |k Pmax -k Pmin | method to reduce the probability of modal coupling flutter.
[0172] Obviously, in the principal stiffness coordinate system {P}, the intersection of the tool tip stiffness ellipsoid and the machining plane can be written as a standard form of the ellipse equation. The difference between its major and minor axes can be used to quantify the principal stiffness difference PSD of the robot on the machining plane. The specific formula is:
[0173]
[0174] Step 4: Define the principal stiffness difference PSD at each machining trajectory point i The sum is the chatter stability index CSI. Taking the maximum chatter stability index CSI as the optimization target, under the constraints of physical interference, flexible space, joint speed and other conditions, the optimal redundancy angle at each point is traversed and calculated, and then the optimal robot milling processing posture is obtained. The specific optimization model is as follows:
[0175]
[0176] Among them, [-90°, 90°] is the limit range of the robot's rotation around the tool axis vector under the interference of physical space after the milling tool is installed at the end of the robot; ξ Flexible It is a flexible workspace where the robot is far away from the singularity; is the rated load speed of the jth joint; μ is the set processing safety factor; ΔT is the processing time between adjacent milling trajectory points; ΔS is the interval distance between adjacent milling trajectory points; f is the set milling feed speed.
[0177] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention and are not intended to limit the present invention.
[0178] First, according to the milling task, the machining trajectory is divided into equal distances, such as Figure 2 As shown. A milling trajectory is discretized into several machining trajectory points [P1, P2, ···, P i ,···,PN-1 ,P N ] T , the spacing between adjacent processing trajectory points is equal.
[0179] Since the tool axis vector of milling needs to be determined by 5 parameters, and the industrial robot has 6 degrees of freedom, it will produce 1 redundant degree of freedom, which is the redundant angle characteristic of the industrial robot. Therefore, for the i-th processing trajectory point P i , the robot's end effector can rotate around the tool axis vector by a redundant angle δ R,i , achieving the same processing task with different spatial configurations.
[0180] At trajectory point P i At, define z M Axis along the tool axis vector direction, x M The axis is along the tool feed direction, and the machining coordinate system {M} can be established according to the right-hand rule. Then the machining coordinate system {M} here is rotated around its z M Axis rotation redundancy angle δ R,i , we can get the tool coordinate system {C}, and the homogeneous transformation matrix of the tool coordinate system {C} relative to the robot base coordinate system {B} is
[0181] The homogeneous transformation matrix As the target matrix, according to the robot inverse kinematics, the robot joint angle configuration ξ corresponding to this condition is solved i The specific formula is:
[0182] ξ i =f IK (δ R,i )
[0183] Among them, f IK (δ R,i ) is based on the trajectory point P i The redundancy angle δ R,i Joint angle configuration ξ i The generalized function to be solved.
[0184] Then, according to the obtained joint angle combination ξ i , and further solve the end tool tip stiffness ellipsoid of the robot in the corresponding posture. First, the robot needs to be modeled for stiffness. The six joints of the robot can be regarded as linear torsion springs with fixed stiffness coefficients, and its joint stiffness matrix can be recorded as K ξ The specific formula is:
[0185]
[0186] in, These are the stiffnesses corresponding to the 6 rotational joints.
[0187] Through the Jacobian matrix J C (ξ), the robot's stiffness characteristics can be mapped from the joint angle space to the end Cartesian space. The specific formula is:
[0188]
[0189] Among them, K c is the tool tip stiffness matrix of the robot. In the internal sub-matrix subscripts, F and M represent the force and torque applied to the tool tip, and D and θ represent the translational displacement and rotational displacement generated by the tool tip. c is the Jacobian matrix at the tool tip of the robot.
[0190] In the dynamic analysis of robot milling, usually only the force on the end and the resulting translational displacement are considered. Therefore, the force / translational displacement submatrix K is taken as FD , solving its eigenvalue and eigenvector, we can get the stiffness ellipsoid equation of the robot tool tip in the corresponding posture, as follows: Figure 3 The specific formula is:
[0191]
[0192] Among them, σ1, σ2, σ3 are the force / translation displacement submatrix K FD The three eigenvalues sorted from large to small correspond to the lengths of the three semi-axes of the stiffness ellipsoid. The corresponding eigenvectors are r1, r2, and r3, and their directions correspond to the directions of the three semi-axes of the stiffness ellipsoid. Define the ellipsoid coordinate system {E}, where x E 、y E 、z E The axis direction is the direction of the three semi-axes, so the standard form of the tool tip stiffness ellipsoid equation can be obtained.
[0193] Then, the rotation matrix of the ellipsoidal coordinate system {E} relative to the base coordinate system {B} is Right now:
[0194]
[0195] From this, the rotation relationship between the ellipsoid coordinate system {E} and the machining coordinate system {M} can be written, and the calculation formula is:
[0196]
[0197] according to The tool tip stiffness ellipsoid equation can be converted from the ellipsoid coordinate system {E} to the machining coordinate system {M} for expression. The specific formula is:
[0198] P p=KM p
[0199] in, E p=[x E ,y E ,z E ] T 、 M p=[x M ,y M ,z M ] T are the coordinate vectors in the ellipsoid coordinate system {E} and the machining coordinate system {M} respectively.
[0200] Thus, the tool tip stiffness ellipsoid equation in the machining coordinate system {M} is obtained. The specific formula is:
[0201]
[0202] Let z M = 0, the elliptical intersection line of the stiffness ellipsoid on the milling plane can be obtained, such as Figure 4 The specific formula is:
[0203]
[0204] Combine similar terms of the above ellipse intersection equation. The specific formula is:
[0205]
[0206] Furthermore, an equivalent substitution is performed. The specific formula is:
[0207]
[0208] Then the above ellipse intersection equation can be rewritten into matrix form. The specific formula is:
[0209]
[0210] Obviously, the matrix For a real symmetric matrix, it can be decomposed into features. The specific formula is:
[0211]
[0212] Among them, λ1 and λ2 are matrices The corresponding eigenvalues, K is the matrix composed of its unit eigenvectors.
[0213] The ellipse intersection equation can be converted from the machining coordinate system {M} to the main stiffness coordinate system {P} through the matrix K. The specific formula is:
[0214] M pK=P p
[0215] in, M p=[x M ,y M ] T 、 P p=[x P ,y P ] T are the coordinate vectors in the ellipsoid coordinate system {M} and the machining coordinate system {E} respectively.
[0216] Thus, the equation of the ellipse intersection line in the main stiffness coordinate system {P} is obtained. The specific formula is:
[0217]
[0218] Obviously, in the principal stiffness coordinate system {P}, the intersection of the tool tip stiffness ellipsoid and the machining plane can be written as a standard form of ellipse equation. Its major and minor axes can quantitatively characterize the principal stiffness characteristics of the robot on the machining plane, such as Figure 5 shown.
[0219] Therefore, according to the major and minor axes of the ellipse, the corresponding principal stiffness difference index PSD can be calculated: i The specific formula is:
[0220]
[0221] Finally, the main stiffness difference index PSD at each processing trajectory point is defined n The sum is the chatter stability index CSI, and the maximum chatter stability index CSI is taken as the optimization target. i and the redundancy angle δ R,i Under these conditions, the corresponding main stiffness difference index PSD i The calculation method traverses all machining trajectory points and redundant angles. Under the constraints of physical interference, flexible space, joint speed, etc., the optimal redundant angle at each point is calculated, and then the optimal robot milling machining posture is obtained. The specific optimization model is as follows:
[0222]
[0223] Among them, [-90°, 90°] is the limit range of the robot's rotation around the tool axis vector under the interference of physical space after the milling tool is installed at the end of the robot; ξ Flexible It is a flexible workspace where the robot is far away from the singularity; is the rated load speed of the jth joint; μ is the set processing safety factor; ΔT is the processing time between adjacent milling trajectory points; ΔS is the interval distance between adjacent milling trajectory points; f is the set milling feed speed.
[0224] The milling posture optimization algorithm process used in the above optimization model is as follows: Figure 6 As shown, the specific calculation steps are:
[0225] S1: Perform equidistant dispersion on the machining trajectory to obtain the spatial coordinate information of each trajectory point.
[0226] S2: Initialize the trajectory point number i to 1.
[0227] S3: Initialize the maximum principal stiffness difference PSD max is 0.
[0228] S4: Calculate trajectory point P i The rotation matrix from the robot base coordinate system {B} to the machining coordinate system {M}
[0229] S5: Initialize the redundancy angle δ R,i -90° as the starting point of the traversal.
[0230] S6: According to the trajectory point P i The spatial coordinates [x i ,y i ,z i ] T and the redundancy angle δ R,i , calculate the homogeneous transformation matrix from the base coordinate system {B} to the tool coordinate system {C}
[0231] S7: According to the target matrix When solving inverse kinematics and there are multiple feasible solutions, we choose the solution in the flexible workspace ξ Flexible The inverse kinematic solution of ξ.
[0232] S8: According to the joint angle configuration ξ, solve the Jacobian matrix J of the tool tip corresponding to the spatial posture C (ξ).
[0233] S9: According to J C (ξ), solve the tool tip point overall stiffness matrix K c , and calculate its force / translation displacement submatrix K FD The eigenvalues σ1, σ2, σ2 are converted into eigenvectors r1, r2, r3, and the ellipsoid equation of the tool tip stiffness in the ellipsoid coordinate system {E} is obtained.
[0234] S10: Rotation matrix obtained according to S4 and the rotation matrix composed of r1, r2, and r3 obtained in S10 Solve Use Transform the stiffness ellipsoid equation to the machining coordinate system {M}.
[0235] S11: Let z M = 0, and solve to obtain the elliptic equation of the stiffness ellipsoid on the machining plane {x M Oy M}. According to the major and minor axes of the elliptic equation, solve the principal stiffness difference PSD.
[0236] S12: Make a determination. If PSD > PSD max and the joint angle configuration ξ meets the joint speed constraint, then PSD max = PSD and the optimal redundant angle δ R,i,opt = δ R,i . If not satisfied, directly execute S13.
[0237] S13: Make a determination. If δ R,i < 90°, then δ R,i = δ R,i + 1°, and return to S6. If not satisfied, the redundant angles at the trajectory point P i have been all traversed, and continue to execute S14.
[0238] S14: Make a determination. If i < N, then i = i + 1, and return to S3. If not satisfied, all the trajectory points have been calculated, and continue to execute S15.
[0239] S15: Output the finally obtained optimal redundant angle combination δ R,opt .
[0240] By determining the optimal redundant angle combination δ R,opt the optimal robot milling machining posture can be obtained accordingly.
[0241] To verify the actual effect of the present invention, a milling machining experiment is designed, and the experimental environment is as Figure 8 shown. The electric spindle system is connected to the robot end through a self-made fixture, and the milling cutter is clamped to the electric spindle. The acceleration sensor is pasted on the electric spindle fixture to collect the vibration acceleration signal of the end part during the machining process, and is transmitted to the computer through a data acquisition instrument for display and storage.
[0242] The machining trajectory of this milling experiment is equally discretized and divided, and the coordinates of the starting point in the base coordinate system {B} are determined. Input the milling machining task data into the aforementioned milling posture optimization algorithm, and the variation trend between the principal stiffness difference at the trajectory point and the corresponding redundant angle can be solved, as Figure 9 shown.
[0243] Determine the optimal redundancy angle δ at each processing trajectory point R,opt , and then get the optimal robot milling posture. Control the robot to perform milling in the optimal spatial posture.
[0244] When milling posture optimization is not performed, the following can be obtained: Figure 10 The experimental results are shown in Figure 2. The overall amplitude of the acceleration time-domain signal is large. Furthermore, in addition to the forced vibration frequency and its harmonic components caused by the periodic cutting of the cutter teeth, the spectrum also exhibits a significant low-frequency component that is very close to the natural frequency of the robot's structure. Therefore, it can be assumed that the milling system is experiencing low-frequency modal coupling chatter.
[0245] After milling posture optimization, the experimental results are as follows Figure 11 As shown in the figure, the overall amplitude of the acceleration time-domain signal is significantly lower than when no posture optimization is performed. Furthermore, the main components of the spectrum are concentrated at the cutter tooth frequency and its multiples, while the low-frequency components are essentially gone. Therefore, it can be concluded that under these conditions, the milling process is stable, with no modal coupling chatter occurring, and the proposed milling spatial posture optimization method is effective.
[0246] Although the embodiments of the present invention have been shown and described above, it will be understood that the above embodiments are illustrative and are not to be construed as limitations on the present invention. A person skilled in the art may change, modify, replace and modify the above embodiments within the scope of the present invention without departing from the principles and purpose of the present invention.
Claims
1. A robot milling posture optimization method based on main stiffness difference, characterized in that The specific steps are as follows: Divide the robot milling machining trajectory into equal distances and number each trajectory point; At each machining trajectory point, the tool tip stiffness ellipsoid of the robot is calculated under the posture corresponding to the specific redundant angle; Project the tool tip point stiffness ellipsoid onto the milling plane to obtain the robot's principal stiffness difference PSD; Taking the maximum sum of the principal stiffness differences (PSDs) corresponding to all machining trajectory points as the optimization goal, under the constraints of physical interference, flexible space, and joint speed, the optimal redundant angles of all trajectory points are traversed to solve the optimal robot milling posture.
2. The robot milling posture optimization method based on principal stiffness difference according to claim 1, characterized in that: The milling machining trajectory is obtained according to the actual milling machining task; the trajectory points are numbered as [P1, P2, ···, P i ,···,PN-1,PN]T.
3. The robot milling posture optimization method based on principal stiffness difference according to claim 2, characterized in that: The process of solving the tool tip stiffness ellipsoid of the robot under the specific redundant angle corresponding posture is as follows: According to the processing trajectory points, the corresponding redundant angle is determined; Determine the robot's spatial posture through coordinate transformation; According to the inverse kinematics of the robot, the joint angle configuration of the redundant angle corresponding to the machining trajectory point is solved; According to the joint angle configuration, solving the end tool tip point stiffness matrix of the robot in the corresponding posture; the end tool tip point stiffness matrix includes a force / translational displacement stiffness submatrix, a force / rotational displacement stiffness submatrix, a moment / translational displacement stiffness submatrix, and a moment / rotational displacement stiffness submatrix; According to the eigenvalues and eigenvectors of the force / translation displacement submatrix, the stiffness ellipsoid of the tool tip of the robot at the corresponding posture is solved.
4. The robot milling posture optimization method based on principal stiffness difference according to claim 3, characterized in that: The coordinate conversion process is to define the machining coordinate system in the robot milling process as the coordinate system {M}, where x M Axis along the tool feed direction, z M Axis along the direction of tool axis vector; the machining coordinate system {M} is rotated around z M The new coordinate system obtained after the axis rotates the redundant angle is defined as the tool coordinate system {C}. The tool coordinate system {C} is then used to reflect the orientation information of the tool axis vector and determine the spatial posture of the robot. The homogeneous transformation matrix of the tool coordinate system {C} relative to the robot base coordinate system {B} is recorded as B C T.
5. The robot milling posture optimization method based on principal stiffness difference according to claim 4, characterized in that: The calculation formula for the joint angle configuration of the redundant angle corresponding to the processing trajectory point is: x i =f IK (d R,i ) Among them, f IK (δ R,i ) is based on the trajectory point P i The redundancy angle δ R,i Joint angle configuration ξ i the generalized function to be solved; According to the joint angle configuration, the stiffness matrix of the robot's end tool tip in the corresponding posture is solved, and the formula is as follows: Among them, K c K is the tool tip stiffness matrix of the robot; FD is the force / translation displacement stiffness submatrix; K Fθ is the force / rotational displacement stiffness submatrix; K MD is the moment / translation displacement stiffness submatrix; K Mθ is the moment / rotational displacement stiffness submatrix; in the subscripts of the internal submatrix, F and M represent the force and moment applied to the tool tip, and D and θ represent the translational displacement and rotational displacement generated by the tool tip; J c is the Jacobian matrix at the robot tool tip, K ξ is the joint stiffness matrix of the robot; According to the force / translation displacement submatrix K FD The eigenvalue and eigenvector of the robot are used to obtain the stiffness ellipsoid equation of the tool tip at the end of the robot under the corresponding posture. The calculation formula is: Among them, σ1, σ2, σ3 are the matrix K FD The three eigenvalues sorted from large to small correspond to the lengths of the three semi-axes of the tool tip stiffness ellipsoid; the corresponding unit eigenvectors are r1, r2, and r3, corresponding to the directions of the three semi-axes of the tool tip stiffness ellipsoid; define the ellipsoid coordinate system {E}, whose x E 、y E 、z E The axis direction is the direction of the three semi-axes, and the ellipsoid equation has a standard form in the ellipsoidal coordinate system {E}.
6. The robot milling posture optimization method based on principal stiffness difference according to claim 5, characterized in that: The process of projecting the tool tip stiffness ellipsoid onto the milling plane to obtain the corresponding ellipse intersection equation is as follows: The rotation matrix of the ellipsoidal coordinate system {E} relative to the base coordinate system {B} Denoted as: The rotation matrix of the machining coordinate system {M} relative to the base coordinate system {B} Expressed as: The relationship between the ellipsoid coordinate system {E} and the machining coordinate system {M} is obtained as follows: according to The tool tip stiffness ellipsoid equation is converted from the ellipsoid coordinate system {E} to the machining coordinate system {M}. The calculation formula is: in, E p=[x E ,y E ,z E ] T 、 M p=[x M ,y M ,z M ] T are the coordinate vectors in the ellipsoid coordinate system {E} and the machining coordinate system {M} respectively; Among them, k ij is the rotation matrix Each element in , i = 1, 2, 3; j = 1, 2, 3); The relationship between the coordinate vectors in the ellipsoid coordinate system {E} and the machining coordinate system {M} is expressed as: The standard equation of the tool tip stiffness ellipsoid is converted from the ellipsoid coordinate system {E} to the machining coordinate system {M} for expression: Let z M =0, the elliptical intersection equation of the rigidity ellipsoid in the milling plane is obtained. The specific formula is:
7. The robot milling posture optimization method based on principal stiffness difference according to claim 6, characterized in that: The calculation process of the robot's main stiffness difference PSD is as follows: By combining like terms and making equivalent substitutions, the ellipse intersection equation can be written in matrix form as follows: make: Transform For the real symmetric matrix Perform feature decomposition and the calculation formula is: Among them, λ1 and λ2 are matrices The corresponding eigenvalue, K is the matrix composed of its unit eigenvector; Through the matrix K, the coordinate vector [x M ,y M ] T By performing a rotation transformation, the ellipse intersection equation can be converted from the machining coordinate system {M} to the principal stiffness coordinate system {P}; the specific formula is as follows: P p=K M p After expansion: in, M p=[x M ,y M ] T 、 P p=[x P ,y P ] T are the coordinate vectors in the machining coordinate system {M} and the principal stiffness coordinate system {P} respectively; In the principal stiffness coordinate system {P}, the intersection of the tool tip stiffness ellipsoid and the machining plane is written as a standard form of the ellipse equation. The principal stiffness difference PSD of the robot on the machining plane is quantified by the difference between the major and minor axes. The formula is as follows: in, 8. The robot milling posture optimization method based on principal stiffness difference according to claim 7, characterized in that: The sum of the principal stiffness differences PSD corresponding to all machining trajectory points is defined as the chatter stability index CSI, and the maximum chatter stability index CSI is used as the optimization target.
9. The robot milling posture optimization method based on principal stiffness difference according to claim 8, characterized in that: The optimization model of the optimal robot milling posture is: Among them, [-90°, 90°] is the limit range of the robot's rotation around the tool axis vector under the interference of physical space after the milling tool is installed at the end of the robot; ξ Flexible It is a flexible workspace where the robot is far away from the singularity; is the rated load speed of the jth joint; μ is the set processing safety factor; ΔT is the processing time between adjacent milling trajectory points; ΔS is the interval distance between adjacent milling trajectory points; f is the set milling feed speed.
10. A robot milling posture optimization system based on primary stiffness difference, characterized by: The method comprises at least one processor and a memory communicatively connected to the at least one processor; wherein the memory stores a computer program executable by the at least one processor, and the computer program is executed by the at least one processor so that the at least one processor can execute the robot milling posture optimization method based on the main stiffness difference as described in any one of claims 1 to 9.
Citation Information
Patent Citations
Milling three-dimensional stability forecasting method of six-freedom-degree series robot
CN108638076A
Elastic deformation and vibration inhibition method for grinding and polishing machining of airplane composite material component by robot
CN111673611A